R實(shí)戰(zhàn) | NGS數(shù)據(jù)時(shí)間序列分析(maSigPro)

R實(shí)戰(zhàn) | NGS數(shù)據(jù)時(shí)間序列分析(maSigPro)

[圖片上傳失敗...(image-5355e7-1653975368931)]

[圖片上傳失敗...(image-533c51-1653975368931)]

跟著Cell學(xué)作圖 | 6.時(shí)間序列分析(Mfuzz包)

一個(gè)答疑教程。

[圖片上傳失敗...(image-95c0e7-1653975368931)]

22

示例數(shù)據(jù)

#BiocManager::install('maSigPro')
library(maSigPro)
# 載入示例數(shù)據(jù)
data(data.abiotic) 
data.abiotic[1:5,1:5]
data(edesign.abiotic)
head(edesign.abiotic)

> data.abiotic[1:5,1:5]
        Control_3H_1 Control_3H_2 Control_3H_3 Control_9H_1 Control_9H_2
STMDF90   0.13735714   -0.3653065  -0.15329448   0.44754535  0.287476796
STMCJ24           NA           NA           NA           NA           NA
STMJH42   0.07864449    0.1002328  -0.17365488  -0.25279484  0.184855409
STMDE66   0.22976991    0.4740975   0.46930716   0.37101059 -0.004992029
STMIX74   0.14407618   -0.4801864  -0.07847999   0.05692331  0.013045420

> head(edesign.abiotic)
             Time Replicate Control Cold Heat Salt
Control_3H_1    3         1       1    0    0    0
Control_3H_2    3         1       1    0    0    0
Control_3H_3    3         1       1    0    0    0
Control_9H_1    9         2       1    0    0    0
Control_9H_2    9         2       1    0    0    0
Control_9H_3    9         2       1    0    0    0

注意：data.abiotic是已經(jīng)標(biāo)準(zhǔn)化過(guò)的基因表達(dá)矩陣。

建立回歸模型

生成回歸矩陣(makeDesignMatrix)

design <- make.design.matrix(edesign.abiotic, degree = 2)
design$groups.vector

示例數(shù)據(jù)有三個(gè)時(shí)間的，故考慮二次回歸模型(degree = 2)。

> design$groups.vector
 [1] "ColdvsControl" "HeatvsControl" "SaltvsControl" "Control"      
 [5] "ColdvsControl" "HeatvsControl" "SaltvsControl" "Control"      
 [9] "ColdvsControl" "HeatvsControl" "SaltvsControl"

尋找重要基因(p.vector)

F檢驗(yàn)確定回歸方程的顯著性，采用BH的校正方式，校正多重假設(shè)檢驗(yàn)的p值。

校正后的p值小于p.vector的參數(shù)Q的基因就作為候選基因，進(jìn)行下一步的分析。通過(guò)fit$SELEC可以得到候選基因的表達(dá)量信息。

fit <- p.vector(data.abiotic, # 標(biāo)準(zhǔn)化的表達(dá)矩陣 
                design, # 實(shí)驗(yàn)設(shè)計(jì)的矩陣 make.design.matrix 生成
                Q = 0.05, # 顯著性水平
                MT.adjust = "BH", 
                min.obs = 20 
                # 最低表達(dá)樣本數(shù) 不應(yīng)小于(degree+1)xGroups+1 
                )

fit$i # 顯著性基因的數(shù)量 
fit$SELEC # 顯著性基因表達(dá)矩陣

尋找顯著性差異(T.fit)

上述的回歸方程是基于所有的自變量的組合構(gòu)建的，接下來(lái)就是通過(guò)逐步回歸法確定最佳的自變量組合。

tstep <- T.fit(fit, # p.vector結(jié)果
               step.method = "backward", 
               alfa = 0.05) # 在逐步回歸中用于變量選擇的顯著性水平

在挑選最佳的自變量組合時(shí)，通過(guò)每種自變量組合對(duì)應(yīng)的回歸模型的擬合優(yōu)度值R-squared來(lái)進(jìn)行判斷，R-squared取值范圍為0到1，數(shù)值越大，越接近1，回歸模型的效果越好。

獲取顯著性基因列表(get.siggenes)

sigs <- get.siggenes(tstep, # T.fit結(jié)果
                     rsq = 0.6, # 逐步回歸中的R-squared截至值
                     vars = "groups")
