### R code from vignette source 'TarSeqQC-vignette.Rnw' ################################################### ### code chunk number 1: General options for R ################################################### options(prompt="> ", continue="+ ", width=78, useFancyQuotes=FALSE, digits=3) ################################################### ### code chunk number 2: TarSeqQC-vignette.Rnw:248-249 ################################################### bedFile<-system.file("extdata", "mybed.bed", package="TarSeqQC", mustWork=TRUE) ################################################### ### code chunk number 3: TarSeqQC-vignette.Rnw:272-273 ################################################### bamFile<-system.file("extdata", "mybam.bam", package="TarSeqQC", mustWork=TRUE) ################################################### ### code chunk number 4: TarSeqQC-vignette.Rnw:286-288 ################################################### fastaFile<-system.file("extdata", "myfasta.fa", package="TarSeqQC", mustWork=TRUE) ################################################### ### code chunk number 5: library call ################################################### suppressMessages(library("TarSeqQC")) suppressMessages(library("BiocParallel")) ################################################### ### code chunk number 6: constructor (eval = FALSE) ################################################### ## library("TarSeqQC") ## library("BiocParallel") ## BPPARAM<-bpparam() ## myPanel<-TargetExperiment(bedFile, bamFile, fastaFile, feature="amplicon", ## attribute="coverage", BPPARAM=BPPARAM) ################################################### ### code chunk number 7: loading TargetExperiment object ################################################### BPPARAM<-bpparam() data(ampliPanel, package="TarSeqQC") myPanel<-ampliPanel rm(ampliPanel) setBamFile(myPanel)<-system.file("extdata", "mybam.bam", package="TarSeqQC", mustWork=TRUE) setFastaFile(myPanel)<-system.file("extdata", "myfasta.fa", package="TarSeqQC", mustWork=TRUE) ################################################### ### code chunk number 8: setters ################################################### # set feature slot value setFeature(myPanel)<-"amplicon" # set attribute slot value setAttribute(myPanel)<-"coverage" ################################################### ### code chunk number 9: setters (eval = FALSE) ################################################### ## # set scanBamP slot value ## scanBamP<-ScanBamParam() ## #set scanBamP which slot ## bamWhich(scanBamP)<-getBedFile(myPanel) ## setScanBamP(myPanel)<-scanBamP ## # set pileupP slot value ## setPileupP(myPanel)<-PileupParam(max_depth=1000) ## # build the featurePanel again ## setFeaturePanel(myPanel)<-buildFeaturePanel(myPanel, BPPARAM) ## # build the genePanel again ## setGenePanel(myPanel)<-summarizePanel(myPanel, BPPARAM) ################################################### ### code chunk number 10: loading TargetExperiment object ################################################### data(ampliPanel, package="TarSeqQC") ################################################### ### code chunk number 11: re-defininf file paths ################################################### # Defining bam file and fasta file names and paths setBamFile(ampliPanel)<-system.file("extdata", "mybam.bam", package="TarSeqQC", mustWork=TRUE) setFastaFile(ampliPanel)<-system.file("extdata", "myfasta.fa", package="TarSeqQC", mustWork=TRUE) ################################################### ### code chunk number 12: exploration ################################################### # show/print myPanel # summary summary(myPanel) #summary at feature level summaryFeatureLev(myPanel) #summary at gene level summaryGeneLev(myPanel) ################################################### ### code chunk number 13: boxplot and density plot (eval = FALSE) ################################################### ## g<-plotAttrExpl(myPanel,level="feature",join=TRUE, log=FALSE, color="blue") ## x11(type="cairo"); ## g ################################################### ### code chunk number 14: length exploration (eval = FALSE) ################################################### ## # explore amplicon length distribution ## plotMetaDataExpl(ampliPanel, "length", log=FALSE, join=FALSE, color= ## "blueviolet") ################################################### ### code chunk number 15: gene exploration (eval = FALSE) ################################################### ## # explore gene's relative frequencies ## plotMetaDataExpl(ampliPanel, "gene", abs=FALSE) ################################################### ### code chunk number 16: read percentages and plot (eval = FALSE) ################################################### ## readFrequencies(ampliPanel) ## plotInOutFeatures(readFrequencies(ampliPanel)) ################################################### ### code chunk number 17: attributeThres ################################################### # definition of the interval extreme values attributeThres<-c(0,1,50,200,500, Inf) ################################################### ### code chunk number 18: plot (eval = FALSE) ################################################### ## # plot panel overview ## g<-plot(myPanel, attributeThres, chrLabels =TRUE) ## g ################################################### ### code chunk number 19: plotFeatPerform (eval = FALSE) ################################################### ## # plot panel overview ## g<-plotFeatPerform(myPanel, attributeThres, complete=TRUE, log=FALSE, ## featureLabs=TRUE, sepChr=TRUE, legend=TRUE) ## g ################################################### ### code chunk number 20: gc exploration (eval = FALSE) ################################################### ## x11(type="cairo") ## biasExploration(myPanel, source="gc", dens=TRUE) ################################################### ### code chunk number 21: frequency table ################################################### # summaryIntervals summaryIntervals(myPanel, attributeThres) ################################################### ### code chunk number 22: coveragePerform (eval = FALSE) ################################################### ## plotAttrPerform(myPanel, attributeThres) ################################################### ### code chunk number 23: low counts features ################################################### getLowCtsFeatures(myPanel, level="gene", threshold=50) ################################################### ### code chunk number 24: low counts features ################################################### getLowCtsFeatures(myPanel, level="feature", threshold=50) ################################################### ### code chunk number 25: geneAttrPerFeat (eval = FALSE) ################################################### ## g<-plotGeneAttrPerFeat(myPanel, geneID="gene4") ## # adjust text size ## g<-g+theme(title=element_text(size=16), axis.title=element_text(size=16), ## legend.text=element_text(size=14)) ## g ################################################### ### code chunk number 26: pileupCounts (eval = FALSE) ################################################### ## # define function parameters ## bed<-getBedFile(myPanel) ## bamFile<-system.file("extdata", "mybam.bam", package="TarSeqQC", mustWork=TRUE) ## fastaFile<-system.file("extdata", "myfasta.fa", package="TarSeqQC", ## mustWork=TRUE) ## scanBamP<-getScanBamP(myPanel) ## pileupP<-getPileupP(myPanel) ## #call pileupCounts function ## myCounts<-pileupCounts(bed=bed, bamFile=bamFile, fastaFile=fastaFile, ## scanBamP=scanBamP, pileupP=pileupP, BPPARAM=BPPARAM) ################################################### ### code chunk number 27: pileupCounts ################################################### data("myCounts", package="TarSeqQC") ################################################### ### code chunk number 28: pileupCounts ################################################### head(myCounts) ################################################### ### code chunk number 29: getRegion ################################################### #complete information for gene7 getRegion(myPanel, level="gene", ID="gene7", collapse=FALSE) #summarized information for gene7 getRegion(myPanel, level="gene", ID="gene7", collapse=TRUE) ################################################### ### code chunk number 30: plotRegion (eval = FALSE) ################################################### ## g<-plotRegion(myPanel, region=c(4500,6800), seqname="chr10", SNPs=TRUE, ## xlab="", title="gene7 amplicons",size=0.5) ## x11(type="cairo") ## g ################################################### ### code chunk number 31: plotFeature (eval = FALSE) ################################################### ## g<-plotFeature(myPanel, featureID="AMPL20") ## x11(type="cairo") ## g ################################################### ### code chunk number 32: plotNtd (eval = FALSE) ################################################### ## g<-plotNtdPercentage(myPanel, featureID="AMPL20") ## g ################################################### ### code chunk number 33: featureInfo ################################################### getFeaturePanel(myPanel)["AMPL20"] ################################################### ### code chunk number 34: readCountsFeat< ################################################### featureCounts<-myCounts[myCounts[, "seqnames"] =="chr10" & myCounts[,"pos"] >= 4866 & myCounts[,"pos"] <= 4928,] ################################################### ### code chunk number 35: readCountsFeat< ################################################### featureCounts[which.min(featureCounts[,"="]),] ################################################### ### code chunk number 36: buildReport (eval = FALSE) ################################################### ## imageFile<-system.file("extdata", "plot.pdf", package="TarSeqQC", ## mustWork=TRUE) ## buildReport(ampliPanel, attributeThres, imageFile ,file="Results.xlsx") ################################################### ### code chunk number 37: loading TarSeqQC example data ################################################### data(ampliPanel, package="TarSeqQC") ampliPanel ################################################### ### code chunk number 38: using TarSeqQC example data (eval = FALSE) ################################################### ## ## setBamFile(ampliPanel)<-system.file("extdata", "mybam.bam", package="TarSeqQC", ## mustWork=TRUE) ## setFastaFile(ampliPanel)<-system.file("extdata", "myfasta.fa", ## package="TarSeqQC", mustWork=TRUE) ################################################### ### code chunk number 39: buildFeaturePanel using TarSeqQC example data (eval = FALSE) ################################################### ## plotFeature(ampliPanel, featureID="AMPL1") ################################################### ### code chunk number 40: Session Info ################################################### sessionInfo()