WGCNA懒人秘籍教程进阶版
大家好,我是风风。今天我们继续我们的WGCNA系列。在上次的推文中,我们对单个数据集进行了WGCNA分析,完整地走完了单个WGCNA的流程,包括:数据清洗、构建网络并进行模块识别、关联表型信息和模块、网络可视化四个方面内容。但是,如果我们现在手上有两份数据集或者多份数据集,我们可以怎么运用这些数据进行WGCNA分析呢?
这就是我们本期推文要来一起探讨的内容——多数据集共识分析。
数据清洗
这里我们适用WGCNA文档提供的两份数据,第一份数据与上期数据一致,是雌性小鼠的肝脏数据,此外,我们再下载一份雄鼠肝脏数据数据在后台回复关键词获取。这两组数据在生物学上非常相似(都是小鼠肝脏数据,只有性别差异),但也存在显着差异。既然有多份数据,那自然就需要用到循环方面的知识。在讲席营的课程中,我给讲席营学员总结了写循环的步骤,这里大家可以反其道而行,看看本期WGCNA的循环,然后 还原出每一个步骤,有助于大家理解代码。接下来我们一样先对数据进行清洗:
rm(list = ls())# BiocManager::install("WGCNA") 安装R包library(WGCNA)options(stringsAsFactors = FALSE)# 加载LiverFemale3600数据集femData <- read.csv("data//LiverFemale3600.csv")# 加载LiverMale3600数据集maleData <- read.csv("data//LiverMale3600.csv");# 查看数据基本信息dim(femData)## [1] 3600 143head(femData)## substanceBXH gene_symbol LocusLinkID ProteomeID cytogeneticLoc CHROMOSOME## 1 MMT00000044 1700007N18Rik 69339 286025 0 16## 2 MMT00000046 Mast2 17776 157466 0 4## 3 MMT00000051 Ankrd32 105377 321939 0 13## 4 MMT00000076 0 383154 0 0 16## 5 MMT00000080 Ldb2 16826 157383 0 5## 6 MMT00000102 Rdhs 216453 0 10_70.0_cM 10## StartPosition EndPosition F2_2 F2_3 F2_14 F2_15 F2_19## 1 50911260 50912491 -0.01810 0.0642 6.44e-05 -0.05800 0.04830## 2 115215318 115372404 -0.07730 -0.0297 1.12e-01 -0.05890 0.04430## 3 74940309 74982847 -0.02260 0.0617 -1.29e-01 0.08710 -0.11500## 4 49345114 49477048 -0.00924 -0.1450 2.87e-02 -0.04390 0.00425## 5 43546124 43613704 -0.04870 0.0582 -4.83e-02 -0.03710 0.02510## 6 1337265 1347607 0.17600 -0.1890 -6.50e-02 -0.00846 -0.00574## F2_20 F2_23 F2_24 F2_26 F2_37 F2_42 F2_43 F2_45## 1 -0.15197410 -0.00129 -0.23600 -0.0307 -0.02610 0.073705890 -0.0466 -0.00673## 2 -0.09380000 0.09340 0.02690 -0.1330 0.07570 -0.009193803 -0.0075 0.01700## 3 -0.06502607 0.00249 -0.10200 0.1420 -0.10200 0.064289290 0.0169 -0.01590## 4 -0.23610000 -0.06900 0.01440 0.0363 -0.01820 0.477874600 0.1440 0.11100## 5 0.08504274 0.04450 0.00167 -0.0680 0.00567 -0.075348680 -0.0673 -0.04720## 6 -0.01807182 -0.12500 -0.06820 0.1250 0.00998 -0.037366600 -0.0402 -0.02190## F2_46 F2_47 F2_48 F2_51 F2_52 F2_54 F2_63 F2_65 F2_66## 1 -0.0193 0.09040 0.0290 0.0356 -0.0388 -0.0360 -0.05600 0.009840 -0.0261## 2 0.0722 -0.08390 0.0273 -0.0784 -0.0178 0.1120 0.12300 0.051700 0.0731## 3 -0.1430 -0.00492 -0.0735 0.0657 -0.0197 -0.1290 -0.14300 -0.061600 0.0419## 4 0.0113 0.11900 0.0225 0.0932 0.1430 0.2640 -0.09280 -0.000635 -0.0126## 5 0.0701 -0.08790 -0.0180 -0.1290 -0.0469 -0.0352 -0.00166 0.058700 -0.0206## 6 0.0269 0.13300 0.0732 0.1070 -0.0362 -0.0696 -0.19400 -0.117000 -0.0400## F2_68 F2_69 F2_70 F2_71 F2_72 F2_78 F2_79 F2_80## 1 0.00856 -0.01180 -0.03350 -0.08310 -0.0471 -0.02820 0.047264410 0.0296## 2 0.08670 0.05710 0.00693 -0.00606 -0.0390 0.01870 0.008471275 -0.0687## 3 -0.29000 -0.10800 -0.09950 -0.00315 0.0975 0.01030 -0.134271000 0.1010## 4 0.06910 0.02260 -0.08630 -0.22900 0.0178 0.00166 0.064096960 0.0103## 5 -0.13000 0.00392 0.05450 -0.11200 0.1070 0.01190 0.008985630 -0.1030## 6 0.06890 0.04320 -0.00338 -0.05270 -0.0416 -0.03040 0.025920240 0.0697## F2_81 F2_83 F2_86 F2_87 F2_88 F2_89 F2_107 F2_108 F2_109## 1 0.0114 0.0498 -0.0249 -0.00264 -0.02050 0.0826 -0.0421 0.0663 0.03620## 2 -0.0114 -0.0262 -0.0215 -0.09580 -0.01930 -0.1140 0.0815 0.0285 0.00299## 3 0.0521 -0.0607 -0.0285 0.02560 -0.01350 0.0796 0.0553 -0.0380 0.12900## 4 -0.0258 -0.0837 0.1880 0.03310 -0.00652 0.1550 0.0458 0.0752 0.12200## 5 -0.1400 -0.0282 -0.1090 0.02070 -0.01370 -0.0288 -0.1220 0.1270 -0.09390## 6 0.1150 0.0953 0.0127 0.05490 0.00311 0.0955 -0.1520 -0.0670 -0.00599## F2_110 F2_111 F2_112 F2_117 F2_119 F2_125 F2_126 F2_127 F2_141## 1 0.0808 -0.0404 0.0877 0.07240 -0.0210 0.04540 -0.03220 -0.00654 0.03490## 2 -0.0407 -0.0657 0.0643 -0.00022 -0.0877 0.00167 0.00321 -0.01260 -0.04530## 3 -0.0361 0.0441 -0.1640 -0.01420 -0.0279 0.00677 0.07360 0.01750 0.10900## 4 -0.0104 0.0914 -0.0355 0.06520 0.1280 0.05940 0.01630 0.00292 0.00714## 5 0.1200 -0.0850 0.1400 0.00867 0.1440 0.08710 -0.03360 0.17300 0.08270## 6 -0.0438 0.0634 0.1380 -0.04010 0.1310 -0.12600 0.00484 -0.00256 -0.06800## F2_142 F2_143 F2_144 F2_145 F2_154 F2_155 F2_156 F2_157## 1 -0.0315 -0.02170 0.00370 0.0322 -0.02150730 -0.000958 -0.0850 0.00462## 2 -0.0579 0.05920 0.00239 -0.0383 0.02457782 -0.030300 -0.1260 -0.06670## 3 -0.0216 -0.01250 0.05460 0.0403 -0.01674888 0.059900 0.0311 -0.05190## 4 -0.0565 0.10200 0.03480 0.0245 0.06776892 0.016500 -0.0382 0.02120## 5 0.0594 -0.00317 -0.06750 0.0495 0.13520570 0.016500 0.0832 0.04350## 6 0.0941 -0.04220 0.12000 0.1080 -0.05128296 -0.005590 0.0136 0.09910## F2_162 F2_163 F2_164 F2_165 F2_166 F2_167 F2_169 F2_180## 1 0.03990 0.0716 -0.0923 0.10900 0.0102 0.0337 0.00911 0.03210## 2 -0.00637 -0.0161 -0.2340 -0.09610 -0.1290 -0.0109 -0.11300 -0.00677## 3 0.01890 0.0207 0.0929 0.00917 0.0874 -0.1260 -0.00949 -0.09900## 4 0.06690 0.0512 -0.2450 1.23000 -0.0402 -0.0635 0.06880 0.03790## 5 0.19300 0.0586 -0.0768 0.04600 0.0484 0.2810 0.07210 -0.00630## 6 0.06770 -0.0520 0.1550 0.07890 0.0336 0.0648 0.14400 0.02770## F2_181 F2_182 F2_187 F2_188 F2_189 F2_190 F2_191 F2_192## 1 0.03144772 0.0543 0.01120 0.01060 0.1130 -0.03960 -0.0504 0.0877## 2 -0.16704700 -0.0239 0.00304 -0.03580 -0.1330 -0.01830 -0.0623 -0.0648## 3 0.02700180 -0.0570 -0.05160 -0.04970 0.1660 0.05000 0.0498 0.0431## 4 -0.02058180 0.0227 0.04180 0.01010 0.2170 0.00206 -0.0155 0.6550## 5 0.37074790 0.0618 0.10800 0.12100 0.0237 0.02960 0.1130 0.0839## 6 0.09297908 0.0601 0.02960 0.00198 0.0251 0.00059 -0.0282 0.0429## F2_194 F2_195 F2_200 F2_201 F2_212 F2_213 F2_214 F2_215 F2_221## 1 -0.0563 -0.00557 -0.0484 -0.0273 -0.10816380 -0.0183 -0.0132 -0.00432 -0.6630## 2 -0.0652 0.05020 -0.0912 -0.0180 0.05682362 -0.0238 0.0721 0.03910 0.1070## 3 -0.0224 -0.10700 0.0715 0.0432 -0.13217820 0.0205 -0.0411 0.07670 -0.0783## 4 0.2820 -0.01310 -0.0387 -0.0667 -0.32395020 -0.0245 0.0865 0.06470 -2.0000## 5 0.1050 0.15500 0.0823 0.1140 0.03542023 -0.2020 0.0822 0.04260 0.1030## 6 0.0697 0.04930 0.0414 -0.0708 -0.10881230 0.0359 -0.0678 -0.11000 -0.1420## F2_222 F2_223 F2_224 F2_225 F2_226 F2_227 F2_228 F2_241 F2_242## 1 0.01440 0.0310 0.00818 -0.00892 -0.08710 0.0129 0.0937 0.0313 0.0821## 2 0.00923 -0.0397 -0.06400 0.06300 -0.00152 0.0555 0.0947 -0.0387 0.0592## 3 -0.06860 -0.0254 -0.05680 -0.13300 -0.07560 -0.0557 -0.0890 -0.1460 -0.0739## 4 0.00874 0.0847 -0.09720 0.00746 -0.55200 0.0415 0.0733 0.0815 0.1100## 5 -0.10100 0.1630 0.07410 -0.01640 0.08700 -0.0557 -0.1910 0.0219 0.0913## 6 0.08430 -0.0610 0.08760 -0.03960 0.10200 0.0190 -0.1190 0.0687 -0.0525## F2_243 F2_244 F2_245 F2_247 F2_248 F2_261 F2_263 F2_264## 1 0.00621 0.0307 -0.13700 0.075300 -0.096881950 -0.01670 -0.0928 -0.00957## 2 -0.00636 0.0614 0.02850 -0.000633 0.001598228 -0.00267 -0.0198 0.16300## 3 -0.01120 -0.0528 0.05050 0.027700 -0.067933370 -0.02220 -0.0684 -0.04930## 4 0.21400 0.0135 -0.13500 -0.003100 0.072318780 0.01030 -0.3150 0.08420## 5 0.01120 0.1190 0.00383 0.041700 -0.038618510 0.11800 0.0123 0.03700## 6 -0.00716 -0.1460 -0.14500 0.029400 0.035281240 -0.05660 0.0917 -0.08080## F2_270 F2_271 F2_272 F2_278 F2_287 F2_288 F2_289 F2_290## 1 0.0287 -0.01300 -0.0292 -0.03810 -0.0488 0.17361240 -0.097900 0.0383## 2 -0.1310 -0.04260 -0.0514 0.07260 -0.0481 -0.16211430 -0.123000 -0.1370## 3 0.0328 0.00537 -0.0259 -0.14400 0.0170 0.25924220 -0.041400 -0.0229## 4 0.0351 NA 0.0730 0.00914 0.0556 0.18311140 0.051700 0.1780## 5 -0.0142 0.00563 -0.0504 -0.05970 -0.0871 0.20897910 -0.000188 -0.0328## 6 0.0362 0.00790 -0.0246 -0.07330 0.0125 -0.04778892 0.082500 0.1360## F2_291 F2_296 F2_298 F2_299 F2_300 F2_302 F2_303 F2_304## 1 0.01850 -0.08937784 0.0230 -0.06250 -0.000142 0.0344 0.0358 -0.0139## 2 -0.05720 -0.07416870 -0.0688 -0.06540 -0.102000 -0.0780 -0.0820 -0.1830## 3 -0.00664 -0.05915232 -0.0134 0.09740 0.015500 -0.0934 0.1780 0.0842## 4 0.05250 -0.21653720 -0.2210 -0.00266 0.545000 0.0127 0.0273 -0.0928## 5 -0.16600 -0.07897525 0.1410 -0.12900 0.090600 -0.1330 -0.2120 -0.0797## 6 0.04620 0.03811979 -0.0346 0.04690 -0.034800 0.0110 0.0323 0.1660## F2_305 F2_306 F2_307 F2_308 F2_309 F2_310 F2_311 F2_312## 1 0.0134 -0.03145069 0.02780 -0.01190 -0.0744 0.00197 -0.0151 -0.0721## 2 -0.0270 -0.09822316 -0.07890 -0.05480 -0.1320 -0.11000 -0.1130 -0.0805## 3 0.0870 0.15520470 0.03410 -0.06830 0.0555 -0.04060 0.0835 0.0514## 4 0.0469 0.10038160 -2.00000 0.05240 0.1260 0.07280 0.0600 -0.0455## 5 -0.0191 -0.11958500 0.00294 -0.10600 -0.0518 -0.13200 0.0494 0.0221## 6 -0.0866 0.05385017 0.09570 -0.00949 0.1120 0.20800 0.0872 -0.0555## F2_320 F2_321 F2_323 F2_324 F2_325 F2_326 F2_327 F2_328 F2_329## 1 -0.0118 0.0200 0.0222 0.047700 -0.0488 0.0168 -0.0309 0.02740 -0.0310## 2 -0.1200 0.0101 -0.1610 -0.049200 -0.0350 -0.0738 -0.1730 -0.07380 -0.2010## 3 0.0713 -0.1130 0.0466 0.000612 0.1210 0.0996 0.1090 0.02730 0.1200## 4 -0.0464 0.0667 -0.1850 -0.270000 0.0803 0.0424 0.1610 0.05120 0.2410## 5 0.0272 -0.0938 0.1020 0.113000 -0.0859 -0.1340 0.0639 0.00731 0.1240## 6 0.0748 -0.1420 0.0590 -0.080000 -0.1200 0.1230 0.1870 0.05410 0.0699## F2_330 F2_332 F2_355 F2_357## 1 0.0660 -0.0199 -0.0146 0.065000## 2 -0.0820 -0.0939 0.0192 -0.049900## 3 -0.0629 -0.0395 0.1090 0.000253## 4 0.3890 0.0251 -0.0348 0.114000## 5 -0.0212 0.0870 0.0512 0.024300## 6 0.0708 0.1450 -0.0399 0.037500dim(maleData)## [1] 3600 132head(maleData)## substanceBXH gene_symbol LocusLinkID ProteomeID cytogeneticLoc CHROMOSOME## 1 MMT00000044 1700007N18Rik 69339 286025 0 16## 2 MMT00000046 Mast2 17776 157466 0 4## 3 MMT00000051 Ankrd32 105377 321939 0 13## 4 MMT00000076 0 383154 0 0 16## 5 MMT00000080 Ldb2 16826 157383 0 5## 6 MMT00000102 Rdhs 216453 0 10_70.0_cM 10## StartPosition EndPosition F2_4 F2_5 F2_6 F2_7 F2_8 F2_9## 1 50911260 50912491 -0.0444 -0.0179 -0.0431 0.03580 0.0263 0.15400## 2 115215318 115372404 0.1250 0.0507 0.1290 0.13900 0.2370 -0.00483## 3 74940309 74982847 -0.1510 -0.0689 -0.0925 0.00353 -0.1610 -0.00932## 4 49345114 49477048 -0.1650 -0.0285 2.0000 0.04570 -0.4550 0.33200## 5 43546124 43613704 -0.0724 -0.0603 -0.0569 0.02610 -0.1130 -0.01210## 6 1337265 1347607 -0.1430 -0.0663 -0.1570 -0.23700 -0.2090 -0.09170## F2_10 F2_13 F2_16 F2_17 F2_18 F2_22 F2_27 F2_28 F2_29## 1 0.000109 0.0254 -0.0294 0.1160 0.0431 -0.0267 -0.2160 -0.12700 0.0377## 2 0.007490 0.0227 0.0355 0.0836 0.1230 0.1180 0.1200 0.16300 0.1570## 3 -0.191000 0.0809 0.0692 -0.1350 -0.0471 -0.0785 -0.0352 0.00584 -0.1070## 4 0.043500 0.0944 0.1640 0.0774 0.0169 -0.1030 -0.2080 -0.25600 0.0204## 5 -0.161000 0.0100 -0.1320 -0.1550 -0.1420 -0.0666 -0.0351 -0.03760 -0.0966## 6 0.060800 -0.1330 -0.0683 -0.2010 -0.2530 -0.2020 -0.1110 -0.12700 -0.0948## F2_30 F2_33 F2_34 F2_35 F2_39 F2_40 F2_41 F2_49 F2_50## 1 -0.07320 -0.0137 0.0434 -0.0277 0.0667 0.0283 0.0541 0.0533 -0.06555326## 2 0.20600 -0.0102 0.1460 0.1890 0.1170 0.2400 0.1560 0.0114 -0.02107601## 3 -0.07020 -0.0273 0.0426 0.0314 0.0751 -0.1070 -0.0586 -0.0698 -0.07634149## 4 -0.04560 -0.8740 -0.8230 0.2260 0.1750 0.0204 0.0801 -0.0481 -0.17293770## 5 0.00728 -0.0629 0.1210 -0.2050 0.0322 -0.0158 -0.0989 -0.0752 -0.03223757## 6 -0.19000 -0.1610 -0.1260 -0.1760 -0.1850 -0.2190 -0.2260 0.0867 -0.08595835## F2_55 F2_56 F2_57 F2_59 F2_60 F2_73 F2_74 F2_75 F2_76## 1 -0.00713 0.0453 0.0256 0.02944015 -0.0459 0.0338 -0.0458 0.0201 0.0300## 2 0.10900 0.1700 0.2540 0.08054645 0.1890 0.1640 0.0728 0.1230 0.1360## 3 -0.03310 -0.0901 -0.0965 -0.11589100 -0.0930 -0.0391 0.0406 -0.0223 -0.0397## 4 0.13600 0.0427 0.0187 0.35591750 0.0437 -0.2150 -0.0366 0.0152 0.0448## 5 -0.06150 0.0164 -0.1050 -0.05905863 -0.1030 0.0122 -0.1220 -0.0603 -0.0907## 6 -0.06300 -0.1770 -0.1320 -0.05455500 -0.2250 -0.1760 -0.0801 -0.1050 -0.1510## F2_84 F2_85 F2_91 F2_92 F2_93 F2_94 F2_104 F2_105 F2_114## 1 -0.0352 -0.1050 0.0259 0.0939 0.04060 0.05805066 -0.0118 0.0143 -0.08070## 2 0.2380 0.1000 0.2040 0.1950 0.06750 -0.09036969 0.2950 -0.0661 -0.02010## 3 -0.0299 -0.0903 -0.2060 -0.1140 -0.01200 -0.04731417 -0.1050 0.0588 0.00895## 4 0.4910 -0.5400 0.0573 -0.0314 0.08910 0.03246458 0.0498 0.0764 -0.07570## 5 -0.0313 -0.0243 -0.2260 0.0257 0.00118 -0.01082061 0.0462 0.0566 0.00530## 6 -0.1560 -0.1650 -0.0885 -0.2140 -0.08690 -0.01983479 -0.2880 -0.0425 -0.10000## F2_115 F2_116 F2_120 F2_121 F2_122 F2_123 F2_124 F2_146## 1 -0.0418 -0.0559 0.00961 0.02130 -0.000128 0.04350 0.01260 0.003750## 2 0.0179 0.0837 0.04040 0.15900 0.004370 0.02910 0.05050 0.049400## 3 0.1190 0.0474 -0.08880 -0.13600 0.052000 -0.00612 0.04040 0.008640## 4 0.0532 -0.1520 0.14000 -0.03820 -0.041300 0.09380 -0.11600 -0.048700## 5 0.0935 -0.0622 0.05640 0.00566 -0.000152 0.07480 -0.00657 -0.000285## 6 -0.1520 -0.1490 -0.03080 -0.10200 -0.093200 -0.04530 -0.16100 -0.085200## F2_147 F2_148 F2_149 F2_151 F2_152 F2_153 F2_158 F2_159 F2_160## 1 0.00994 -0.0225 0.0593 -0.00857 0.0288 0.0761 0.000479 -0.0189 0.0438## 2 0.17200 -0.0412 0.0968 0.04930 -0.0367 -0.1340 0.138000 -0.0126 0.0757## 3 0.02550 -0.0475 0.0802 0.04530 0.0184 0.0162 -0.052900 0.0576 -0.0076## 4 0.07400 0.0380 0.0568 -0.00238 -0.0396 0.0121 0.026400 0.0114 0.0108## 5 0.13500 0.1200 -0.0286 0.15700 -0.0247 0.1090 0.004630 -0.1240 -0.0387## 6 -0.18200 -0.0417 -0.1450 -0.04530 -0.0119 0.0662 -0.063400 0.0423 -0.0895## F2_170 F2_171 F2_172 F2_173 F2_174 F2_176 F2_178 F2_179 F2_183## 1 0.0149 0.02290 0.0812 -0.0100 0.0492 0.03220 0.07230 -0.0196 -0.05150## 2 0.0853 0.14800 -0.0538 0.1300 0.1850 0.02230 0.00528 0.0265 0.03850## 3 -0.0349 -0.03930 0.0696 0.0564 -0.0620 0.02440 0.00459 -0.0327 0.00872## 4 0.0861 0.01890 0.0772 0.0169 0.0694 0.00808 0.15500 -0.1810 -0.03080## 5 0.0269 0.00419 -0.0258 -0.1100 0.0790 0.08090 -0.02610 -0.0216 -0.08210## 6 -0.1090 -0.11600 0.0621 -0.1820 -0.1480 -0.09400 0.00701 -0.0180 0.06090## F2_184 F2_185 F2_186 F2_197 F2_198 F2_199 F2_207 F2_208## 1 0.00377 0.03590 0.02331811 0.08710 0.00320 -0.0152 0.0919 0.0745## 2 0.19300 0.06140 0.05443614 -0.09730 0.02270 0.0731 0.1870 0.1540## 3 -0.04460 -0.07370 -0.16528400 0.00276 0.00964 -0.0403 -0.0760 -0.0429## 4 -0.01700 -0.12100 -0.04767130 -0.06740 0.00838 0.0253 0.2100 -0.3510## 5 0.03000 0.00615 0.05199314 0.04700 0.04130 -0.0335 0.1610 0.1570## 6 -0.18000 0.00157 -0.05937405 -0.04100 -0.04790 -0.1440 -0.2910 -0.2530## F2_209 F2_210 F2_216 F2_217 F2_218 F2_219 F2_220 F2_230## 1 -0.07960 0.0848 -0.093800 -0.0898 0.0472 0.00513 0.0578 0.05616089## 2 0.14400 0.0594 0.109000 0.0791 0.2110 0.08110 0.1580 0.19241050## 3 -0.12000 -0.0627 -0.029200 0.1090 -0.0459 -0.06390 -0.1700 -0.09710876## 4 0.09110 0.0349 -0.024900 -0.0165 0.7450 0.04310 0.0427 0.38320980## 5 0.00777 0.0935 0.000275 -0.0371 0.0980 0.07460 0.2250 -0.11742250## 6 -0.11300 -0.0358 -0.042800 -0.1930 -0.1750 -0.02980 -0.1190 -0.15757000## F2_231 F2_232 F2_233 F2_234 F2_235 F2_236 F2_237 F2_238## 1 0.1470 0.018600 0.0976 0.0160 0.05150205 0.0394 0.00542 0.000242## 2 0.1410 0.056600 0.2570 0.2590 0.14049010 0.0965 0.04190 0.009570## 3 -0.0163 -0.000807 -0.1110 -0.1750 -0.09649123 0.0154 -0.00482 0.014500## 4 0.1750 -0.040400 0.0284 -0.1630 0.02090355 0.0610 0.04090 0.004970## 5 -0.0112 0.007410 0.2130 0.0578 0.06377663 -0.0739 -0.03110 0.019900## 6 -0.0319 -0.046300 -0.2130 -0.2990 -0.10599170 -0.0209 -0.14300 0.069700## F2_239 F2_249 F2_250 F2_251 F2_252 F2_254 F2_256 F2_257## 1 -0.01540 -0.02430 -0.1010 0.0626 -0.060100 0.11600 0.03889860 0.07270702## 2 0.11900 0.08050 0.1460 0.0296 0.243000 0.18900 0.13016450 0.03534575## 3 -0.00822 0.00863 -0.0533 -0.0225 0.011700 -0.19800 -0.06286667 -0.13364770## 4 0.19500 0.04790 -0.2420 0.1500 -0.000738 0.21100 0.06825731 0.04275748## 5 -0.02510 0.03110 -0.0222 NA 0.133000 -0.00411 -0.08267811 0.08027854## 6 -0.08810 -0.13200 -0.1830 -0.1090 -0.237000 -0.19800 -0.15300000 0.00877483## F2_265 F2_266 F2_268 F2_274 F2_275 F2_276 F2_279 F2_280## 1 -0.0290 0.0550 -0.0312 -0.02870776 0.05570 -0.0859 0.01570 0.1010## 2 0.0221 0.1020 0.1030 0.07293987 0.00983 0.0640 0.05220 0.2420## 3 -0.0235 -0.0451 -0.0247 -0.68900000 0.02710 -0.0721 0.00623 -0.1590## 4 0.2240 0.1280 0.0340 0.12850620 -0.09060 0.3490 -0.04130 0.0187## 5 -0.0183 -0.0851 -0.0846 -0.19800000 -0.02600 -0.1410 0.00820 -0.0193## 6 -0.0432 -0.0188 -0.1010 0.03046819 -0.05890 -0.0467 -0.10800 -0.2750## F2_281 F2_282 F2_284 F2_285 F2_286 F2_292 F2_294## 1 -0.02040 -0.00133 0.0414 0.020115580 -0.00453 0.1898726 0.04873549## 2 -0.01090 0.04050 0.0824 0.013043140 0.12100 0.0674650 -0.02203408## 3 0.00717 0.03830 0.0193 0.007803106 -0.06740 0.1602482 -0.03922225## 4 0.01140 0.05380 1.9100 -0.088830460 -0.00285 0.1820795 -0.14910580## 5 -0.12600 -0.06070 -0.0211 0.206402900 -0.01670 0.1148936 -0.02899761## 6 0.00944 -0.04300 -0.1100 -0.099250960 -0.12500 -0.1783375 -0.08796206## F2_295 F2_313 F2_314 F2_315 F2_316 F2_317 F2_318 F2_343## 1 0.01950 0.00240 -0.09950 -0.0872 -0.103662100 0.0242 0.00536 0.1340## 2 -0.01470 0.19700 0.09810 0.0618 0.098719220 0.0104 0.09670 -0.0248## 3 0.11700 -0.00744 0.00862 0.0130 -0.002592110 0.0946 0.01590 -0.0934## 4 0.14100 0.04860 -0.03720 0.7800 0.280451100 -0.0560 0.02180 0.2100## 5 0.00608 0.05360 -0.04540 -0.1290 0.001011547 0.0877 -0.07280 -0.0284## 6 -0.02930 -0.17800 -0.09560 -0.0600 -0.067627370 -0.0127 -0.07340 0.0180## F2_353## 1 0.15584910## 2 0.11533460## 3 -0.13519600## 4 0.24050990## 5 -0.13719800## 6 -0.06457439# 设置数据集数量,这里为2nSets = 2# 定义数据集名称setLabels = c("Female liver", "Male liver")shortLabels = c("Female", "Male")# 提取表达数据multiExpr <- vector(mode = "list", length = nSets) # 创建List存储数据multiExpr[[1]] <- list(data = as.data.frame(t(femData[-c(1:8)]))) # 提取femData表达数据names(multiExpr[[1]]$data) <- femData$substanceBXHrownames(multiExpr[[1]]$data) <- names(femData)[-c(1:8)] multiExpr[[2]] <- list(data = as.data.frame(t(maleData[-c(1:8)]))) # 提取maleData表达数据names(multiExpr[[2]]$data) <- maleData$substanceBXHrownames(multiExpr[[2]]$data) <- names(maleData)[-c(1:8)]# 检查数据集格式exprSize <- checkSets(multiExpr)exprSize## $nSets## [1] 2## ## $nGenes## [1] 3600## ## $nSamples## [1] 135 124## ## $structureOK## [1] TRUE# 检查基因和样本是否合格gsg <- goodSamplesGenesMS(multiExpr, verbose = 3)## Flagging genes and samples with too many missing values...## ..step 1## ..bad gene count: 0, bad sample counts: 0, 0gsg$allOK## [1] TRUE# 如果为FALSE,则运行下面的代码if (!gsg$allOK){ # 打印被移除的基因的信息: if (sum(!gsg$goodGenes) > 0) printFlush(paste("Removing genes:", paste(names(multiExpr[[1]]$data)[!gsg$goodGenes], collapse = ", "))) for (set in 1:exprSize$nSets) { if (sum(!gsg$goodSamples[[set]])) printFlush(paste("In set", setLabels[set], "removing samples", paste(rownames(multiExpr[[set]]$data)[!gsg$goodSamples[[set]]], collapse = ", "))) # 去除质量差的基因和样本 multiExpr[[set]]$data = multiExpr[[set]]$data[gsg$goodSamples[[set]], gsg$goodGenes]; } exprSize <- checkSets(multiExpr)}# 样本聚类观察离群样本sampleTrees <- list()for (set in 1:nSets){ sampleTrees[[set]] = hclust(dist(multiExpr[[set]]$data), method = "average")}# 绘制聚类树par(mfrow = c(2,1))par(mar = c(0, 4, 2, 0))for (set in 1:nSets) plot(sampleTrees[[set]], main = paste("Sample clustering on all genes in", setLabels[set]), xlab="", sub="", cex = 0.7)
# 从上图可以发现雌性小鼠数据集的最左侧有一个游离的样本,需要进行剪切处理baseHeight = 16# 设置当雌性数据集剪切时,雄性数据集对应的剪切的高度cutHeights <- c(baseHeight, baseHeight*exprSize$nSamples[2]/exprSize$nSamples[1])# 重新绘图par(mfrow = c(2,1))par(mar = c(0, 4, 2, 0))for (setin 1:nSets){ plot(sampleTrees[[set]], main = paste("Sample clustering on all genes in", setLabels[set]), xlab="", sub="", cex = 0.7); abline(h = cutHeights[set], col = "red");}
for (set in1:nSets){# 寻找被剪切的样本 labels <- cutreeStatic(sampleTrees[[set]], cutHeight = cutHeights[set])# 保留1的样本 keep <- (labels == 1) multiExpr[[set]]$data <- multiExpr[[set]]$data[keep, ]}# 再一次检查数据exprSize <- checkSets(multiExpr)exprSize## $nSets## [1] 2## ## $nGenes## [1] 3600## ## $nSamples## [1] 134 124## ## $structureOK## [1] TRUE# 读取临床表型数据traitData <- read.csv("data//ClinicalTraits.csv")dim(traitData)## [1] 361 38names(traitData)## [1] "X""Mice""Number"## [4] "Mouse_ID""Strain""sex"## [7] "DOB""parents""Western_Diet"## [10] "Sac_Date""weight_g""length_cm"## [13] "ab_fat""other_fat""total_fat"## [16] "comments""X100xfat_weight""Trigly"## [19] "Total_Chol""HDL_Chol""UC"## [22] "FFA""Glucose""LDL_plus_VLDL"## [25] "MCP_1_phys""Insulin_ug_l""Glucose_Insulin"## [28] "Leptin_pg_ml""Adiponectin""Aortic.lesions"## [31] "Note""Aneurysm""Aortic_cal_M"## [34] "Aortic_cal_L""CoronaryArtery_Cal""Myocardial_cal"## [37] "BMD_all_limbs""BMD_femurs_only"# 去除无关的列allTraits <- traitData[, - c(31, 16)]allTraits <- allTraits[, c(2, 11:36)]dim(allTraits)## [1] 361 27names(allTraits)## [1] "Mice""weight_g""length_cm"## [4] "ab_fat""other_fat""total_fat"## [7] "X100xfat_weight""Trigly""Total_Chol"## [10] "HDL_Chol""UC""FFA"## [13] "Glucose""LDL_plus_VLDL""MCP_1_phys"## [16] "Insulin_ug_l""Glucose_Insulin""Leptin_pg_ml"## [19] "Adiponectin""Aortic.lesions""Aneurysm"## [22] "Aortic_cal_M""Aortic_cal_L""CoronaryArtery_Cal"## [25] "Myocardial_cal""BMD_all_limbs""BMD_femurs_only"# 构建具有表型特征的多数据集结构Traits <- vector(mode="list", length = nSets)for (set in1:nSets){ setSamples <- rownames(multiExpr[[set]]$data) traitRows <- match(setSamples, allTraits$Mice) Traits[[set]] <- list(data = allTraits[traitRows, -1]) rownames(Traits[[set]]$data) <- allTraits[traitRows, 1]}dim(Traits[[1]][["data"]])## [1] 134 26dim(Traits[[2]][["data"]])## [1] 124 26nGenes <- exprSize$nGenesnSamples <- exprSize$nSamples# 保存数据save(multiExpr, Traits, nGenes, nSamples, setLabels, shortLabels, exprSize, file = "Consensus-dataInput.RData")
构建网络并进行模块识别
对于WGCNA来说,最重要的就是构建网络并识别出相应的模块,将基因按照相似性聚集成相应的模块,这里我们一样适用自动识别的方面来构建网络和识别模块:
# 开启多线程并行运算enableWGCNAThreads()## Allowing parallel execution with up to 11 working processes.# 加载数据lnames <- load(file = "Consensus-dataInput.RData");lnames## [1] "multiExpr" "Traits" "nGenes" "nSamples" "setLabels" ## [6] "shortLabels" "exprSize"# 获取multiExpr中的数据集数量nSets <- checkSets(multiExpr)$nSets# 设置软阈值powers <- c(seq(4,10,by=1), seq(12,20, by=2))powerTables <- vector(mode = "list", length = nSets)# 拓扑异构分析for (setin 1:nSets) powerTables[[set]] <- list(data = pickSoftThreshold(multiExpr[[set]]$data, powerVector=powers, verbose = 2)[[2]])## pickSoftThreshold: will use block size 3600.## pickSoftThreshold: calculating connectivity for given powers...## ..working on genes 1 through 3600 of 3600## Power SFT.R.sq slope truncated.R.sq mean.k. median.k. max.k.## 1 4 0.506 -1.42 0.973 56.50 47.2000 202.0## 2 5 0.681 -1.72 0.940 32.20 25.1000 134.0## 3 6 0.902 -1.50 0.962 19.90 14.5000 94.8## 4 7 0.921 -1.67 0.917 13.20 8.6800 84.1## 5 8 0.904 -1.72 0.876 9.25 5.3900 76.3## 6 9 0.859 -1.70 0.836 6.80 3.5600 70.5## 7 10 0.833 -1.66 0.831 5.19 2.3800 65.8## 8 12 0.853 -1.48 0.911 3.33 1.1500 58.1## 9 14 0.876 -1.38 0.949 2.35 0.5740 51.9## 10 16 0.907 -1.30 0.970 1.77 0.3090 46.8## 11 18 0.912 -1.24 0.973 1.39 0.1670 42.5## 12 20 0.931 -1.21 0.977 1.14 0.0951 38.7## pickSoftThreshold: will use block size 3600.## pickSoftThreshold: calculating connectivity for given powers...## ..working on genes 1 through 3600 of 3600## Power SFT.R.sq slope truncated.R.sq mean.k. median.k. max.k.## 1 4 0.577 -1.58 0.938 51.30 45.1000 196.0## 2 5 0.630 -1.68 0.926 29.10 23.4000 132.0## 3 6 0.861 -1.36 0.968 17.90 12.9000 91.8## 4 7 0.903 -1.52 0.932 11.90 7.5200 79.6## 5 8 0.908 -1.61 0.889 8.28 4.5600 72.0## 6 9 0.895 -1.62 0.867 6.07 2.8800 66.2## 7 10 0.882 -1.59 0.868 4.63 1.9000 61.3## 8 12 0.890 -1.49 0.904 2.96 0.8760 53.5## 9 14 0.931 -1.38 0.958 2.09 0.4370 47.3## 10 16 0.954 -1.31 0.980 1.57 0.2270 42.1## 11 18 0.970 -1.26 0.984 1.24 0.1190 37.7## 12 20 0.984 -1.22 0.995 1.01 0.0642 33.9# 结果可视化colors <- c("black", "red")plotCols = c(2,5,6,7)colNames <- c("Scale Free Topology Model Fit", "Mean connectivity", "Median connectivity","Max connectivity")ylim <- matrix(NA, nrow = 2, ncol = 4)for (setin 1:nSets){for (col in 1:length(plotCols)) { ylim[1, col] = min(ylim[1, col], powerTables[[set]]$data[, plotCols[col]], na.rm = TRUE) ylim[2, col] = max(ylim[2, col], powerTables[[set]]$data[, plotCols[col]], na.rm = TRUE) }}par(mfcol = c(2,2));par(mar = c(4.2, 4.2 , 2.2, 0.5))cex1 = 0.7;for (col in 1:length(plotCols)) for (setin 1:nSets){if (set==1) { plot(powerTables[[set]]$data[,1], -sign(powerTables[[set]]$data[,3])*powerTables[[set]]$data[,2], xlab="Soft Threshold (power)",ylab=colNames[col],type="n", ylim = ylim[, col], main = colNames[col]); addGrid(); }if (col==1) { text(powerTables[[set]]$data[,1], -sign(powerTables[[set]]$data[,3])*powerTables[[set]]$data[,2], labels=powers,cex=cex1,col=colors[set]); } else text(powerTables[[set]]$data[,1], powerTables[[set]]$data[,plotCols[col]], labels=powers,cex=cex1,col=colors[set]);if (col==1) { legend("bottomright", legend = setLabels, col = colors, pch = 20) ; } else legend("topright", legend = setLabels, col = colors, pch = 20) ;}
# 根据图片结果,我们选择6为两个数据集的最佳软阈值,接下来进行网络构建net <- blockwiseConsensusModules( multiExpr, power = 6, minModuleSize = 30, deepSplit = 2, pamRespectsDendro = FALSE, mergeCutHeight = 0.25, numericLabels = TRUE, minKMEtoStay = 0, saveTOMs = TRUE, verbose = 5)## Calculating consensus modules and module eigengenes block-wise from all genes## Calculating topological overlaps block-wise from all genes## Flagging genes and samples with too many missing values...## ..step 1## ..bad gene count: 0, bad sample counts: 0, 0## ....Working on set 1## TOM calculation: adjacency..## ..will not use multithreading.## Fraction of slow calculations: 0.396405## ..connectivity..## ..matrix multiplication (system BLAS)..## ..normalization..## ..done.## ....Working on set 2## TOM calculation: adjacency..## ..will not use multithreading.## Fraction of slow calculations: 0.438385## ..connectivity..## ..matrix multiplication (system BLAS)..## ..normalization..## ..done.## Cluster size 3600 broken into 2108 1492 ## Cluster size 2108 broken into 1126 982 ## Done cluster 1126 ## Done cluster 982 ## Done cluster 2108 ## Done cluster 1492 ## Cluster size 3600 broken into 1576 2024 ## Cluster size 1576 broken into 694 882 ## Done cluster 694 ## Done cluster 882 ## Done cluster 1576 ## Cluster size 2024 broken into 1045 979 ## Done cluster 1045 ## Done cluster 979 ## Done cluster 2024 ## ..Working on block 1 .## ....Working on set 1## ....Working on set 2## ....Calculating consensus network## ..Working on block 1 .## ....clustering and detecting modules..## ..done.## ....calculating eigengenes..## multiSetMEs: Calculating module MEs.## Working on set 1 ...## Working on set 2 ...## ....checking consensus modules for statistical meaningfulness..## ....checking for genes that should be reassigned..## ......module membership p-values..## ......module membership scores..## ......individual modules.. ## ..merging consensus modules that are too close..## mergeCloseModules: Merging modules whose distance is less than 0.25## multiSetMEs: Calculating module MEs.## Working on set 1 ...## moduleEigengenes: Calculating 20 module eigengenes in given set.## Working on set 2 ...## moduleEigengenes: Calculating 20 module eigengenes in given set.## multiSetMEs: Calculating module MEs.## Working on set 1 ...## moduleEigengenes: Calculating 18 module eigengenes in given set.## Working on set 2 ...## moduleEigengenes: Calculating 18 module eigengenes in given set.## Calculating new MEs...## multiSetMEs: Calculating module MEs.## Working on set 1 ...## moduleEigengenes: Calculating 18 module eigengenes in given set.## Working on set 2 ...## moduleEigengenes: Calculating 18 module eigengenes in given set.# 保留模块信息consMEs <- net$multiMEsmoduleLabels <- net$colors# 将模块标签转为颜色标签moduleColors <- labels2colors(moduleLabels)consTree <- net$dendrograms[[1]]# 可视化plotDendroAndColors(consTree, moduleColors, "Module colors", dendroLabels = FALSE, hang = 0.03, addGuide = TRUE, guideHang = 0.05, main = "Consensus gene dendrogram and module colors")
# 保存数据save(consMEs,moduleLabels, moduleColors,consTree, file = "Consensus-NetworkConstruction-auto.RData")
将共识模块与雌性小鼠特定模块相关联
# 加载上期推文得到的雌性小鼠的WGCNA模块数据数据lnames <- load("FemaleLiver-02-networkConstruction-auto.RData")lnames## [1] "MEs""moduleLabels""moduleColors""geneTree"# 重命名变量femaleLabels <- moduleLabelsfemaleColors <- moduleColorsfemaleTree <- geneTreefemaleMEs <- orderMEs(MEs,                      greyName = "ME0")# 加载上一步分析得到的共识性模块lnames <- load("Consensus-NetworkConstruction-auto.RData")lnames## [1] "consMEs""moduleLabels""moduleColors""consTree"# 按照模块标签在有序模块特征基因中出现的顺序进行分离femModuleLabels <- substring(names(femaleMEs), 3)femModuleLabels## [1] "9""7""1""5""16""2""3""6""4""13""18""11""17""8""10"## [16] "12""14""15""0"consModuleLabels <- substring(names(consMEs[[1]]$data), 3)consModuleLabels## [1] "13""7""8""11""16""15""17""12""9""3""14""4""6""1""10"## [16] "2""5""0"# 转换颜色标签femModules <- labels2colors(as.numeric(femModuleLabels))femModules## [1] "magenta""black""turquoise""green""lightcyan"## [6] "blue""brown""red""yellow""salmon"## [11] "lightgreen""greenyellow""grey60""pink""purple"## [16] "tan""cyan""midnightblue""grey"consModules <- labels2colors(as.numeric(consModuleLabels))consModules## [1] "salmon""black""pink""greenyellow""lightcyan"## [6] "midnightblue""grey60""tan""magenta""brown"## [11] "cyan""yellow""red""turquoise""purple"## [16] "blue""green""grey"# 雌性小鼠模块和共识模块的数量nFemMods <- length(femModules)nFemMods## [1] 19nConsMods <- length(consModules)nConsMods## [1] 18# 初始化p值表和相应计数表pTable <- matrix(0, nrow = nFemMods, ncol = nConsMods);CountTbl <- matrix(0, nrow = nFemMods, ncol = nConsMods);# 进行成对比较for (fmod in1:nFemMods)for (cmod in1:nConsMods) { femMembers <- (femaleColors == femModules[fmod]) consMembers <- (moduleColors == consModules[cmod]) pTable[fmod, cmod] <- -log10(fisher.test(femMembers, consMembers, alternative = "greater")$p.value) CountTbl[fmod, cmod] <- sum(femaleColors == femModules[fmod] & moduleColors == consModules[cmod]) }# 截断p值pTable[is.infinite(pTable)] <- 1.3*max(pTable[is.finite(pTable)])pTable[pTable > 50] = 50# 计算各模块大小femModTotals <- apply(CountTbl, 1, sum)femModTotals## [1] 123 211 609 312 58 460 409 221 316 91 34 100 47 157 106 94 77 76 99consModTotals <- apply(CountTbl, 2, sum)consModTotals## [1] 82 245 240 105 44 62 37 89 236 292 63 279 251 423 183 353 266 350# 可视化分析par(mfrow=c(1,1))par(cex = 1.0)par(mar=c(8, 10.4, 2.7, 1)+0.3)# Use function labeledHeatmap to produce the color-coded table with all the trimmingslabeledHeatmap(Matrix = pTable, xLabels = paste(" ", consModules), yLabels = paste(" ", femModules), colorLabels = TRUE, xSymbols = paste("Cons ", consModules, ": ", consModTotals, sep=""), ySymbols = paste("Fem ", femModules, ": ", femModTotals, sep=""), textMatrix = CountTbl, colors = blueWhiteRed(100)[50:100], main = "Correspondence of Female set-specific and Female-Male consensus modules", cex.text = 1.0, cex.lab = 1.0, setStdMargins = FALSE)
# 从热图可以看出,上期推文我们选择的感兴趣模块——蓝色模块,其对应的共识模块是黄色和红色两个模块,那么共识模块中的黄色和红色两个模块的基因可以作为我们后续关注的感兴趣模块
关联共识模块和表型数据
如果前面没有对单独数据集先进行WGCNA分析,那我们也可以和单个数据集的WGCNA分析一样,把识别到的共识模块和临床表型信息进行关联,得到和感兴趣的表型信息相关性最大的共识模块,选择其作为感兴趣模块
# 加载前面的网络信息lnames <- load(file = "Consensus-NetworkConstruction-auto.RData")lnames## [1] "consMEs" "moduleLabels" "moduleColors" "consTree"exprSize <- checkSets(multiExpr)exprSize## $nSets## [1] 2## ## $nGenes## [1] 3600## ## $nSamples## [1] 134 124## ## $structureOK## [1] TRUEnSets <- exprSize$nSets# 建立变量以包含模块-特质相关性moduleTraitCor = list()moduleTraitPvalue = list()# 计算相关性for (setin 1:nSets){ moduleTraitCor[[set]] = cor(consMEs[[set]]$data, Traits[[set]]$data, use = "p") moduleTraitPvalue[[set]] = corPvalueFisher(moduleTraitCor[[set]], exprSize$nSamples[set])}MEColors <- labels2colors(as.numeric(substring(names(consMEs[[1]]$data), 3)))MEColorNames <- paste("ME", MEColors, sep="")# 可视化分析# 第一个数据集set = 1textMatrix <- paste(signif(moduleTraitCor[[set]], 2), "\n(", signif(moduleTraitPvalue[[set]], 1), ")", sep = "");dim(textMatrix) <- dim(moduleTraitCor[[set]])par(mar = c(6, 8.8, 3, 2.2));labeledHeatmap(Matrix = moduleTraitCor[[set]], xLabels = names(Traits[[set]]$data), yLabels = MEColorNames, ySymbols = MEColorNames, colorLabels = FALSE, colors = blueWhiteRed(50), textMatrix = textMatrix, setStdMargins = FALSE, cex.text = 0.5, zlim = c(-1,1), main = paste("Module--trait relationships in", setLabels[set]))
# 第二个数据集set = 2textMatrix <- paste(signif(moduleTraitCor[[set]], 2), "\n(", signif(moduleTraitPvalue[[set]], 1), ")", sep = "");dim(textMatrix) <- dim(moduleTraitCor[[set]])par(mar = c(6, 8.8, 3, 2.2));labeledHeatmap(Matrix = moduleTraitCor[[set]], xLabels = names(Traits[[set]]$data), yLabels = MEColorNames, ySymbols = MEColorNames, colorLabels = FALSE, colors = blueWhiteRed(50), textMatrix = textMatrix, setStdMargins = FALSE, cex.text = 0.5, zlim = c(-1,1), main = paste("Module--trait relationships in", setLabels[set]))
# 初始化矩阵以保存相关性和p值consensusCor<- matrix(NA, nrow(moduleTraitCor[[1]]), ncol(moduleTraitCor[[1]]))consensusPvalue<- matrix(NA, nrow(moduleTraitCor[[1]]), ncol(moduleTraitCor[[1]]))# 找到两个数据集一致的负相关信息negative<- moduleTraitCor[[1]] < 0 & moduleTraitCor[[2]] < 0consensusCor[negative]<- pmax(moduleTraitCor[[1]][negative], moduleTraitCor[[2]][negative])consensusPvalue[negative]<- pmax(moduleTraitPvalue[[1]][negative],                                 moduleTraitPvalue[[2]][negative])# 找到两个数据集一致的正相关信息positive<- moduleTraitCor[[1]] > 0 & moduleTraitCor[[2]] > 0consensusCor[positive]<- pmin(moduleTraitCor[[1]][positive], moduleTraitCor[[2]][positive])consensusPvalue[positive]<- pmax(moduleTraitPvalue[[1]][positive],                                 moduleTraitPvalue[[2]][positive])textMatrix<- paste(signif(consensusCor, 2), "\n(",signif(consensusPvalue,1), ")", sep = "");dim(textMatrix) <- dim(moduleTraitCor[[set]])# 关联两个数据集的模块相关性结果进行可视化par(mar = c(6, 8.8, 3, 2.2));labeledHeatmap(Matrix = consensusCor,xLabels = names(Traits[[set]]$data),yLabels = MEColorNames,ySymbols = MEColorNames,colorLabels = FALSE,colors = blueWhiteRed(50),textMatrix = textMatrix,setStdMargins = FALSE,cex.text = 0.5,zlim = c(-1,1),main = paste("Consensus module--trait relationships across\n",paste(setLabels,collapse = " and ")))
# 这样我们通过热图就可以直接获得我们感兴趣的模块信息了
好了,今天的内容到此结束!今天我们一起学习的是两个数据集的的WGCNA共识性分析,这比单个数据集的WGCNA更难,也提供了更多的分析可能性。
比如我们手上有两份数据集,一份来自TCGA,另一份来自GEO,那么我们就可以利用本期使用的方法进行两个数据集的共识性分析,找到两个数据集的共识性模块并关联感兴趣表型,从而筛选感兴趣模块得到目标基因集,这样比我们使用单个数据集的分析更加准确;再比如,如果我们要进行"关于不同性别患者患肺癌"的WGCNA分析,我们可以仿照本推文的做法,将肺癌数据集按照性别一分为二,然后进行WGCNA分析,找到男性患者和女性患者的共识性模块,得到共识性基因,找到男性和女性肺癌患者中跟吸烟最相关的基因。
当然还可以有更多的玩法,大家可以自行挖掘,也欢迎留言讨论。
后台回复“feng 61”获取代码和数据,我们下期见!
单细胞分析专栏传送门
碎碎念专栏传送门(完结)
风之美图系列传送门(完结)
END

撰文丨风   风
排版丨四金兄
主编丨小雪球
欢迎大家关注解螺旋生信频道-挑圈联靠公号~

继续阅读
阅读原文