解螺旋公众号·陪伴你科研的第2618天

非肿瘤多数据集联合分析

在进行生信分析的时候，如果搜索到好几个GEO数据集，不知道大家有没有这些困惑：

多个GEO数据集，应该分别分析，还是合并以后再一起分析呢？

哪些数据集能合并哪些不能合并呢？

合并数据需要怎样特殊的处理吗？

这些问题困扰我很久，在雪球老师和阿琛老师的指导下，终于有了一点心得。今天就让我们通过这篇2021年1月发表在Journal of Translational Medicine（IF：5.531）上的文章，一起来看看多个数据集的分析思路吧。

期刊简介

文章标题

Three hematologic/immune system-specific expressed genes are considered as the potential biomarkers for the diagnosis of early rheumatoid arthritis through bioinformatics analysis

要素拆解

题目：通过生物信息学分析指出3种血液学/免疫系统特异性表达基因是诊断早期类风湿性关节炎的潜在生物标志物

疾病：类风湿性关节炎（rheumatoid arthritis, RA），骨关节炎（osteoarthritis, OA）

数据来源：GEO （GSE77298，GSE82107，GSE55584，GSE55457，GSE55235，GSE89408）

分析策略：表达差异 + 功能聚类 + 交互网络 + 临床意义

背景知识

GEO数据来源非常丰富，既有芯片数据也有测序数据，即使是芯片数据也有多个芯片平台，所以一般来说，数据集不能直接合并分析。

小贴士：想了解更多，推荐学习解螺旋生信全书上篇段位一的内容

在合并数据集的时候最好遵循以下原则：

1. 同一平台的数据合并分析的效果更好。GEO的数据集有不同的平台，比如本文的GEO数据就涉及GPL570,GPL96,GPL11154这个3个平台，一般来说，同一个平台的数据合并分析的效果较好，所以作者选择合并GSE77298和GSE82107这两个数据集，因为它们都来自GPL570平台。此外，作者还合并来自GPL96平台的3个数据集作为验证集。当然，不同的平台也能合并，但是程序更复杂，需要前期的数据处理也更多，如果有条件，还是尽可能合并同一平台的数据。

2. 尽可能合并来自同一组织的样本。本文采用的都是来自关节滑膜的数据，那么看看其他文章。比如下面这篇文章中，虽然前2个数据集GSE8294和GSE22255均来自GPL570平台，但其样本一个是来自全血，一个是来自外周血单核细胞，所以更推荐合并后2个数据集GSE16561和GSE37587,这2个数据集不仅来自同一平台，而且样本都是来自全血组织。

3. 合并数据集时，需要去批次效应。由于各种GEO数据集可能有不同的实验仪器、不同的实验试剂、不同的实验人员等，会导致不同数据集的批次差异，这也是为什么尽可能选择同一平台同一组织来源的数据集合并效果会更好。无论合并什么样的数据，都需要去批次效应，以减少外部因素对数据的影响。

小贴士：想了解更多，推荐学习解螺旋生信全书下篇段位一的内容

4. 可以合并分析，也可以分别分析以后取交集。多个数据集分析时，其实非常灵活，没有严格规定必须合并分析或者分开分析，可以尝试不同的分析方法，看看哪种结果更利于把文章的故事说清楚。如果多个数据集合并分析的结果不理想或者不能满足以上条件时，分别分析以后取交集也是不错的选择。在过往的推文（←点击链接查看）中我们就遇到过因数据集来自不同平台，故作者分别分析了各个数据集的差异基因，然后取了交集的情况。

这些是我的一点心得，大家还有什么合并数据集的秘诀吗？期待在评论区看到你的分享哦〜

数据解读

本文一共有8个图3个表。作者合并来自GPL570平台的2个数据集GSE77298和GSE82107作为测试集，其中包括16个RA样本和10个OA样本；将4个数据集GSE55584、GSE55457、GSE55235、GSE89408数据集作为验证集(表1)。分析差异表达基因（图1），利用BioGPS数据库鉴定的组织/器官的特异性表达基因分布（表2），进行GSEA（图2）和GO/KEGG（图3）富集分析、PPI网络分析（图4）、筛选hub基因（表3和图4），构建mRNA－miRNA共表达网络（图5），利用验证集来验证hub基因的表达差异（图6）和ROC曲线（图7），最终筛选出3个关键基因的ceRNA（图8）。

复现工具

◆ 仙桃学术工具

（https://www.xiantao.love/products）

◆ BioGPS数据库

（http://biogps.org/）

◆ NetworkAnalyst数据库

（https://www.networkanalyst.ca/）

