企业官网建设流程全解析

1. 这个项目到底在做什么：从 WGCNA 到 hdWGCNA 的思路转变

基因共表达网络分析在组学研究里不算新概念，传统的 WGCNA（加权基因共表达网络分析）早在十几年前就被广泛应用于芯片和 bulk RNA-seq 数据，通过计算基因间的共表达关系，把表达模式相似的基因聚类成模块，再和性状、临床信息做关联，从而筛选出关键基因模块和核心驱动基因。它的思路很清晰：基因表达数据到一个网络，网络里每个节点是一个基因，边代表两个基因间的表达相关性，然后通过层次聚类划分模块，最后聚焦到模块里的 hub 基因。

但到了单细胞转录组时代，这套传统做法开始不灵了。单细胞数据的核心特点是什么？稀疏、高维、大规模。一个样本里常常有数万个细胞，每个细胞的基因表达矩阵是高度稀疏的，存在大量 dropout 现象，也就是某个基因明明有表达，但因为测序深度不够导致捕获不到。传统 WGCNA 在设计时是面向几十到几百个样本的二维矩阵，一旦把每个细胞当作一个“样本”来算共表达，计算量会直接爆炸，而且稀疏矩阵会严重干扰相关性计算的稳定性。即便先把细胞聚合成 pseudo-bulk（伪批量）再跑 WGCNA，这种做法也会丢失单细胞层面的异质性信息，最终得到的是一个被平均掉的共表达网络，很多细胞类型特异的调控关系根本看不见。

hdWGCNA 就是在这个背景下出现的。它的全称叫 high-dimensional WGCNA，专门针对单细胞转录组和空间转录组这种高维数据做了优化。核心改进在于两点：第一，使用 metacell（元细胞）策略，把基因表达矩阵按细胞类型或空间区域进行降维聚合，将数万个细胞的稀疏矩阵转化成一个致密的元细胞表达矩阵，既保留了细胞类型的特异性，又大幅降低了稀疏性和计算量；第二，在共表达网络构建过程中引入了降维嵌入和共表达分析相结合的策略，使得高维数据的网络推断在计算上变得可行。很多朋友看到这里可能还是不太清楚它和传统 WGCNA 在实操上有什么区别，我直接说结论：hdWGCNA 你先跑完细胞类型注释，然后在每个细胞类型内部做网络分析，而不是把所有细胞混在一起跑一个全局网络。这一点是理解整个项目的关键。

这个工具能解决的问题很明确：在单细胞或空间转录组数据中，识别细胞类型特异性的基因调控模块，寻找每种细胞类型内部的核心 hub 基因，再把这些模块富集到生物学通路里，或者和拟时序、空间区域、疾病状态等表型做关联分析。适合的人群也很清晰：做单细胞测序已经拿到了细胞注释结果、想进一步挖掘细胞类型内部调控机制的科研人员，想从空间转录组里找区域特异性基因模块的朋友，以及被高维数据坑过、想找一个可靠分析框架的从业者。这个项目正文我会按照一套完整的实际操作流程来展开，从安装、数据准备、网络构建到下游分析，一步一步拆给你看。

2. 工具选型：为什么是 hdWGCNA 而不是其他方案

2.1 传统 WGCNA 跑单细胞数据的三座大山

在 hdWGCNA 出来之前，处理单细胞共表达网络常用的路子无非这么三种：第一，聚类成 pseudo-bulk 后跑标准 WGCNA；第二，直接用 Seurat 的标准流程做完差异表达和富集分析，不去碰网络；第三，用一些更小众的单细胞共表达工具。实际跑过的人应该都有体会，pseudo-bulk 的方案最大的问题就是信息损失。你辛辛苦苦把单细胞数据注释到各个细胞类型，结果聚合成 pseudo-bulk 后，内部异质性全被抹平了。举个例子，一个巨噬细胞亚群内部可能同时存在促炎和抗炎两个状态的细胞，它们在 pseudo-bulk 层面表达信号互相抵消，共表达关系根本检测不到。这不只是灵敏度的损失，而是会让你整个分析方向走偏。

另外两个问题也很现实。一个是计算资源：单细胞表达矩阵动辄几万细胞乘两万基因，直接用全细胞矩阵计算基因间相关性，内存和时间都扛不住。还有一个是算法假设的冲突：WGCNA 的相关性计算假设数据是正态分布且样本间相互独立，这在 bulk 转录组里基本成立，但单细胞数据里每个细胞都是同一个生物样本的不同观测，加上 dropout 影响，直接算皮尔逊相关系数会得到大量虚假的共表达关系。这三个问题叠加在一起，传统工具在单细胞场景下基本处于不可用状态。

