################################################### ### chunk number 1: load library ################################################### #line 52 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" library(lumi) ################################################### ### chunk number 2: load example dataset ################################################### #line 74 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## load example data (a methyLumiM object) data(example.lumiMethy) ## summary of the example data example.lumiMethy ## print sample Names sampleNames(example.lumiMethy) ################################################### ### chunk number 3: sampleRelation ################################################### #line 87 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotSampleRelation(example.lumiMethy, method='mds', cv.Th=0) ################################################### ### chunk number 4: sampleRelationTree ################################################### #line 96 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotSampleRelation(example.lumiMethy, method='cluster', cv.Th=0) ################################################### ### chunk number 5: load example titration dataset ################################################### #line 107 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## load the tritration data (a methyLumiM object) data(example.methyTitration) ## summary of the example data example.methyTitration ## print sample Names sampleNames(example.methyTitration) ################################################### ### chunk number 6: sampleRelationTitration ################################################### #line 118 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotSampleRelation(example.methyTitration, method='mds', cv.Th=0) ################################################### ### chunk number 7: densityMTitration ################################################### #line 135 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot the density density(example.methyTitration, xlab="M-value") ################################################### ### chunk number 8: densityM ################################################### #line 146 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## specify the colors of control and treatment samples sampleColor <- rep(1, ncol(example.lumiMethy)) sampleColor[grep("Treat", sampleNames(example.lumiMethy))] <- 2 density(example.lumiMethy, col=sampleColor, xlab="M-value") ## plot the density ################################################### ### chunk number 9: boxplotM ################################################### #line 160 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## Because the distribution of M-value has two modes, we use a boxplot different from regular ones boxplot(example.lumiMethy) ################################################### ### chunk number 10: densityColorBiasBoth ################################################### #line 178 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotColorBias1D(example.lumiMethy) ################################################### ### chunk number 11: densityColorBiasMethy ################################################### #line 187 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotColorBias1D(example.lumiMethy, channel='methy') ################################################### ### chunk number 12: densityColorBiasUnmethy ################################################### #line 196 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotColorBias1D(example.lumiMethy, channel='unmethy') ################################################### ### chunk number 13: boxplotColorBiasMethy ################################################### #line 208 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" boxplotColorBias(example.lumiMethy, channel='methy') ################################################### ### chunk number 14: boxplotColorBiasUnmethy ################################################### #line 218 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" boxplotColorBias(example.lumiMethy, channel='unmethy') ################################################### ### chunk number 15: color balance summary ################################################### #line 225 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## summary of color balance information of individual samples colorBiasSummary(example.lumiMethy[,1:8], channel='methy') ################################################### ### chunk number 16: scatterColorBias1 ################################################### #line 236 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotColorBias2D(example.lumiMethy, selSample=1, cex=2) ################################################### ### chunk number 17: boxplotColorBiasSum ################################################### #line 252 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" boxplotColorBias(example.lumiMethy, channel='sum') ################################################### ### chunk number 18: densityColorBiasSum ################################################### #line 263 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotColorBias1D(example.lumiMethy, channel='sum') ################################################### ### chunk number 19: densityIntensity ################################################### #line 274 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" density(estimateIntensity(example.lumiMethy), xlab="log2(CpG-site Intensity)") ################################################### ### chunk number 20: boxplotIntensity ################################################### #line 284 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" boxplot(estimateIntensity(example.lumiMethy)) ################################################### ### chunk number 21: pairsColor ################################################### #line 297 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## get the color channel information colorChannel <- as.character(pData(featureData(example.lumiMethy))[, "COLOR_CHANNEL"]) ## replace the "Red" and "Grn" as color names defined in R colorChannel[colorChannel == 'Red'] <- 'red' colorChannel[colorChannel == 'Grn'] <- 'green' ## select a subet of sample for pair plot selSample <- c( "Ctrl1", "Ctrl1.rep", "Treat1", "Treat1.rep") ## plot pair plot with the dots in scatter plot colored based on the color channels pairs(estimateIntensity(example.lumiMethy[, selSample]), dotColor= colorChannel, main="Pair plot of CpG-site Intensity") ################################################### ### chunk number 22: color balance adjustment ################################################### #line 321 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## summary of color balance information of individual samples lumiMethy.c.adj <- lumiMethyC(example.lumiMethy) ################################################### ### chunk number 23: densityColorBiasSumAdj ################################################### #line 328 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotColorBias1D(lumiMethy.c.adj, channel='sum') ################################################### ### chunk number 24: boxplotColorBiasSumAdj ################################################### #line 338 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" boxplotColorBias(lumiMethy.c.adj, channel='sum') ################################################### ### chunk number 25: scatterColorBias1Adj ################################################### #line 348 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot the color balance adjusted scatter plot of two color channels plotColorBias2D(lumiMethy.c.adj, selSample=1, cex=2) ################################################### ### chunk number 26: pairsColorAdj ################################################### #line 359 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot pairwise plot after color balance adjustment pairs(estimateIntensity(lumiMethy.c.adj[, selSample]), dotColor= colorChannel, main="Pair plot of CpG-site Intensity") ################################################### ### chunk number 27: background adjustment ################################################### #line 374 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ##separately adjust backgrounds of two color channels lumiMethy.b.adj <- lumiMethyB(example.lumiMethy, method="bgAdjust2C", separateColor=TRUE) ##background adjustment of individual samples lumiMethy.bc.adj <- lumiMethyB(lumiMethy.c.adj, method="bgAdjust2C") ################################################### ### chunk number 28: bgDensityMethy ################################################### #line 385 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot the background mode of methylated probe data of first five example samples plotColorBias1D(example.lumiMethy[,1:5], channel='methy', xlim=c(-1000,5000), logMode=FALSE) ################################################### ### chunk number 29: bgAdjDensityMethy ################################################### #line 396 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot the background mode of methylated probe data of first five example samples plotColorBias1D(lumiMethy.b.adj [,1:5], channel='methy', xlim=c(-1000,5000), logMode=FALSE) ################################################### ### chunk number 30: bcAdjDensityMethy ################################################### #line 407 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot the background mode of methylated probe data of first five example samples plotColorBias1D(lumiMethy.bc.adj [,1:5], channel='methy', xlim=c(-1000,5000), logMode=FALSE) ################################################### ### chunk number 31: Normalization ################################################### #line 427 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## Perform SSN normalization based on color balance adjusted data lumiMethy.c.ssn <- lumiMethyN(lumiMethy.c.adj, method='ssn') ## Perform quantile normalization based on color balance adjusted data lumiMethy.c.quantile <- lumiMethyN(lumiMethy.c.adj, method='quantile') ################################################### ### chunk number 32: sampleRelationTreeNormalized ################################################### #line 437 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" plotSampleRelation(lumiMethy.c.ssn, method='cluster', cv.Th=0) ################################################### ### chunk number 33: densitySSN ################################################### #line 447 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot the density of M-values after SSN normalization density(lumiMethy.c.ssn, col= sampleColor, main="Density plot after SSN normalization") ################################################### ### chunk number 34: densityQuantile ################################################### #line 458 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## plot the density of M-values after quantile normalization density(lumiMethy.c.quantile, col= sampleColor, main="Density plot after quantile normalization") ################################################### ### chunk number 35: densityIntensityNormalizedSSN ################################################### #line 469 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" density(estimateIntensity(lumiMethy.c.ssn), col= sampleColor, xlab="log2(CpG-site Intensity)") ################################################### ### chunk number 36: densityIntensityNormalizedQ ################################################### #line 479 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" density(estimateIntensity(lumiMethy.c.quantile), col= sampleColor, xlab="log2(CpG-site Intensity)") ################################################### ### chunk number 37: boxplotColorBiasNormalized ################################################### #line 489 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" boxplotColorBias(lumiMethy.c.quantile, channel='sum') ################################################### ### chunk number 38: pairMNormalize ################################################### #line 499 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## select a subet of sample for pair plot selSample <- c( "Ctrl1", "Ctrl1.rep", "Treat1", "Treat1.rep") ## plot pair plot with the dots in scatter plot colored based on the color channels pairs(lumiMethy.c.quantile[, selSample], dotColor= colorChannel, main="Pair plot of M-value after normalization") ################################################### ### chunk number 39: gammaFit ################################################### #line 515 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## Fit the two component gamma mixture model of the first sample of example.lumiMethy fittedGamma <- gammaFitEM(exprs(example.lumiMethy)[,1], plotMode=TRUE) ################################################### ### chunk number 40: estimate methylation status based on M-value ################################################### #line 526 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## estimate the methylation status based on the results of gammaFitEM methyCall <- methylationCall(fittedGamma) table(methyCall) ################################################### ### chunk number 41: estimate methylation status of a LumiMethyM object ################################################### #line 534 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## estimate the methylation status of a LumiMethyM object methyCall <- lumiMethyStatus(example.lumiMethy[,1:4]) head(methyCall) ## retrieve the methylation probability matrix methyProb <- attr(methyCall, "probability") head(methyProb) ################################################### ### chunk number 42: user defined preprocessing functions eval=FALSE ################################################### ## #line 548 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## ## suppose "userB" is a user defined background adjustment method ## lumiMethy.b.adj <- lumiMethyB(example.lumiMethy, method=userB, separateColor=TRUE) # not Run ## ## ## suppose "userC" is a user defined color balance adjustment method ## lumiMethy.c.adj <- lumiMethyC(example.lumiMethy, method=userC, separateColor=TRUE) # not Run ## ## ## suppose "userN" is a user defined probe level normalization method ## lumiMethy.c.n <- lumiMethyN(lumiMethy.c.adj, method= userN, separateColor=TRUE) # not Run ################################################### ### chunk number 43: Color channel information ################################################### #line 563 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## retrieve the featureData information ff <- pData(featureData(example.lumiMethy)) ## show the color channel information head(ff) ## add user provided color channel information if it is not existed in the featureData # example.lumiMethy <- addColorChannelInfo(example.lumiMethy, lib="IlluminaHumanMethylation27k.db") ################################################### ### chunk number 44: options to separately process each color channel eval=FALSE ################################################### ## #line 577 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## ## suppose "userB" is a user defined background adjustment method ## lumiMethy.b.adj <- lumiMethyB(example.lumiMethy, separateColor=TRUE) ## ## ## suppose "userC" is a user defined color balance adjustment method ## lumiMethy.c.adj <- lumiMethyC(example.lumiMethy, separateColor=TRUE) ## ## ## suppose "userN" is a user defined probe level normalization method ## lumiMethy.c.n <- lumiMethyN(lumiMethy.c.adj, separateColor=TRUE) ################################################### ### chunk number 45: Estimate detection call of a CpG site eval=FALSE ################################################### ## #line 592 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## ## Estimate the detection call of a CpG site ## presentCount <- detectionCall(example.lumiMethy) ################################################### ### chunk number 46: Preprocessing methylation microarray eval=FALSE ################################################### ## #line 611 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" ## ## library(lumi) ## ## specify the file name ## # fileName <- 'Example_Illumina_Methylation_profile.txt' ## ## load the data ## # example.lumiMethy <- lumiMethyR(fileName, lib="IlluminaHumanMethylation27k.db") # Not Run ## ## ## Quality and color balance assessment ## data(example.lumiMethy) ## ## summary of the example data ## example.lumiMethy ## ## ## preprocessing and quality control after normalization ## plotColorBias1D(example.lumiMethy, channel='sum') ## boxplotColorBias(example.lumiMethy, channel='sum') ## ## select interested sample to further check color balance in 2D scatter plot ## plotColorBias2D(example.lumiMethy, selSample=1) ## ## ## Color balance adjustment between two color channels ## lumiMethy.c.adj <- lumiMethyC(example.lumiMethy) ## ## Check color balance after color balance adjustment ## boxplotColorBias(lumiMethy.c.adj, channel='sum') ## ## ## Background adjustment is skipped because the SSN normalization includes background adjustment ## ## ## Normalization ## ## Perform SSN normalization based on color balance adjusted data ## lumiMethy.c.ssn <- lumiMethyN(lumiMethy.c.adj, method='ssn') ## ## ## Or we can perform quantile normalization based on color balance adjusted data ## # lumiMethy.c.q <- lumiMethyN(lumiMethy.c.adj, method='quantile') ## ## ## plot the density of M-values after SSN normalization ## density(lumiMethy.c.ssn, main="Density plot of M-value after SSN normalization") ## ## comparing with the density of M-values before normalization ## density(example.lumiMethy, main="Density plot of M-value of the raw data") ## ## ## output the normlized M-value as a Tab-separated text file ## write.exprs(lumiMethy.c.ssn, file='processedMethylationExampleData.txt') ################################################### ### chunk number 47: sessionInfo ################################################### #line 671 "vignettes/lumi/inst/doc/methylationAnalysis.Rnw" toLatex(sessionInfo())