◆ Cytoscape软件＋ cytoHubba插件 + MCODE插件

学术造假

因表1为从GEO平台GEO datasets平台中获得的数据集信息整理所得，故不做具体复现展示。我们从figure1开始此次复现之旅。

Figure 1

比较RA和OA的差异表达基因

仙桃学术（https://www.xiantao.love/），点击【数据集检索】

小贴士：免费版/基础版/高级版每日可提交的分析次数不同，根据需要提升等级吧~

输入数据集名称【GSE77298】 → 【检索】 → 【选择样本】

勾选RA样本 → 【添加到样本库】

输入数据集名称【GSE82107】 → 【检索】 → 【选择样本】

勾选OA样本 → 【添加到样本库】 → 【进入我的样本库】

将【OA】勾选【加入参考组】 → 将【RA】勾选【加入实验组】 → 【提交分析】

待状态变为【完成】时即可下载结果

很神奇的是，仙桃已经默默的完成了2个数据集的合并、数据校正、去批次效应、差异基因分析、热图火山图的绘制。

这里直接【热图下载】和【火山图下载】即为图1

【CSV表格下载】可下载所有分析结果，保存为“差异分析.csv”

用Excel打开下载的CSV表格，筛选log2 FC > 1或 < -1，并且p<0.05的差异基因，另存为新文件“DEG.xlsx”以供后续分析。

小贴士：一般建议是用adj.p.value<0.05来寻找差异基因，如果使用adj.p.value的差异分析的结果不是那么好，这里我们用p值来替代继续进行后续的分析

Table 2

BioGPS鉴定的组织/器官特异性表达基因的分布

