### R code from vignette source 'TitanCNA.Rnw' ### Encoding: UTF-8 ################################################### ### code chunk number 1: TitanCNA.Rnw:75-80 ################################################### library(TitanCNA) infile <- system.file("extdata", "test_alleleCounts_chr2.txt", package = "TitanCNA") data <- loadAlleleCounts(infile, genomeStyle = "NCBI") names(data) ################################################### ### code chunk number 2: TitanCNA.Rnw:85-89 ################################################### data$chr <- setGenomeStyle(data$chr, genomeStyle = "UCSC") head(data$chr) data$chr <- setGenomeStyle(data$chr, genomeStyle = "NCBI") head(data$chr) ################################################### ### code chunk number 3: TitanCNA.Rnw:96-102 ################################################### tumWig <- system.file("extdata", "test_tum_chr2.wig", package = "TitanCNA") normWig <- system.file("extdata", "test_norm_chr2.wig", package = "TitanCNA") gc <- system.file("extdata", "gc_chr2.wig", package = "TitanCNA") map <- system.file("extdata", "map_chr2.wig", package = "TitanCNA") cnData <- correctReadDepth(tumWig, normWig, gc, map, genomeStyle = "NCBI") head(cnData) ################################################### ### code chunk number 4: TitanCNA.Rnw:109-111 ################################################### logR <- getPositionOverlap(data$chr, data$posn, cnData) data$logR <- log(2^logR) #transform the log ratio to natural logs ################################################### ### code chunk number 5: TitanCNA.Rnw:121-123 ################################################### data <- filterData(data, c(1:22, "X", "Y"), minDepth = 10, maxDepth = 200, positionList = NULL, centromere = NULL, centromere.flankLength = 10000) ################################################### ### code chunk number 6: TitanCNA.Rnw:131-137 ################################################### numClusters <- 2 params <- loadDefaultParameters(copyNumber = 5, numberClonalClusters = numClusters, symmetric = TRUE, hetBaselineSkew = 0, alleleEmissionModel = "binomial", data = data) params ################################################### ### code chunk number 7: TitanCNA.Rnw:151-153 (eval = FALSE) ################################################### ## params$ploidyParams$phi_0 <- 2 # for diploid or ## params$ploidyParams$phi_0 <- 4 # for tetraploid/ployploid ################################################### ### code chunk number 8: TitanCNA.Rnw:158-161 ################################################### K <- length(params$genotypeParams$alphaKHyper) params$genotypeParams$alphaKHyper <- rep(500, K) params$ploidyParams$phi_0 <- 1.5 ################################################### ### code chunk number 9: TitanCNA.Rnw:169-171 (eval = FALSE) ################################################### ## library(doMC) ## registerDoMC(cores = 4) #use 4 cores on a single machine ################################################### ### code chunk number 10: TitanCNA.Rnw:186-193 ################################################### convergeParams <- runEMclonalCN(data, params, maxiter = 3, maxiterUpdate = 50, useOutlierState = FALSE, txnExpLen = 1e15, txnZstrength = 5e5, normalEstimateMethod = "map", estimateS = TRUE, estimatePloidy = TRUE) names(convergeParams) ################################################### ### code chunk number 11: TitanCNA.Rnw:213-215 ################################################### optimalPath <- viterbiClonalCN(data, convergeParams) head(optimalPath) ################################################### ### code chunk number 12: TitanCNA.Rnw:226-236 ################################################### results <- outputTitanResults(data, convergeParams, optimalPath, filename = NULL, posteriorProbs = FALSE, subcloneProfiles = TRUE, correctResults = TRUE, proportionThreshold = 0.05, proportionThresholdClonal = 0.05, is.haplotypeData = FALSE) names(results) head(results$corrResults) ## corrected results convergeParams <- results$convergeParam ## use corrected parameters results <- results$corrResults ## use corrected results ################################################### ### code chunk number 13: eval ################################################### segs <- outputTitanSegments(results, id = "test", convergeParams, filename = NULL, igvfilename = NULL) head(segs) ################################################### ### code chunk number 14: TitanCNA.Rnw:302-304 (eval = FALSE) ################################################### ## outparam <- paste("test_cluster02_params.txt", sep = "") ## outputModelParameters(convergeParams, results, outparam, S_Dbw.scale = 1) ################################################### ### code chunk number 15: TitanCNA.Rnw:324-325 (eval = FALSE) ################################################### ## outputModelParameters(convergeParams, results, outparam, S_Dbw.scale = 10) ################################################### ### code chunk number 16: TitanCNA.Rnw:335-341 ################################################### ploidy <- tail(convergeParams$phi, 1) ploidy normal <- tail(convergeParams$n, 1) normal plotCNlogRByChr(results, segs = segs, chr = 2, ploidy = ploidy, normal = normal, ylim = c(-2, 2), cex = 0.25, xlab = "", main = "Chr 2") ################################################### ### code chunk number 17: TitanCNA.Rnw:354-356 ################################################### plotAllelicRatio(results, chr = 2, ylim = c(0, 1), cex = 0.25, xlab = "", main = "Chr 2") ################################################### ### code chunk number 18: TitanCNA.Rnw:370-375 ################################################### norm <- tail(convergeParams$n, 1) norm # estimated normal contamination 1 - convergeParams$s[, ncol(convergeParams$s)] # estimated cellular prevalence plotClonalFrequency(results, chr = 2, normal = norm, ylim = c(0, 1), cex = 0.25, xlab = "", main = "Chr 2") ################################################### ### code chunk number 19: TitanCNA.Rnw:384-385 ################################################### plotSubcloneProfiles(results, chr = 2, cex = 1, spacing = 2, main = "Chr 2") ################################################### ### code chunk number 20: TitanCNA.Rnw:393-395 ################################################### plotSegmentMedians(segs, chr=2, resultType = "LogRatio", plotType = "CopyNumber", plot.new = TRUE, ylim = c(0, 4), main="Chr 2") ################################################### ### code chunk number 21: TitanCNA.Rnw:402-403 ################################################### toLatex(sessionInfo())