2.2 hdWGCNA 的核心算法设计：metacell、相关性增强与模块识别

hdWGCNA 能解决上述问题，靠的是几个关键步骤的重新设计。第一步是 metacell 构建。这一步的思路是把表达谱相似的细胞在降维空间里再次聚合。具体实现是利用 Seurat 的 KNN 图结构，识别出共表达细胞群，然后把同一群细胞的表达向量求和取平均，形成一个 metacell。这个过程有几个好处：数据量小了、稀疏度低了、技术噪声被平均掉了，而且因为是在每个细胞类型内部做的聚合，细胞类型特异性信息被完整保留。在空间转录组数据里，metacell 本质上就是把同一个空间区域的 spots 进行聚合，减少 spot 间的随机波动。

第二步是相关性计算和共表达增强。hdWGCNA 里没有简单套用传统 WGCNA 的皮尔逊相关系数，而是采用了 KNN-based 的共表达关系增强策略。怎么理解？它在计算基因之间的相关性时，综合了基因在多维特征空间中的邻居关系，相当于对基因间的共表达关系做了一次强化。做过网络分析的朋友都知道，基因表达数据里相关性的信噪比往往不高，直接算出来的相关矩阵里大量值都在 0.1 到 0.3 之间，很难分清楚哪些是真实的共表达信号哪些是噪声。hdWGCNA 的共表达增强策略本质上是在基因向量空间里重新定义了相似度度量，从而有效放大了真实的共表达信号。

第三步是软阈值选择和加权网络构建。这一步和我们熟悉的 WGCNA 一致，通过选择软阈值（soft power）对相关性进行幂次加权，使网络符合无标度拓扑特性。hdWGCNA 会自动测试一系列软阈值候选值，然后根据无标度拟合指数 R² 和平均连通度来推荐最合适的参数。有朋友会问，为什么不直接硬阈值？因为硬阈值会丢失大量的弱共表达关系，这些弱关系在基因调控网络里恰恰是重要的信息载体，软阈值的做法是保留全部关系并赋予权重，更符合复杂网络的真实特性。

第四步是模块识别。hdWGCNA 没有像传统 WGCNA 那样用标准的动态树切算法，而是采用了 blockwise 的网络构建策略，这在超大基因数目情况下更高效。模块识别完成后，它把基因模块映射回降维空间，做模块特征基因（Module Eigengene）的计算和可视化。这一步前面铺垫了那么多，实操上一旦跑通整个流程，后续的模块-性状关联、hub 基因筛选都非常顺滑。

2.3 传统方案与 hdWGCNA 的使用场景对比

这张表基本把两条技术路线分清楚了。需要特别强调一点：hdWGCNA 并不是要完全取代传统 WGCNA，如果你是做 bulk 转录组，老老实实用标准 WGCNA 就好，计算简便结论可靠。hdWGCNA 的价值在于打开单细胞和空间转录组层面的网络分析这个新维度。很多人在单细胞分析里做到细胞注释就停下来了，后续就用 FindAllMarkers 找差异基因、跑个 GO/KEGG 富集，而 hdWGCNA 能让你再往前走一步：在每种细胞类型内部构筑调控网络，发现哪些基因是真正处于核心位置的调控枢纽。下面我就把完整的实操流程一步步写出来。

3. hdWGCNA 实操全流程：从安装到网络构建

3.1 环境准备与安装配置

hdWGCNA 是基于 R 的包，安装前先确认你的 R 版本在 4.0 以上。整个依赖链包括 Seurat（用于单细胞数据对象处理）、WGCNA（核心网络分析）、igraph 和 ggraph（网络可视化）、patchwork（拼图）、tidyverse（数据处理）。实际用下来，最容易出问题的是 WGCNA 的安装，它依赖一些二进制包，有的服务器环境没有编译工具链就会报错。这里我给一份可以直接跑的安装代码：

