R/createAnnotationPackage.hg19.R

In AdrianS85/RnBeadsAnnotationCreator_rn6: RnBeads annotation creator

########################################################################################################################
## createAnnotationPackage.hg19.R
## created: 2014-02-13
## creator: Yassen Assenov
## ---------------------------------------------------------------------------------------------------------------------
## Annotation package creation for HG19.
########################################################################################################################

## F U N C T I O N S ###################################################################################################

#' createAnnotationPackage.hg19
#'
#' Helper function to create RnBeads annotation package for genome assembly hg19.
#'
#' @return None (invisible \code{NULL}).
#' @author Fabian Mueller
#' @noRd
createAnnotationPackage.hg19 <- function() {

	suppressPackageStartupMessages(library(BSgenome.Hsapiens.UCSC.hg19))
	suppressPackageStartupMessages(library(FDb.InfiniumMethylation.hg19))

	## Genomic sequence and supported chromosomes
	GENOME <- 'BSgenome.Hsapiens.UCSC.hg19'
	assign('GENOME', GENOME, .globals)
	CHROMOSOMES <- c(1:22, "X", "Y")
	names(CHROMOSOMES) <- paste0("chr", CHROMOSOMES)
	assign('CHROMOSOMES', CHROMOSOMES, .globals)
	rm(GENOME, CHROMOSOMES)

	## Download SNP annotation
	logger.start("SNP Annotation")
	vcf.files <- paste0(DBSNP.FTP.BASE, "human_9606_b147_GRCh37p13/VCF/00-All.vcf.gz")
	update.annot("snps", "polymorphism information", rnb.update.dbsnp, ftp.files = vcf.files)
	logger.info(paste("Using:", attr(.globals[['snps']], "version")))
	rm(vcf.files)
	logger.completed()

	## Define genomic regions
	biomart.parameters <- list(
		database.name = "ENSEMBL_MART_ENSEMBL",
		dataset.name = "hsapiens_gene_ensembl",
		required.columns = c(
			"id" = "ensembl_gene_id",
			"chromosome" = "chromosome_name",
			"start" = "start_position",
			"end" = "end_position",
			"strand" = "strand",
			"symbol" = "hgnc_symbol",
			"entrezID" = "entrezgene"),
		host = "feb2014.archive.ensembl.org")
	logger.start("Region Annotation")
	update.annot("regions", "region annotation", rnb.update.region.annotation,
		biomart.parameters = biomart.parameters)
	rm(biomart.parameters)
	logger.completed()

	## Define genomic sites
	logger.start("Genomic Sites")
	update.annot("sites", "CpG annotation", rnb.update.sites)
	logger.completed()

	## Define Infinium 27k probe annotations
	logger.start("Infinium 27k")
	ftp.table <- paste0(GEO.FTP.BASE, "GPL8490/GPL8490_HumanMethylation27_270596_v.1.2.csv.gz")
	table.columns <- rnb.get.illumina.annotation.columns("27k")
	update.annot("probes27", "Infinium 27K annotation", rnb.update.probe27k.annotation,
		ftp.table = ftp.table, table.columns = table.columns)
	.globals[['sites']][['probes27']] <- .globals[['probes27']][['probes']]
	logger.completed()

	## Define Infinium 450k probe annotations
	logger.start("Infinium 450k")
	ftp.table <- paste0(GEO.FTP.BASE, "GPL13534/GPL13534_HumanMethylation450_15017482_v.1.1.csv.gz")
	table.columns <- rnb.get.illumina.annotation.columns("450k")
	update.annot("probes450", "Infinium 450K annotation", rnb.update.probe450k.annotation,
		ftp.table = ftp.table, table.columns = table.columns)
	.globals[['sites']][["probes450"]] <- .globals[['probes450']][["probes"]]
	rm(ftp.table, table.columns)
	logger.completed()

	## Define MethylationEPIC probe annotations
	logger.start("MethylationEPIC")
	table.columns <- rnb.get.illumina.annotation.columns("EPIC")
	update.annot("probesEPIC", "MethylationEPIC annotation", rnb.update.probeEPIC.annotation,
		table.columns = table.columns)
	.globals[['sites']][["probesEPIC"]] <- .globals[['probesEPIC']][["probes"]]
	logger.completed()

	## Add annotation columns to the context probes, showing if they are covered by an assay
	logger.start("Updating Site Annotation with Probes")
	.globals[['sites']] <- rnb.update.site.annotation.with.probes(sites = .globals[['sites']],
		query.probes = c("probes27", "probes450", "probesEPIC"),
		platform.names = c("HumanMethylation27", "HumanMethylation450", "MethylationEPIC"))
	logger.completed()

	## Create all possible mappings from regions to sites
	logger.start("Mappings")
	update.annot("mappings", "mappings", rnb.create.mappings)
	logger.completed()

	## Export the annotation tables
	rnb.export.annotations.to.data.files()
}