前两天有个学生跑来问我，手里握着四十多个差异表达基因，问我怎么从中找出真正在肿瘤里起核心作用的那几个。这个问题几乎每个做癌症组学的人都会遇到，而答案往往不是再盯着单个基因死磕，而是把这些基因放进一张网络里去看。所谓的癌症基因网络，就是把基因当作节点、把基因之间的蛋白互作或调控关系当作边，用开源工具构建出一张可分析、可验证的关系图。这两年我帮团队跑了不少类似的分析，从TCGA表达数据到STRING互作网络，再到Cytoscape的可视化和模块识别，整套链路基本都能用免费开源工具走通，而且每一步的结论都可复现。这篇文章就把我实际操作中验证过的方法、工具选型以及那些容易翻车的地方完整写出来。
1. 从"单基因"到"基因网络"：癌症研究为什么必须换一种思路
1.1 单基因研究的天花板
以前我们做癌症分子机制，最常见的一套逻辑是：先通过差异表达筛选出一批候选基因，然后挨个做功能实验，验证哪个基因能促进增殖、哪个能抑制凋亡。这种"一个基因一个坑"的打法在二十年前是主流，现在也仍然有用，但问题很明显——肿瘤不是一个基因失控，而是整个调控系统出了问题。你筛出来的那批差异基因，单独看每个都挺有道理，放在一起却不知道他们之间谁管谁、谁和谁抱团、哪几个属于同一条通路。
我自己接过不少类似的求助，很多人的分析卡在同一个地方：手上有差异基因列表，有表达矩阵，也有临床信息，但就是不知道下一步该用什么方法把这堆基因串成线索。这时候把基因放到网络里去看，往往是性价比最高的突破方式。
1.2 网络视角能解决什么实际问题
构建基因网络的核心价值有三个。第一是找关键节点，也就是网络里连接度特别高的hub基因，这类基因往往在功能上更重要，也更容易成为潜在的标志物或治疗靶点；第二是找功能模块，网络里会自然形成一些联系紧密的小团体，这些模块常常对应的就是特定生物学过程或信号通路；第三是做降维和优先级排序，几千个差异基因不可能逐个做实验，但通过网络分析能把范围缩到几十个基因，直接给下游湿实验指明方向。
举一个常见的例子。你从TCGA的肺癌数据里筛出800个差异基因，直接看KEGG富集结果会觉得什么通路都有，但如果你先用PPI网络把这些基因串起来，再用MCODE跑一遍模块识别，很可能会发现最大的那个模块集中在细胞周期和DNA复制上，而模块里的核心基因就那么七八个。这个信息对后续实验设计是决定性的。
1.3 为什么这个领域几乎被开源工具主导
市面上做基因网络分析的商业软件不是没有，但生信圈的主流实践基本都跑在开源工具上。原因很直接：癌症基因网络分析这种工作，数据格式五花八门，分析流程经常要按具体研究问题反复调整，闭源软件很难跟上这种灵活性。更重要的是，审稿人和合作方都要求结果可复现，开源工具从数据下载到每一步参数设置都能完整记录，别人照着跑一遍能得到同样结果。再加上R语言生态和Cytoscape插件体系已经把常用功能覆盖得很完整，实在没有必要去用那些又贵又封闭的方案。
2. 数据起跑线：表达谱、互作组、临床注释数据源怎么选
2.1 表达数据：TCGA和GEO的使用场景完全不同
做基因网络分析，表达数据是第一层输入。TCGA和GEO是目前最常用的两个来源，但它们服务于不同场景。
TCGA的优势在于泛癌、多组学、临床信息齐全。如果你想做"某个癌种的共表达网络 + 模块与生存预后关联"，首选TCGA。TCGA的转录组数据是HTSeq-Counts或者HTSeq-FPKM，建议用STAR处理后的Counts做差异分析，或者用log2(TPM+1)做共表达网络。下载方式我推荐用R的TCGAbiolinks包，它比直接在GDC官网手动点要顺手得多，而且能直接拿到整理好的临床信息。
GEO则更适合做外部验证和跨平台复现。比如你用TCGA数据构建了一个共表达模块，需要在独立队列里确认这个模块仍然成立，这时候去GEO下载一个相同癌种的数据集，用GEOquery包读取表达矩阵，对齐相同基因后重新计算模块特征基因表达，就能做验证。需要注意GEO里不同平台的数据格式差异很大，有affy芯片的CEL文件、有illumina的txt矩阵、也有已经标准化好的Series Matrix，拿到手之后务必确认基因注释的版本，否则后续ID转换会给你添很多麻烦。
2.2 互作数据：STRING、BioGRID、GeneMANIA各有各的脾气
构建蛋白互作网络时，互作数据来源决定了网络的基础质量。我日常最常用的是STRING数据库。STRING的特点是把实验验证、数据库注释、共表达、文本挖掘等多种证据整合成一条combined score，分数越高代表互作越可靠。在线版操作简单，也提供API和R包，适合快速构建某个基因集的互作网络。
BioGRID的定位和STRING不一样。它偏向收录经过实验验证的互作关系，包括酵母双杂交、亲和纯化质谱等方法得到的直接物理互作。如果你要构建一个高可信度的核心网络，BioGRID是更好的底库，它不会把"文本挖掘出来的可能关联"也算进去，代价是覆盖范围比STRING小。GeneMANIA则偏功能关联网络，适合回答"这批基因是否参与同一生物学功能"这类问题，做出来的网络边代表功能关联而非物理结合，解读时需要区分开。
实际分析中建议以STRING为主、BioGRID为辅，先看两个来源之间的交集互作——那些在两个数据库里同时出现的互作关系通常更可靠，也更容易说服审稿人。
2.3 下载和预处理里容易被忽略的细节
数据下载看似没技术含量，但恰恰是翻车高发区。TCGA的样本名有TCGA-XX-XXXX-XX这种格式，不同组织来源的样本后缀不一样，如果你不做过滤，把实体瘤和血液样本混在一起跑WGCNA，模块结构会直接被干扰。GEO数据里则经常存在同一基因对应多个探针的情况，需要按最大表达量或者按方差最大做collapse处理后再进入分析。
另外还有一个细节：互作数据下载后先检查基因ID类型是Symbol还是Ensembl ID，STRING导出的TSV里有时候混着两种ID格式，如果你不做统一，后面在Cytoscape里会出现两个节点长得一模一样但ID不同的情况，网络结构直接失真。我一般会在进入分析前写一段脚本，统一把ID转成官方Symbol，并检查是否有重复项。
3. 核心开源工具实战：从基因列表到一张可分析的网络图
3.1 用STRING构建PPI网络：在线版也能做出科研级结果
构建PPI网络最快的方式是直接用STRING在线版。操作路径是：打开string-db.org，选择Multiple Proteins，把差异基因列表粘贴进去，物种选Homo sapiens，然后点Search。出来结果后先别急着用默认参数，点Settings把confidence score设为0.700，这是大多数文献里公认的高置信度阈值。如果互作数量太少，可以先降到0.400跑一个探索性网络看看整体结构，但最终用于模块分析和hub基因筛选的网络建议保持在0.700这一档。
在线版最容易被忽略的功能是Exports。在网络结果页右边选择"as simple tabular text output"，能下载到包含node1、node2、combined_score三列的TSV文件。这个文件是后续所有下游分析的基础。我习惯同时下载一份network layout文件，就是带坐标信息的格式，后面在Cytoscape里加载可以直接沿用STRING的布局，省去重新排版的时间。
3.2 Cytoscape导入与可视化：节点大小和颜色别只靠审美
拿到TSV之后，接下来的工作基本都在Cytoscape里完成。Cytoscape是目前基因网络可视化的事实标准，它本身是一个开源桌面软件，再加上MCODE、CytoHubba、clusterMaker这些插件，功能已经覆盖了从网络拓扑分析到模块识别的完整链路。
导入方式很简单：File → Import → Network from File，选中刚才的TSV文件。导入时Cytoscape会自动识别交互列，但你还是需要在预览界面确认source和target列指的是node1和node2，interaction列选combined_score，否则有时候它会把score识别为主题，导致边的属性错乱。
导入完成后，默认的展示方式所有节点一样大，没法看。这时候需要把网络拓扑参数计算出来：Tools → Analyze Network，勾选Treat the network as undirected。计算完成后节点表里会出现Degree等属性。接下来在Style面板里，把Node Size映射到Degree列，映射方式选Continuous Mapping，颜色也按Degree或者模块归属做映射。这个映射的本质含义是：在网络中连接越多的节点越大越醒目，你一眼就能看到哪些基因处于网络的核心位置。
3.3 MCODE模块识别：参数不是随便填的
网络图好看不等于有生物学意义，关键一步是做模块识别。MODE目前最常用的插件，全称是Molecular Complex Detection，算法逻辑是找网络里连接紧密的区域。插件从Cytoscape的Apps菜单里安装，装好后按Cluster → Open MCODE打开面板。
MCODE的几个核心参数值得认真对待。Degree Cutoff默认是2，意思是节点度至少为2才参与聚类；Node Score Cutoff默认0.2，控制子图打分门槛；K-Core默认2，要求最终模块里的每个节点至少有2个邻居；Max. Depth默认100，限制种子节点向外的搜索深度。对于癌症基因网络，我通常保留默认参数先跑一轮，看结果里模块数量和大小是否合理。如果模块太大（比如超过100个基因），说明网络太稠密，可以适度提高degree cutoff；如果模块太碎，说明网络太稀疏，可能需要降低confidence阈值重新构建网络。
MCODE跑完后会在网络下方显示模块列表，按Score排序，每个模块有对应的基因集合。这些模块就是后续下游分析的富集对象。
3.4 hub基因筛选：CytoHubba的12种算法怎么选
模块识别出来之后，下一个任务是从模块里挑核心基因。Cytoscape的CytoHubba插件提供了12种拓扑算法，从最简单的Degree、Betweenness、Closeness，到更复杂的MCC、DMNC、EPC都有。很多人在这里犯的错是只看Degree一种指标，把连接数最多的基因当成hub基因。
我的做法是至少看三个维度。先用Degree看连接度，再用Betweenness看网络中的桥接能力，然后用MCC看核心子网络的重要性。MCC这个算法的名字叫Maximal Clique Centrality，它对网络中的紧密小团体特别敏感，适合找模块内部真正的核心节点。实际操作时，CytoHubba面板里同时勾选这几种算法，分别输出Top 10基因，然后取交集。交集基因的数量通常在3到6个，这些才是真正经得起多维度检验的hub基因。
3.5 没有现成互作数据时：用WGCNA做共表达网络
有时候你面对的情况是：手上只有表达矩阵，没有对应的互作数据库覆盖，或者研究的是稀有肿瘤，已知互作信息极少。这时候还有一条独立路径——共表达网络。WGCNA是这个领域最经典的开源R包，全称Weighted Correlation Network Analysis，核心思想是用加权相关系数矩阵描述基因之间的共表达关系，然后划分出共表达模块。
WGCNA的关键参数是软阈值β。它不是一个随意定的值，WGCNA会通过pickSoftThreshold函数计算一系列候选幂次，然后选择能让网络达到无标度拓扑特性的那个值，通常要求scale-free topology fit index R²达到0.85以上。我实际跑下来，大部分癌症转录组数据取β=6到β=10之间就能达到标准。
得到模块之后，每个模块的特征基因（module eigengene）可以和临床性状做相关分析，比如肿瘤分期、生存状态、分子亚型，从而找到"和表型最相关的模块"。这一思路和PPI网络分析形成互补：PPI告诉你蛋白层面谁和谁能结合，共表达告诉你转录层面谁和谁协同变化。
4. 容易被低估的坑：ID转换、阈值设定、hub基因筛选的实操经验
4.1 基因ID不一致是网络分析最隐蔽的翻车点
我在实际项目中踩过最深的坑，就是基因ID的类型不统一。STRING默认显示的是官方Symbol，但TCGA下载的表达矩阵里，基因标识符可能是Ensembl Gene ID，也可能是Entrez ID，而富集分析时很多R包又只认Entrez ID。如果你在构建网络时用的是Symbol，富集时直接拿同一批Symbol去跑clusterProfiler，看着好像没问题，实际会有相当比例的基因匹配不上，尤其是一些已经有别名变更的基因，比如"KIT"在不同版本里可能有不同写法。
解决方式很简单：进R用biomaRt或者AnnotationDbi做一次批量映射。我一般用clusterProfiler自带的bitr函数，比如把Symbol转成Entrez ID：
library(clusterProfiler) library(org.Hs.eg.db) symbols <- c("TP53", "EGFR", "VEGFA", "MYC") entrez <- bitr(symbols, fromType = "SYMBOL", toType = "ENTREZID", OrgDb = org.Hs.eg.db)跑完可以检查一下na比例。如果超过10%找不到对应关系，优先怀疑数据源注释版本不一致，不要直接跳过。
4.2 STRING置信度阈值不是越高越好
很多人以为confidence设得越高网络越可靠，但实际并非如此。0.700能让网络保留实验验证级别的互作，但如果你的关键词集里包含大量研究不够充分的基因，高阈值会把它们全部孤立出来，这些孤立节点既进不了模块，也做不了功能富集，等于白分析了。
比较合理的做法是分两步走。先用0.400构建一个宽松网络，观察基因集整体是否形成大的连通区域，了解大概的网络规模；然后再用0.700重新构建，用于正式的模块和hub筛选。如果两个阈值下的核心模块高度一致，那这个发现就非常稳健，可以在文章里直接呈现；如果差异巨大，说明网络基础不够牢，先回数据源头检查差异基因的筛选流程更稳妥。
4.3 hub基因只是统计分析结果，不等于功能结论
这是我在评审和搭档合作中最常强调的一点。CytoHubba选出来的hub基因，本质上是网络拓扑分析的结果，它在计算层面是核心节点，但不等于它一定是这个癌症的关键驱动基因。hub基因是否真实参与肿瘤发生发展，最终需要表达验证、生存分析和功能实验来作答。
所以模块分析完成后，别急着写结论。我通常会把hub基因拎出来，在TCGA里看肿瘤和正常组织之间的表达差异，再跑一次Kaplan-Meier生存分析，看看高表达组的生存是否显著不同。至少要让网络预测的结果和独立临床数据对上，再进入湿实验环节。这一步虽然只是数据层面的验证，但能筛掉大量虚假的阳性信号。
5. 从"好看的图"到"站得住的结论"：网络结果的验证与解读
5.1 模块或子网络的功能富集：用clusterProfiler跑GO和KEGG
模块识别完成，hub基因筛选结束，接下来一定要做功能富集。否则你只是给审稿人看了一张密密麻麻的网络图，他们完全无法判断这个模块到底在做什么。富集分析我固定用R的clusterProfiler包，它支持GO、KEGG、Reactome等多种数据库，输出可以直接取数据框做可视化。
以MCODE跑出的一个包含56个基因的模块为例，先用bitr把Symbol转成Entrez ID，然后：
library(clusterProfiler) library(org.Hs.eg.db) ego <- enrichGO(gene = module_entrez$ENTREZID, OrgDb = org.Hs.eg.db, ont = "BP", pAdjustMethod = "BH", pvalueCutoff = 0.05, qvalueCutoff = 0.2)跑完之后看前几条显著富集的生物学过程。如果一个模块富集到"DNA replication"和"cell cycle"，那这个模块的生物学身份就清楚了，后续和你关注的两表型关联顺理成章。KEGG富集同理，enrichKEGG需要联网读取最新通路注释。
5.2 模块与临床表型的关联：从网络分析到临床价值
构建共表达网络时，WGCNA天然支持模块-表型关联分析，这是它相比单纯PPI网络的突出优势。每个模块计算出一个module eigengene，相当于这个模块在所有样本中的综合表达水平，然后与肿瘤分期、TNM分期、生存状态等临床变量做相关。
PPI网络分析也可以做类似的关联。把识别出的模块基因拿回来，在表达矩阵中提取这些基因的表达值，用平均表达量或者第一主成分作为模块分数，同样可以跟临床表型做关联。形式上不一定要正规的统计模型，简单的分组比较加箱线图就能看出趋势。需要提醒的是，多重比较时要关注P值校正问题，我在实际项目里会直接用BH方法，避免单个模块P值显著但全模块整体无效的情况。
5.3 外部验证：最低成本也能做得像模像样
网络分析最容易被质疑的一点是：结果是不是只在你这套数据里成立？所以无论如何，外部验证不可省略。最经济的方式是找一个独立的GEO数据集，相同癌种、相近样本量，下载表达矩阵后提取你hub基因的对应表达值，重复一次分组比较或生存分析。
如果外部队列的结果方向和原数据一致，那么这个网络分析的结论说服力就上去了。如果方向不一致，先不用急着否定自己。检查几件事：外部数据集的组织类型是否和TCGA一致、样本量和分组是否均衡、基因注释版本是否匹配。很多时候方向不一致不是因为网络分析错了，而是因为临床队列的异质性导致差异不显著。
网络分析这条路我前后跑了很多轮，最大的体会是：不要贪快，每一个步骤的输入和参数都要能说出为什么。开源工具的好处在于每个环节都透明、可回溯，你做的每一次阈值调整都有记录，这对发表和后续合作来说是最宝贵的资产。最后再分享一个小习惯：每次分析完，把基因列表、网络文件、参数设置和代码脚本打包归档，哪怕半年后再翻出来，也能立刻复现当时的完整结果。