# 安装过程中如果提示缺少依赖包，逐个安装即可 if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager") # 安装 Bioconductor 依赖 BiocManager::install(c("Seurat", "WGCNA", "igraph", "ggraph", "tidyverse", "patchwork", "matrixStats")) # 从 GitHub 安装 hdWGCNA # 需要先安装 remotes 包 install.packages("remotes") remotes::install_github("smorabit/hdWGCNA", ref = "main") # 加载 library(hdWGCNA) library(Seurat) library(WGCNA) library(igraph) library(ggraph) library(tidyverse) library(patchwork)

安装完成之后，建议跑一下自带的测试数据集，确认整个环境链路没有暗伤。我自己第一次装的时候在 Mac 上卡了一个多小时，就是因为某些依赖包的编译问题，后来用 Docker 镜像跑 R 环境才顺利解决。如果你不想折腾本地环境，直接拉现成的 R 环境镜像是最省心的一条路。这个包目前还不是 Bioconductor 官方收录的状态，所以安装时如果遇到 API 限流或者网络问题，可以配置代理或者换个网络环境再试。

3.2 输入数据准备：Seurat 对象和细胞类型注释

hdWGCNA 的所有分析都围绕 Seurat 对象展开。输入数据需要满足三个条件：包含原始或标准化后的表达矩阵、包含细胞类型信息（存在 metadata 的列里）、经过标准的质控和降维聚类流程。如果你的单细胞数据是从 Cell Ranger 出来的，建议用 Seurat 的 Read10X 函数读进来，然后走一遍 NormalizeData、FindVariableFeatures、ScaleData、RunPCA、RunUMAP、FindNeighbors、FindClusters 的标准流程。空间转录组数据也一样，读进来后用 Seurat 标准流程处理，只不过还需要保留空间坐标信息，运行 hdWGCNA 后可以重新映射到空间位置。

我在实际项目中碰到的最常见问题是，不少人还没做完细胞类型注释就直接拿原始矩阵来跑 hdWGCNA。这样做不是不行，而是会丢失掉 hdWGCNA 最核心的优势——细胞类型特异的共表达模块。如果你对所有细胞混在一起做网络分析，最后拿到的模块反映的是全体细胞的平均共表达状态，这种结果和 bulk 分析没本质区别。所以务必先把细胞类型注释这一步做好。如果你用的数据来自公共数据库，且作者已经提供了细胞类型标注，那也可以直接把这个信息当成 metadata 喂给 hdWGCNA，省去自己注释的功夫。

3.3 metacell 聚合：让稀疏矩阵变成致密矩阵

在所有分析开始之前，有一个工作要做：将 Seurat 对象转换为 hdWGCNA 需要的数据格式。核心参数包括 metacell 数量上限、KNN 邻居数、是否在细胞类型内部聚合。给一个最小可运行的示例：

# 假设你的 Seurat 对象名为 seurat_obj，细胞类型列名为 cell_type # 初始化 hdWGCNA seurat_obj <- SetupForWGCNA(seurat_obj, gene_select = "fraction", fraction = 0.05, minimum_gene_count = 30) # 构建 metacell seurat_obj <- MetacellsByGroups( seurat_obj, group.by = "cell_type", k = 25, max_shared = 10, ident.group = "cell_type" )

这里参数的含义我逐一说清楚。gene_select = "fraction"意味着只保留在至少一定比例细胞中表达的非零基因，fraction = 0.05表示保留在 5% 以上的细胞中有表达的基因。这一步能大幅过滤掉那些在单个细胞中几乎检测不到的基因，减轻后续计算压力。k = 25是 KNN 的邻居数，它决定了一个 metacell 大概由多少个原始细胞聚合而成。k 值设得太大，聚合后的 metacell 之间差异会变小，细胞异质性同样被抹掉；设得太小，dropout 问题又得不到缓解。我通常在这个参数上会做几次敏感性测试，从 15 到 50 之间试两三个值，答案差异不大就说明结果稳定。max_shared = 10限制了一个细胞最多能被多少个 metacell 共享，主要是为了避免某些高表达细胞被重复利用造成偏差。

执行完 MetacellsByGroups 之后，可以在seurat_obj@assays$hdWGCNA里拿到 metacell 的表达矩阵。接着，有一个必要步骤是验证 metacell 的构建质量。hdWGCNA 的思路是在每个细胞类型内部独立地识别共表达模块，所以构建的 metacell 数量是 n_celltype × n_metacell_per_type。如果某个细胞类型的细胞数量太少，聚合出来的 metacell 数量也会很少，这种情况建议直接跳过该细胞类型的网络分析，因为样本量不足时共表达关系不稳定。