# vars參數(shù)有3種
# all  每個(gè)基因直接給出一個(gè)最佳的回歸模型
# groups  只給出不同實(shí)驗(yàn)條件下相比control組中的差異基因
# each 會(huì)給出時(shí)間點(diǎn)和實(shí)驗(yàn)條件的所有組合對(duì)應(yīng)差異基因列表

names(sigs)
names(sigs$sig.genes$ColdvsControl)
sigs$sig.genes$ColdvsControl$sig.profiles # 查看cold vs control的結(jié)果

結(jié)果可視化

韋恩圖(suma2Venn)

suma2Venn(sigs$summary[, c(2:4)]) # 左圖
suma2Venn(sigs$summary[, c(1:4)]) # 右圖
# 這個(gè)韋恩圖面積大小與數(shù)量不成比例 較普通

[圖片上傳失敗...(image-53f699-1653975368931)]

see.genes()

see.genes(sigs$sig.genes$ColdvsControl, # 差異基因表達(dá)矩陣
          show.fit = T, # 是否顯示回歸擬合線(虛線)
          dis =design$dis, # 回歸設(shè)計(jì)矩陣
          cluster.method="hclust" , # 聚類方法
          cluster.data = 1, 
          k = 9) # 聚類數(shù)目

[圖片上傳失敗...(image-e14690-1653975368931)]

這一步生成兩個(gè)圖，如圖可分別查看。注意調(diào)整圖片顯示區(qū)域大小，以免報(bào)錯(cuò)。

[圖片上傳失敗...(image-d192c2-1653975368931)]

[圖片上傳失敗...(image-6a30f8-1653975368931)]

PlotGroups()

選擇某一特定genes的表達(dá)進(jìn)行可視化。

# 選取STMDE66基因
STMDE66 <- data.abiotic[rownames(data.abiotic)=="STMDE66", ]
PlotGroups (STMDE66, 
            edesign = edesign.abiotic)

[圖片上傳失敗...(image-ba61b6-1653975368931)]

# 添加回歸擬合線
PlotGroups (STMDE66, 
            edesign = edesign.abiotic, 
            show.fit = T, 
            dis = design$dis, 
            groups.vector = design$groups.vector)

[圖片上傳失敗...(image-cb99f-1653975368931)]

示例數(shù)據(jù)和代碼領(lǐng)取

詳見(jiàn)：https://mp.weixin.qq.com/s/MFdPrXrD_gErao-Jtsu6Fg

參考

