## ----style, eval=TRUE, echo=FALSE, results='asis'------------------------ options(max.print=1000) stopifnot(BiocInstaller::biocVersion() == "2.13") BiocStyle::latex() library(knitr) opts_chunk$set(cache=TRUE, tidy=FALSE) ## ----results='hide'------------------------------------------------------ library(parathyroidSE) ## ----eval=FALSE---------------------------------------------------------- ## data(package="parathyroidSE") ## ------------------------------------------------------------------------ data(parathyroidGenesSE) parathyroidGenesSE class(parathyroidGenesSE) dim(parathyroidGenesSE) ## ------------------------------------------------------------------------ class(assay(parathyroidGenesSE)) dim(assay(parathyroidGenesSE)) assay(parathyroidGenesSE)[1:8, 1:8] ## ------------------------------------------------------------------------ colnames(parathyroidGenesSE) ## ------------------------------------------------------------------------ colSums(assay(parathyroidGenesSE)) ## ------------------------------------------------------------------------ sum(rowSums(assay(parathyroidGenesSE)) != 0) ## ------------------------------------------------------------------------ rowData(parathyroidGenesSE) ## ------------------------------------------------------------------------ colData(parathyroidGenesSE) ## ------------------------------------------------------------------------ table(colData(parathyroidGenesSE)$treatment) ## ----bampaths------------------------------------------------------------ bamdir <- system.file("extdata", package="parathyroidSE") bampaths <- list.files(bamdir, pattern="bam$", full.names=TRUE) bampaths ## ------------------------------------------------------------------------ library(Rsamtools) gal0 <- readGAlignmentsFromBam(bampaths[1]) gal0 ## ------------------------------------------------------------------------ galp0 <- readGAlignmentPairsFromBam(bampaths[1]) galp0 ## ------------------------------------------------------------------------ length(galp0) ## ------------------------------------------------------------------------ first(galp0) last(galp0) ## ------------------------------------------------------------------------ param <- ScanBamParam(flag=scanBamFlag(isDuplicate=FALSE, isNotPrimaryRead=FALSE)) readGAlignmentPairsFromBam(bampaths[1], param=param) ## ----makeTranscriptDbFromBiomart, eval=FALSE----------------------------- ## ## Requires INTERNET ACCESS and takes about 6 min. Please don't try to run this! ## library(GenomicFeatures) ## txdb <- makeTranscriptDbFromBiomart(biomart="ensembl", ## dataset="hsapiens_gene_ensembl") ## ----exonsBy, eval=FALSE------------------------------------------------- ## ex_by_gene <- exonsBy(txdb, by="gene") # GRangesList object ## ------------------------------------------------------------------------ data(exonsByGene) exonsByGene ## ------------------------------------------------------------------------ length(exonsByGene) ## ----BamFileList--------------------------------------------------------- library(Rsamtools) bamfile_list <- BamFileList(bampaths, index=character()) ## ------------------------------------------------------------------------ read_count0 <- summarizeOverlaps(exonsByGene, bamfile_list, ignore.strand=TRUE, singleEnd=FALSE) read_count0 ## ------------------------------------------------------------------------ read_count <- summarizeOverlaps(exonsByGene, bamfile_list, ignore.strand=TRUE, singleEnd=FALSE, param=param) read_count ## ------------------------------------------------------------------------ colSums(assay(read_count0)) colSums(assay(read_count))