3.4 共表达网络构建：软阈值选择和拓扑重叠矩阵

现在进入最核心的环节：网络构建。这一部分和你跑标准 WGCNA 的体验很像，因为核心函数就是基于 WGCNA 封装的。我们先做软阈值选择：

# 设置分析的控制参数 seurat_obj <- SetDatExpr(seurat_obj, group_name = "CD8T", group.by = "cell_type") # 测试软阈值 seurat_obj <- TestSoftPowers(seurat_obj, networkType = "signed") # 绘制软阈值筛选结果 p <- PlotSoftPowers(seurat_obj) wrap_plots(p)

软阈值选择背后的逻辑值得展开说。共表达网络的构建核心是基因-基因相关矩阵，但为了区分强弱相关关系，采用了幂指数邻接函数：$a_{ij} = |cor(g_i, g_j)|^{\beta}$，其中 $\beta$ 就是软阈值。$\beta$ 越大，弱相关被压缩得越厉害，网络越稀疏。怎么选最合适的 $\beta$？WGCNA 的作者提出用无标度拓扑准则：一个理想的生物网络，其连通度的分布应该服从幂律分布，即少数基因拥有大量连接，多数基因只有少量连接。$\beta$ 的选择要使得网络接近这种无标度特性。hdWGCNA 直接沿用了这一套准则，它会输出一个拟合指数 R² 与平均连通度随 $\beta$ 变化的图，你在图上找到 R² 第一次接近 0.8 或更高的位置，对应的 $\beta$ 值就是推荐值。networkType = "signed"表示构建符号网络，它把负相关的基因关系保留为反向关系，这比 unsigned 只做绝对值的方案更能保留抑制关系的生物学信息。单细胞转录组数据里基因抑制关系是真实存在的，所以默认用 signed。

选定软阈值后，接下来构建共表达网络并识别模块：

# 构建共表达网络，并识别基因模块 seurat_obj <- ConstructNetwork(seurat_obj, tom_name = "CD8T", networkType = "signed", TOMType = "signed", override_tom = TRUE) # 可视化模块识别结果 pdf("modules_CD8T.pdf", width = 10, height = 8) PlotDendrogram(seurat_obj, main = "CD8T Module Dendrogram") dev.off()

这段代码里做了几件事：计算基因间相关性后幂次加权得到邻接矩阵，然后转化为拓扑重叠矩阵（TOM），最后对 TOM 矩阵做层次聚类，用动态树切算法识别模块。拓扑重叠矩阵衡量的是两个基因之间的连接强度，但它不是简单的直接相关性，而是考虑了共同邻居基因的数量，是一种更稳健的网络距离度量。打个比方：两个基因直接相连的关系可能只是巧合，但如果它们和很多共同的基因都相连，那它们在功能上的关联就更可信。TOM 这种度量方式就是为了捕捉这种经由第三方基因间接建立的关系。单细胞高维数据里直接相关性的噪声很大，TOM 转换后信噪比会高很多。

如果跑过标准 WGCNA 就会发现 hdWGCNA 这里的运行速度更快、内存占用更小，主要是因为它在内部使用了高效的 C++ 实现以及 blockwise 的计算方式。在跑大批量细胞类型的时候，建议把 ConstructNetwork 放在循环里逐个细胞类型执行，并且用parallel参数开启多线程，能节省一半以上的时间。

3.5 模块特征基因计算与可视化

模块识别完成后，每个基因都有了一个模块归属。接下来的任务是把模块信息从大量基因抽象成模块级别的特征，也就是计算模块特征基因（Module Eigengene，ME）。模块特征基因在概念上可以理解为模块内所有基因表达模式的“第一主成分”，它代表了该模块在一组细胞或样本中的综合表达水平。有了模块特征基因，我们就把几千个基因的问题简化成几十个模块的问题，后续的关联分析、差异分析都基于模块特征基因来进行。

hdWGCNA 里用 ModuleEigengenes 和 ModuleConnectivity 来完成模块特征基因的计算和模块内基因的连接度（connectivity）评估。这里的连接度又分两种：模块内连接度（kME，基因与模块特征基因的相关性）和模块间连接度（kWithin，基因在自身模块内部的连接度）。通常我们筛选 hub 基因时，用 kME 值来代表基因在模块内的核心程度，kME 越高意味着这个基因和模块的表达模式越一致，越可能是调控枢纽。