Bioconductor - maSigPro(https://bioconductor.org/packages/release/bioc/html/maSigPro.html)

往期內(nèi)容

  1. (免費(fèi)教程+代碼領(lǐng)取)|跟著Cell學(xué)作圖系列合集
  2. Q&A | 如何在論文中畫出漂亮的插圖？
  3. Front Immunol 復(fù)現(xiàn) | 4. 使用estimate包評(píng)估腫瘤純度
  4. R繪圖 | 氣泡散點(diǎn)圖+擬合曲線
  5. 跟著 Cell 學(xué)作圖 | 桑葚圖(ggalluvial)
  6. R繪圖 | 對(duì)比條形圖+連線
  7. R繪圖 | 一幅小提琴圖的美化之旅
  8. R實(shí)戰(zhàn) | 給聚類加個(gè)圈圈(ggunchull)
  9. R繪圖 | 描述性統(tǒng)計(jì)常用圖(散點(diǎn)圖+柱狀圖+餅圖)

[圖片上傳失敗...(image-d6ce5c-1653975368931)]# R實(shí)戰(zhàn) | NGS數(shù)據(jù)時(shí)間序列分析(maSigPro)

[圖片上傳失敗...(image-50f56f-1653975368931)]

[圖片上傳失敗...(image-3758c5-1653975368931)]

跟著Cell學(xué)作圖 | 6.時(shí)間序列分析(Mfuzz包)

一個(gè)答疑教程。

[圖片上傳失敗...(image-1297de-1653975368931)]

22

示例數(shù)據(jù)

#BiocManager::install('maSigPro')
library(maSigPro)
# 載入示例數(shù)據(jù)
data(data.abiotic) 
data.abiotic[1:5,1:5]
data(edesign.abiotic)
head(edesign.abiotic)

> data.abiotic[1:5,1:5]
        Control_3H_1 Control_3H_2 Control_3H_3 Control_9H_1 Control_9H_2
STMDF90   0.13735714   -0.3653065  -0.15329448   0.44754535  0.287476796
STMCJ24           NA           NA           NA           NA           NA
STMJH42   0.07864449    0.1002328  -0.17365488  -0.25279484  0.184855409
STMDE66   0.22976991    0.4740975   0.46930716   0.37101059 -0.004992029
STMIX74   0.14407618   -0.4801864  -0.07847999   0.05692331  0.013045420

> head(edesign.abiotic)
             Time Replicate Control Cold Heat Salt
Control_3H_1    3         1       1    0    0    0
Control_3H_2    3         1       1    0    0    0
Control_3H_3    3         1       1    0    0    0
Control_9H_1    9         2       1    0    0    0
Control_9H_2    9         2       1    0    0    0
Control_9H_3    9         2       1    0    0    0

注意：data.abiotic是已經(jīng)標(biāo)準(zhǔn)化過(guò)的基因表達(dá)矩陣。

建立回歸模型

生成回歸矩陣(makeDesignMatrix)

design <- make.design.matrix(edesign.abiotic, degree = 2)
design$groups.vector

示例數(shù)據(jù)有三個(gè)時(shí)間的，故考慮二次回歸模型(degree = 2)。

> design$groups.vector
 [1] "ColdvsControl" "HeatvsControl" "SaltvsControl" "Control"      
 [5] "ColdvsControl" "HeatvsControl" "SaltvsControl" "Control"      
 [9] "ColdvsControl" "HeatvsControl" "SaltvsControl"

尋找重要基因(p.vector)

F檢驗(yàn)確定回歸方程的顯著性，采用BH的校正方式，校正多重假設(shè)檢驗(yàn)的p值。

校正后的p值小于p.vector的參數(shù)Q的基因就作為候選基因，進(jìn)行下一步的分析。通過(guò)fit$SELEC可以得到候選基因的表達(dá)量信息。

fit <- p.vector(data.abiotic, # 標(biāo)準(zhǔn)化的表達(dá)矩陣 
                design, # 實(shí)驗(yàn)設(shè)計(jì)的矩陣 make.design.matrix 生成
                Q = 0.05, # 顯著性水平
                MT.adjust = "BH", 
                min.obs = 20 
                # 最低表達(dá)樣本數(shù) 不應(yīng)小于(degree+1)xGroups+1 
                )

fit$i # 顯著性基因的數(shù)量 
fit$SELEC # 顯著性基因表達(dá)矩陣

尋找顯著性差異(T.fit)

上述的回歸方程是基于所有的自變量的組合構(gòu)建的，接下來(lái)就是通過(guò)逐步回歸法確定最佳的自變量組合。

tstep <- T.fit(fit, # p.vector結(jié)果
               step.method = "backward", 
               alfa = 0.05) # 在逐步回歸中用于變量選擇的顯著性水平

在挑選最佳的自變量組合時(shí)，通過(guò)每種自變量組合對(duì)應(yīng)的回歸模型的擬合優(yōu)度值R-squared來(lái)進(jìn)行判斷，R-squared取值范圍為0到1，數(shù)值越大，越接近1，回歸模型的效果越好。

獲取顯著性基因列表(get.siggenes)

sigs <- get.siggenes(tstep, # T.fit結(jié)果
                     rsq = 0.6, # 逐步回歸中的R-squared截至值
                     vars = "groups")
# vars參數(shù)有3種
# all  每個(gè)基因直接給出一個(gè)最佳的回歸模型
# groups  只給出不同實(shí)驗(yàn)條件下相比control組中的差異基因
# each 會(huì)給出時(shí)間點(diǎn)和實(shí)驗(yàn)條件的所有組合對(duì)應(yīng)差異基因列表

names(sigs)
names(sigs$sig.genes$ColdvsControl)
sigs$sig.genes$ColdvsControl$sig.profiles # 查看cold vs control的結(jié)果

結(jié)果可視化

韋恩圖(suma2Venn)

suma2Venn(sigs$summary[, c(2:4)]) # 左圖
suma2Venn(sigs$summary[, c(1:4)]) # 右圖
# 這個(gè)韋恩圖面積大小與數(shù)量不成比例 較普通

[圖片上傳失敗...(image-70965a-1653975368931)]

see.genes()

see.genes(sigs$sig.genes$ColdvsControl, # 差異基因表達(dá)矩陣
          show.fit = T, # 是否顯示回歸擬合線(虛線)
          dis =design$dis, # 回歸設(shè)計(jì)矩陣
          cluster.method="hclust" , # 聚類方法
          cluster.data = 1, 
          k = 9) # 聚類數(shù)目

[圖片上傳失敗...(image-a2527d-1653975368931)]

這一步生成兩個(gè)圖，如圖可分別查看。注意調(diào)整圖片顯示區(qū)域大小，以免報(bào)錯(cuò)。

[圖片上傳失敗...(image-a09102-1653975368931)]

[圖片上傳失敗...(image-a4c051-1653975368931)]

PlotGroups()

選擇某一特定genes的表達(dá)進(jìn)行可視化。

# 選取STMDE66基因
STMDE66 <- data.abiotic[rownames(data.abiotic)=="STMDE66", ]
PlotGroups (STMDE66, 
            edesign = edesign.abiotic)

[圖片上傳失敗...(image-1045d2-1653975368931)]

# 添加回歸擬合線
PlotGroups (STMDE66, 
            edesign = edesign.abiotic, 
            show.fit = T, 
            dis = design$dis, 
            groups.vector = design$groups.vector)

[圖片上傳失敗...(image-2d7f07-1653975368931)]

示例數(shù)據(jù)和代碼領(lǐng)取

點(diǎn)贊、在看 本文，分享至朋友圈集贊20個(gè)并保留30分鐘，截圖發(fā)至微信mzbj0002領(lǐng)取。

木舟筆記2022年度VIP可免費(fèi)領(lǐng)取。

木舟筆記2022年度VIP企劃

權(quán)益：

  1. 2022年度木舟筆記所有推文示例數(shù)據(jù)及代碼(在VIP群里實(shí)時(shí)更新)。

    <img src="https://mzbj-1312199408.cos.ap-shanghai.myqcloud.com/image-20220529221704900.png" alt="image-20220529221704900" style="zoom:67%;" />

  2. 木舟筆記科研交流群。

  3. 半價(jià)購(gòu)買跟著Cell學(xué)作圖系列合集(免費(fèi)教程+代碼領(lǐng)取)|跟著Cell學(xué)作圖系列合集。

收費(fèi)：

**99￥/人**?？商砑游⑿牛篳mzbj0002` 轉(zhuǎn)賬，或直接在文末打賞。

<img src="https://mzbj-1312199408.cos.ap-shanghai.myqcloud.com/image-20220529221602758.png" alt="image-20220529221602758" style="zoom:50%;" />

參考

Bioconductor - maSigPro(https://bioconductor.org/packages/release/bioc/html/maSigPro.html)

往期內(nèi)容

  1. (免費(fèi)教程+代碼領(lǐng)取)|跟著Cell學(xué)作圖系列合集
  2. Q&A | 如何在論文中畫出漂亮的插圖？
  3. Front Immunol 復(fù)現(xiàn) | 4. 使用estimate包評(píng)估腫瘤純度
  4. R繪圖 | 氣泡散點(diǎn)圖+擬合曲線
  5. 跟著 Cell 學(xué)作圖 | 桑葚圖(ggalluvial)
  6. R繪圖 | 對(duì)比條形圖+連線
  7. R繪圖 | 一幅小提琴圖的美化之旅
  8. R實(shí)戰(zhàn) | 給聚類加個(gè)圈圈(ggunchull)
  9. R繪圖 | 描述性統(tǒng)計(jì)常用圖(散點(diǎn)圖+柱狀圖+餅圖)

[圖片上傳失敗...(image-fc3c55-1653975368931)]

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請(qǐng)聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成，瀏覽時(shí)請(qǐng)結(jié)合常識(shí)與多方信息審慎甄別。
平臺(tái)聲明：文章內(nèi)容（如有圖片或視頻亦包括在內(nèi)）由作者上傳并發(fā)布，文章內(nèi)容僅代表作者本人觀點(diǎn)，簡(jiǎn)書系信息發(fā)布平臺(tái)，僅提供信息存儲(chǔ)服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容