經(jīng)過預處理之后的數(shù)據(jù),就可以進行差異分析了。對于甲基化芯片而言,有兩個方面的差異分析 DMP 差異甲基化探針 DMR 差異甲基化區(qū)域
在ChAMP 包中,champ.DMP 函數(shù)用于分析差異甲基化探針,champ.DMR 函數(shù)用于分析差異甲基化區(qū)域。本章我們先看下差異探針的分析 champ.DMP 函數(shù)的用法示例
myDMP <- champ.DMP()
在差異分析時,我們需要兩個輸入數(shù)據(jù),一個就是探針的表達譜數(shù)據(jù),beta matrix, 另外一個就是樣本的分組信息。 在champ.DMP 函數(shù)中,默認myNorm 作為歸一化之后的beta matrix,對于樣本的分組信息,ChAMP 默認從Samplesheet.csv文件中讀取,在數(shù)據(jù)導入成功后,myLoad$pd 代表的就是SampleSheet.csv文件的信息,所以myLoad$pd$Sample_Group 代表樣本的分組信息。 在差異分析時,最關(guān)鍵的就是差異分組問題。在實驗設計階段,有很多類型的分組設計,比如最常見的case_vs_control, 兩個group的分組;多個組織,比如3個組織,共3個group; 時間序列,比如藥物處理后的幾個時間點。不同的實驗設計,在差異分析時,想要關(guān)注的差異點自然不同,在分析時也要采取不同的分析策略。 對于ChAMP 來說,上述的幾種分組設計都是支持。 champ.DMP 計算過程分為以下3步:
1. 檢測樣本分組,確定分組比較測試數(shù)據(jù)分成T和C兩組,每組各4個樣本、 在這一步,需要確定兩個因素: 分組的類型 主要識別分組變量是離散型還是連續(xù)型,如果是字符串型character , 就是離散型,如果是數(shù)值型numeric , 就是連續(xù)型。不同的類型,采取的差異分析策略不同。 測試數(shù)據(jù)是字符型的兩個group, 具體的輸入信息如下 分組的個數(shù) 確定group的個數(shù),2個group 肯定是兩者之間進行差異分析,但是當group 個數(shù)3個或以上時,就需要確定如何分組比較。 默認情況下兩兩之間都進行差異分析,如果你不需要這么多的差異結(jié)果,可以通過compare.group 參數(shù)指定, compare.group 參數(shù)的值是一個list, list 中的每個元素是一個長度為2的向量,指定了用于差異分析的兩個group
[ Section 1: Check Input Pheno Start ] You pheno is character type. Your pheno information contains following groups. >> C:4 samples. T:4 samples. [The power of statistics analysis on groups contain very few samples may not strong.] pheno contains only 2 phenotypes compare.group parameter is NULL, two pheno types will be added into Compare List. C_to_T compare group : C, T [ Section 1: Check Input Pheno Done ]
進行差異探針分析通過調(diào)用limma 函數(shù)進行差異分析,默認通過BH 方法進行多重建設檢驗的校正,p.adjust < 0.05 的認為是差異探針 可以通過adjPVal 參數(shù)修改p.adjust的閾值,當然也可以修改adjust.method 參數(shù)的值,調(diào)整多重假設檢驗校正的算法,默認值為BH , 可選值包括 “none”, “BH”, “BY”, “holm”。 [ Section 2: Find Differential Methylated CpGs Start ] Start to Compare : C, T Contrast Matrix Contrasts Levels pT-pC pC -1 pT 1 You have found 4283 significant MVPs with a BH adjusted P-value below 0.05. Calculate DMP for C and T done. [ Section 2: Find Numeric Vector Related CpGs Done ]
添加差異探針的注釋信息之前的分析都是針對探針的beta matrix 進行的分析,找的差異探針之后,我們肯定希望知道這個探針對應的基因,染色體位置等注釋信息。這一步實際就是在已有的差異結(jié)果的基礎上,追加探針的注釋信息。 [ Section 3: Match Annotation Start ] [ Section 3: Match Annotation Done ]
結(jié)果展示str(myDMP) List of 1 $ C_to_T:’data.frame’: 4283 obs. of 23 variables: ..$ logFC : num [1:4283] 0.724 … ..$ AveExpr : num [1:4283] 0.398 … ..$ t : num [1:4283] 28.2 … ..$ P.Value : num [1:4283] 1.74e-08… ..$ adj.P.Val : num [1:4283] 0.00703… ..$ B : num [1:4283] 7.6 7.05 … ..$ C_AVG : num [1:4283] 0.0358 … ..$ T_AVG : num [1:4283] 0.759 … ..$ deltaBeta : num [1:4283] 0.724 … ..$ CHR : Factor w/ 25 levels “”,”1”,”10”,”11”,..: ..$ MAPINFO : int [1:4283] 141516291 … ..$ Strand : Factor w/ 3 levels “”,”F”,”R”: 3 2 … ..$ Type : Factor w/ 2 levels “I”,”II”: 2 2 … ..$ gene : Factor w/ 20622 levels “”,”A1BG” … ..$ feature : chr [1:4283] “Body” “Body” … ..$ cgi : Factor w/ 4 levels “island”,”opensea”,..: … ..$ feat.cgi : Factor w/ 28 levels “1stExon-island”,..: 13 … ..$ UCSC_CpG_Islands_Name: Factor w/ 27177 levels “”,”chr1:10003165-10003585”,..: 18278 … ..$ DHS : logi [1:4283] NA NA NA NA NA NA … ..$ Enhancer : logi [1:4283] TRUE TRUE NA NA NA NA … ..$ Phantom : Factor w/ 11912 levels “”,”high-CpG:100009860- … ..$ Probe_SNPs : Factor w/ 56004 levels “”,”rs10000615”,..: ..$ Probe_SNPs_10 : Factor w/ 35790 levels “”,”rs10000804” …
myDMP 就是最終的差異分析結(jié)果,是一個list 對象,list中的每個元素是兩個group之間差異分析的結(jié)果。
測試數(shù)據(jù)只有兩個分組,所以list 中只有一個元素。差異分析的結(jié)果是一個data.frame 對象,可以分成3個部分。 從logFC 到B 的部分是limma 差異輸出結(jié)果, C_AVG 到deltaBeta 是每組表達量的均值,deltaBate 是兩組均值的差,CHR 到Probe_SNPs_10 是探針的注釋信息。
|