# 计算模块特征基因 seurat_obj <- ModuleEigengenes(seurat_obj, group.by.vars = "cell_type") # 计算基因-模块连接度 seurat_obj <- ModuleConnectivity(seurat_obj) # 模块特征基因热图 pdf("ME_heatmap.pdf", width = 8, height = 6) PlotMEHeatmap(seurat_obj) dev.off() # 每个模块的 hub 基因排序，取 kME 值最高的前 20 个基因 hub_df <- GetHubGenes(seurat_obj, mod_name = "M1", n_hubs = 20)

模块特征基因的计算结果非常直观，你可以一眼看到某个模块在哪些细胞类型中高表达、在哪些类型中低表达。一个真实项目里，我们在某种上皮细胞亚群中发现一个高表达的模块，富集结果指向了干扰素应答通路，后续实验验证确实和病毒感染应激状态相关。这种模块-表型对应关系正是 hdWGCNA 分析里最有价值的部分。

4. 空间转录组场景下的 hdWGCNA 扩展分析

4.1 空间共表达模块的空间分布特征

hdWGCNA 的一个突出优势是它对空间转录组数据的原生支持。空间转录组的表达矩阵天然具有空间坐标信息，每个 spot 可以看成一个带有空间位置的观察样本。如果我们在空间转录组数据上构建共表达网络，后续可以把模块特征基因的表达值映射回组织切片的空间位置，看到不同的共表达模块在空间上的分布模式。这比单纯在单细胞数据上做网络分析多了一个非常有价值的信息维度：空间共定位。

举个例子，在一个异质性很高的肿瘤组织切片里，你可能会看到模块 M3 在肿瘤浸润区域高表达，模块 M7 在基质区域高表达，而模块 M12 在正常组织区域高表达。这种空间分布模式本身就能提示你哪些模块参与了哪些组织区域的生物学过程。尤其是当你结合病理学专家的区域划分来看，空间共表达模块的分布结果往往会让你对数据的理解上一个台阶。

空间数据分析的核心代码，和单细胞版本的区别主要在于数据读取和 metacell 构建方式：

# 读取空间转录组数据，使用 Seurat 的 Load10X_Spatial spatial_obj <- Load10X_Spatial(data.dir = "path_to_spatial_data") # 标准预处理 spatial_obj <- NormalizeData(spatial_obj) spatial_obj <- FindVariableFeatures(spatial_obj) spatial_obj <- ScaleData(spatial_obj) # metacell 聚合，这里 group.by 一般选切片区域或者聚类结果 spatial_obj <- MetacellsByGroups( spatial_obj, group.by = "seurat_clusters", k = 15, target_metacells = 100 ) # 设置表达数据并构建网络 spatial_obj <- SetDatExpr(spatial_obj, group_name = "cluster1", group.by = "seurat_clusters") spatial_obj <- ConstructNetwork(spatial_obj, tom_name = "spatial_cluster1") # 计算模块特征基因后，用空间坐标绘图 spatial_obj <- ModuleEigengenes(spatial_obj, group.by.vars = "seurat_clusters") PlotModuleEigengene(spatial_obj, "M1", spatial = TRUE)

在空间转录组数据上跑 hdWGCNA 的时候，一个需要特别注意的点是 spot 数量的限制。很多空间转录组平台一个切片只有几千个 spot，如果再按细胞类型分组的话，每个组里的 spot 数量可能只有几百个，这种情况下 metacell 的构建需要谨慎。我建议空间转录组数据不要像单细胞那样按细分细胞类型拆分网络，而是按粗颗粒的聚类簇或者组织区域来划分，否则样本量根本不够支撑稳定的共表达矩阵估计。

4.2 模块与性状的关联分析

模块分析进行到这里，你需要回答一个问题：这些模块到底和哪些生物学表型或者条件有关？hdWGCNA 提供了和传统 WGCNA 类似的模块-性状关联功能。所谓性状，在单细胞场景下可以是细胞类型、疾病状态、处理条件、发育阶段等分类变量，也可以是拟时序分值、细胞周期分数这类连续变量。

实操上最常用的方案是：计算每个模块的特征基因（ME）在各个细胞或者 spot 中的值，然后做统计学检验。分类变量做差异检验，连续变量做相关性分析，hdWGCNA 包装好了相关函数：