进入BioGPS数据库(http://biogps.org)，分别将“DEG.xlsx”中的基因输入数据库，然后查找在人体组织中，同样出现特异表达的基因。我们以CXCL13为例。

输入【CXCL13】 → 【Search】

点击【CXCL13】

可以看到，CXCL13在扁桃体和淋巴结中高表达。

同样的操作，将其他基因输入，最终鉴定了23个在特定组织或器官系统中表达的基因，并整理成表2。

Figure 2

GSEA 图显示了 RA 组和 OA 组中与免疫相关基因集

回到仙桃学术（https://www.xiantao.love/） → 【生信工具】

【高级版】 → 【立即使用】

注：免费版和基础版都可以进行统计和可视化，由于高级版功能最全，这里选择高级版作为范例

将“差异分析.csv”整理成【Gene.Symbol】和【logFC】两列，另存为“GSEAinput.xlsx”

【功能聚类（圈）】 → 【GSEA富集】 → 【GSEA分析】 → 上传“GSEAinput.xlsx” → 【确认】

在【历史记录】中状态为【完成】时即可下载结果

【分析工具】 → 【功能聚类（圈）】 → 【GSEA富集】 → 【GSEA可视化】 → 选择刚刚的【GSEA富集分析结果】 → 【确认】 → 即可【保存结果】或【下载图片】

Figure 3

GO/KEGG分析

【功能聚类（圈）】 → 【GO|KEGG】 → 【GO|KEGG富集分析】 → 上传“GSEAinput.xlsx” → 选择“GO:BP” → 【确认】 → 【保存结果】

小贴士：【GO|KEGG富集分析】需要上传文件的标准格式只需要基因名一列，这里偷了个小懒，直接上传“GSEAinput.xlsx”也是可以的；这里一定要记得【保存结果】，这样就能直接进入下一步【GO|KEGG可视化】了~

【GO|KEGG可视化】 → 选择刚刚分析的【BP】 → 【确认】 → 【保存结果】或【下载图片】即为图3A

将【GO|KEGG富集分析】中选择“KEGG”并【GO|KEGG可视化】，即为图3B

Table 3

cytoHubba识别出的15个hub基因

【交互网络（联）】 → 【STRING蛋白互作】 → 【PPI蛋白互作分析】 → 上传“GSEAinput.xlsx” → 【确认】 → 【Excel表格下载】

将基因对应的上调or下调属性列为一个属性列表,保存为“attribute.xlsx”

打开Cytoscape软件 → 【File】 → 【Import】 → 【Network from File】 → 上传刚刚下载的Excel表格“STRING蛋白互作.xlsx”

【File】 → 【Import】 → 【Table from File】 → 上传刚刚保存的Excel表格“attribute.xlsx”

【cytoHubba】 → 【Calculate】 → 【Top15】 → 【Check the first-stage nodes】 → 【Submit】 → 【Save Current Rank】

根据这15个基因，结合前面保存的“DEG.xlsx”，可整理出表3。

Figure 4

由 MCODE 提取的 DEG 和四个集群模块的 PPI 网络

仍然在cytoHubba中，【Apps】 → 【MCODE】

调整好分析参数后【Analyze Current Network】

保存图片即为图4A

这时形成了4个cluster，选择第一个cluster → 【Select】 → 【Nodes】 → 【Hide Unselected Nodes】，保存图片即为图4B

图4C-4E制作方法相似。将cytoHubba得到的hub基因与MCODE得到的hub基因取交集，便得到8个hub基因：GZMA，TTK，DLGAP5，CD27，PRC1,CD52,CXCL13,CEP55。

Figure 5

mRNA 和miRNA 的共表达网络

作者使用了五个在线 miRNA 数据库，即 RNA22、DIANA-micro T、miRWalk、miRDB 和 miRcode，来预测hub基因的靶向miRNA。但因某些数据库目前连不上，这里介绍一种简便方法。

进入NetworkAnalyst数据库（https://www.networkanalyst.ca/） → 【Gene List Input】

选择【H.sapiens(human)】 → 【Official Gene Symbol】 → 输入8个hub基因名 → 【Upload】 → 【Proceed】

选择【Gene-miRNA Interactions】 → 【miRTarBase v8.0】 → 【OK】 → 【Proceed】

【Proceed】

调整图片格式即为图5

Figure 6

GEO数据库的4个数据集验证了8个特异表达的hub基因

按图1的方法，合并GSE55584,GSE55457,GSE55235,进行差异分析后，下载【表达谱.csv】

打开【表达谱.csv】，查找CD27，复制此行

新建一个Excel表格，【粘贴】，【转置】

根据样本名，将标本分为“RA”和“OA”，保存为“ROC.xlsx”,可用于图7的分析。

将“OA”整理为一列，“RA”整理为一列，保存文件。

回到仙桃学术生信工具（https://www.xiantao.love/products/） → 【基础绘图】 → 【分组比较图】 → 上传刚刚保存的文件 → 【确认】 → 【保存结果】或【下载图片】

图6的其他小图做法类似。

Figure 7

8个特异表达的hub基因的ROC曲线

打开“ROC.xlsx”,整理为如下格式，保存。

回到仙桃学术生信工具（https://www.xiantao.love/products/） → 【基础绘图】 → 【诊断性ROC－独立指标】 → 上传刚刚保存的“ROC.xlsx”文件 → 【确认】 → 【保存结果】或【下载图片】

图７的其他小图做法类似。

Figure 8

PRC1、TTK 和 GZMA 的三个 ceRNA 网络以及潜在的 RNA 调控途径

根据图６和图７的分析，最终筛选出３个关键基因PRC1, TTK, GZMA。

利用StarBase（3.0版）（http://starbase.sysu.edu.cn/index.php）用于预测与选定miRNA相互作用的lncRNA和circRNA。最近StarBase正在进行数据维护，没办法打开。等到StarBase再次开放时，可参考一本八道老师的推文“starBase3.0王者归来，研究RNA必备神器”，来制作图8。

本文亮点

1. 本文一共采用了6个数据集，在数据集的选择上也是花了心思的

2. 在合并数据集方面，本文做得非常规范，严格按照对同一平台同一组织来源的数据集进行合并

3. 对于不同平台的数据集，本文将其分别当作测试集和验证集，而且将2个不同平台的数据集当作2个验证集，让验证更可靠

文章扩展

1. 仙桃学术的数据集检索功能已经非常强大了，目前已经收录了近2万个表达谱数据，对GEO数据的一站式分析非常友好，特别是非肿瘤研究者的福音。

2. 在数据集合并的时候，最好在附图中展示合并校正前后的箱式图、PCA图、UMAP图，这样能证明数据集合并时去除了批次效应，并且尽可能的减小了数据集合并的差异。贴心的技术小哥哥已经为大家想到了。

细心的小伙伴也已经发现，在差异分析的结果中已经有校正前后的箱式图、PCA图、UMAP图，直接下载就能使用了。

—END—

撰文丨dodo

排版丨豨莶

临床论文发表的必修技：Clinicaltrials.gov 注册

Cancer Cell：明星分子p53调控肿瘤发生的新机制

本周直播预告

领悟科研优人一步

直播主题：新影响因子发布，有哪些值得关注的地方？

直播时间：7月8日晚18点-20点

本次直播将在哔哩哔哩和微信视频号同步播出，同时也会在双平台进行直播抽奖。大家可以去微信视频号进行直播预约，更多直播精彩内容不要错过奥。

戳按钮预约直播

本站仅提供存储服务，所有内容均由用户发布，如发现有害或侵权内容，请点击举报。