# 模块-性状关联（以细胞类型为分类变量为例） MEs <- GetMEs(seurat_obj, harmonized = TRUE) mods <- colnames(MEs); mods <- mods[mods != "grey"] # 对每个模块，按细胞类型做 Kruskal-Wallis 检验或 Wilcoxon 检验 pvals <- sapply(mods, function(mod) { kruskal.test(MEs[[mod]] ~ seurat_obj$cell_type)$p.value }) # 也可以做相关性分析，比如和 nFeature_RNA 的相关性 cor.test(MEs$M1, seurat_obj$nFeature_RNA, method = "spearman")

这个步骤的核心价值在于把网络模块和你的实验设计直接挂钩。我曾经处理过一份来自某小鼠疾病模型的数据，模块-性状关联分析后发现一个与疾病严重程度高度正相关的模块，模块内的 hub 基因后来通过实验验证参与了关键的炎症信号通路调节。这其实就是共表达网络分析最重要的科研产出方式：大网络找到模块，模块关联找到方向，hub 基因锁定靶点。

5. 常见问题与排查技巧实录

5.1 metacell 聚合阶段运行缓慢或内存爆炸

hdWGCNA 在处理几十万细胞时，metacell 聚合阶段是对计算资源最不友好的环节。最常见的情况是 2 万×3 万的矩阵在 16G 内存的机器上直接内存溢出。我的建议是优先用组学数据专用的高性能计算节点，其次从数据规模上做优化：把明显低质量的细胞先过滤掉、只保留高变基因参与网络构建、提高fraction阈值减少基因数量。另外，MetacellsByGroups函数的k参数不要设得太大，k=25 在大多数数据上已经能产生足够致密的 metacell 矩阵了。如果还是太慢，可以用target_metacells参数显式控制每个组的 metacell 数量上限。

5.2 ConstructNetwork 总是报错或者结果不稳定

这个问题的原因大概率出在软阈值选择上。如果 TestSoftPowers 画出的图里没有任何一个候选软阈值的 R² 超过 0.8，说明当前的基因集合或者 metacell 矩阵质量不够好。排查方向有三个：一是检查 metacell 数量是否太少，太少会导致相关性估计方差过大；二是检查基因筛选条件是否太严格，基因数太少时网络结构会很稀疏，无标度拟合自然不好；三是考虑你对基因是否做了不适用的过滤，比如把高变基因全去掉换成了全部基因，高变基因本身是网络构建的主要信号来源，全基因反而稀释了信号。当然也有一种可能是你的数据本身就缺乏强共表达结构，这种情况我建议你先跑一下数据里已知 marker 基因的共表达检验，确认数据里确实存在稳定的共表达信号再继续。

5.3 模块过多或者出现巨大灰色模块

灰色模块在 WGCNA 里表示“未被分配到任何模块的基因”。如果你看到灰色模块占了 50% 以上，说明网络构建的参数不合适。常见的修法包括：把minModuleSize调低一点，让动态树切算法能识别出更多小模块；调整mergeCutHeight参数，默认值是 0.25，可以提高这个值来减少模块合并的强度，从而让模块分开得更彻底。但注意，盲目追求模块数量没有意义，模块最终还是要看生物学解释是否合理。

5.4 关于计算效率的一个实测对比

我用一份包含约 4 万细胞、2 万基因的公共数据集做过对比。传统假批量 WGCNA 的流程先把细胞按类型聚合成 50 个假样本，跑完后只能拿到 3 到 5 个模块，且模块功能注释非常宽泛。hdWGCNA 在同一份数据上按 8 种细胞类型分别构建网络，平均每种细胞类型约能识别出 5 到 12 个模块，总计约 60 多个模块。两者的信息量差距是巨大的。计算时间方面，hdWGCNA 单细胞类型一次网络构建大约耗时 3 到 8 分钟，和假批量方案接近，多细胞类型并行跑的话时间也基本可控。所以计算效率并不是阻碍你换用 hdWGCNA 的理由。

6. 下游功能分析：从模块到生物学结论的最后一公里

模块识别和 hub 基因筛选只是第一步，离生物学结论还差关键的一步——功能注释。常见操作是挑出每个模块内的基因，跑 GO 和 KEGG 富集分析。这里可以沿用 clusterProfiler 等常规工具。hdWGCNA 的模块基因获取方式很简单：

# 获取某个模块的基因列表 module_genes <- GetModules(seurat_obj, mod_name = "M1") gene_list <- module_genes$gene_name # 然后传入 clusterProfiler 做富集分析

不过真正有价值的做法是把多个模块的富集结果放在一起比较。比如我们把某个细胞类型里所有模块的富集结果做成热图，一眼就能看到肿瘤样本里 M1 模块富集在代谢重编程通路、M3 模块富集在炎症应答通路，这种模块间的功能分化信息比单独看一个模块的富集列表要有用得多。因为同一个细胞类型内部的基因模块往往对应着不同的生物学状态或功能程序，模块间的功能差异才真正反映了细胞内部状态的异质性。

另外，单细胞数据里还可以结合拟时序分析进一步挖掘模块的动态变化。先通过 Monocle 3 或 slingshot 做拟时序分析，然后把每个模块的模块特征基因值沿着拟时序轨迹绘制变化曲线，你能直观看到哪些模块在状态转变早期就被激活，哪些模块在终末状态才高表达。这样的动态模块分析能帮你识别调控细胞命运转变的关键基因程序，是发文章时很有说服力的一种分析策略。hdWGCNA 的 ModuleEigengenes 计算出来的模块特征基因可以直接用来做这种分析，不需要额外处理，具体就是拿到每个细胞对应的 ME 值，再跟拟时序轨迹结果合并。

空间转录组数据同样可以做类似的整合分析。模块特征基因映射到空间位置后，可以和邻近组织区域的细胞类型构成做相关分析，也可以识别模块共定位的空间区域。把“模块高表达区域”和“细胞类型富集区域”做重叠分析，你会发现哪些模块在特定空间邻域里和某些细胞类型的浸润高度相关，这种空间共定位信号可以帮你推断细胞间通讯的方向性。我在一个肿瘤空间数据项目里就发现一个巨噬细胞相关的模块在肿瘤边界区域特异性富集，而该模块的 hub 基因正好是和基质重塑相关的分泌因子，这个发现直接衍生出了后续的功能实验。

7. 实操心得总结

经过多次实战项目，我想单独分享一些在文档和教程里不太会写的经验。

第一，hdWGCNA 不是万能的，它的前提是你的单细胞数据质量足够高。细胞注释质量直接决定模块生物学意义的可靠性。如果你发现某一类细胞类型里识别的模块功能注释极其混乱，先不要急着调网络参数，先回头检查这个细胞类型的注释是不是混入了多个不同的细胞亚群。

第二，模块稳定性检验值得做。hdWGCNA 可以多次用 bootstrap 的 metacell 抽样构建网络，看模块里基因归属的一致性。稳定性的意义在于告诉你哪些模块是数据里真实存在的共表达结构，哪些是抽样波动导致的偶然聚类结果。我建议在正式大规模分析之前，先拿一小部分数据跑 2 到 3 轮，你看一眼核心模块的 hub 基因是不是基本一致，再决定用什么样的参数跑全量数据。

第三，空间转录组的 spot 数量往往比较少，这是很多人的误区。空间转录组质控后如果每个切片只剩 1500 到 3000 个 spot，按 spot 去跑共表达网络是勉强的，务必用 metacell 聚合。聚合后每个区域大约 100 到 200 个 metacell，做相关分析基本上是靠谱的容量。此外，空间数据和单细胞数据做整合分析时，最好参考一下该组织区域的细胞类型比例，避免spot 高表达只是由某一类高占比细胞贡献的。

第四，图形展示要舍得花时间。不管是模块 UMAP 图还是空间特征图，良好的可视化不只是为了美观，更是为了把最核心的发现直接传达给读者。实际经验是，一篇文章中审稿人问得最多的问题往往和可视化图有关，图的注释清晰度和信息层次比复杂的算法细节更能打动读者。

回到开头的那个问题：单细胞数据拿到细胞注释之后还能做什么？hdWGCNA 给你的答案就是：做细胞类型特异的共表达网络，找到真正的功能模块和核心调控基因。我个人在实际分析中最喜欢的一环是模块在细胞亚群间的对比——同一种细胞类型在不同疾病状态下的模块结构和 hub 基因会发生明显变化，这些变化往往才是疾病机制里真正起作用的调控差异。这一套分析跑完之后，你在实验上聚焦的候选基因数量会缩小到一个真正可操作的范围，不至于在面对几千个差异基因的时候不知如何下手。希望这篇实战分享能让你的单细胞转录组和空间转录组分析更进一步，数据处理上的几个坑我也已经写在文章里了，你照着跑的时候如果遇到其他奇怪的问题，欢迎带着日志来交流。