1. 项目背景与核心价值
基因共表达网络分析是转录组研究中的重要工具，而hdWGCNA（high-dimensional Weighted Gene Co-expression Network Analysis）作为WGCNA的升级版，特别适合处理单细胞和空间转录组这类高维数据。我在处理10x Genomics单细胞数据和Visium空间数据时发现，传统WGCNA在细胞亚群分析中会出现模块识别模糊的问题，而hdWGCNA通过改进的邻接矩阵计算和模块检测算法，能更精准地捕捉细胞类型特异的共表达模式。
这个方法的突破性在于：
- 处理单细胞数据时，通过智能降维保留生物变异信号
- 整合空间坐标信息时，能同时分析基因共表达和空间共定位
- 输出的基因模块可直接关联到细胞类型注释和空间热点区域
2. 核心算法原理拆解
2.1 改进的邻接矩阵计算
传统WGCNA使用Pearson相关系数计算基因间相似性，而hdWGCNA采用三步优化：
- 稳健相关系数计算：对稀疏单细胞数据使用Sparse Cor算法
adjacency = hdWGCNA::SparseCor(exprs_mat, method="pearson")- 幂律变换：β值通过scale-free拓扑拟合自动选择
- 拓扑重叠矩阵(TOM)计算：引入kNN图约束降低噪声
2.2 分层模块检测
采用迭代式分层聚类：
- 第一层：基于TOM的动态树切割
- 第二层：对大型模块进行子模块划分
- 模块合并：使用eigengene相关性阈值0.8
关键参数：deepSplit参数建议设为2-4，空间数据需要更高分辨率时可设为4
3. 单细胞数据实操流程
3.1 数据预处理
library(Seurat) library(hdWGCNA) # 标准单细胞预处理流程 sc_data <- NormalizeData(sc_data) %>% FindVariableFeatures() %>% ScaleData() # hdWGCNA专用处理 sc_data <- SetupForWGCNA( sc_data, features = VariableFeatures(sc_data), metacell_size = 20 # 元细胞聚合参数 )3.2 网络构建与模块识别
# 分细胞类型构建网络 sc_data <- ConstructNetwork( sc_data, soft_power = 6, # 通过pickSoftThreshold确定 tom_name = "TOM" ) # 模块识别 sc_data <- FindModules( sc_data, min_module_size = 30, merge_thresh = 0.8 )3.3 结果可视化
# 模块特征基因热图 ModuleHeatmap( sc_data, plot_modules = TRUE, cell_subset = "cell_type1" ) # 模块-性状关联分析 ModuleTraitCorrelation( sc_data, traits = c("nFeature_RNA", "percent.mt") )4. 空间转录组整合分析
4.1 空间权重矩阵构建
# 创建空间邻接矩阵 spatial_graph <- FindSpatialNeighbors( visium_data, k = 6, # 近邻数 weight_type = "gaussian" # 距离衰减权重 ) # 整合空间信息的TOM计算 spatial_tom <- SpatialTOM( exprs = GetAssayData(visium_data), spatial_graph = spatial_graph, alpha = 0.5 # 表达相似性与空间相似性的平衡参数 )4.2 空间模块可视化
# 模块的空间分布 SpatialFeaturePlot( visium_data, features = GetModules(sc_data)$genes, alpha = c(0.1, 1) ) # 共表达网络的空间布局 PlotSpatialNetwork( visium_data, feature1 = "GeneA", feature2 = "GeneB", edge_threshold = 0.7 )5. 实战经验与避坑指南
5.1 参数优化策略
- 软阈值选择：单细胞数据通常需要较低soft power(4-6)，空间数据需要更高(6-8)
- 模块最小基因数：单细胞建议30-50，空间数据建议50-100
- 处理批次效应：建议先用Harmony整合后再运行hdWGCNA
5.2 常见报错解决
- 内存不足问题：
options(future.globals.maxSize = 8000 * 1024^2) # 增加内存限制- 模块数量过多：调整deepSplit参数从4降到2
- 基因过滤过严：检查SetupForWGCNA中的min_mean_quantile参数
5.3 高级技巧
- 跨数据集比较：使用modulePreservation函数
- 时间序列分析：对伪时间分箱后分别构建网络
- 多组学整合：用WGCNA模块特征基因关联ATAC-seq峰
6. 下游分析思路
6.1 功能富集分析
# 使用clusterProfiler进行GO分析 module_genes <- GetModules(sc_data)$module1 ego <- enrichGO( gene = module_genes, OrgDb = org.Hs.eg.db, ont = "BP" ) dotplot(ego, showCategory=10)6.2 驱动基因识别
# 基于模块内连接度和基因重要性评分 hub_genes <- IdentifyHubGenes( sc_data, module = "module1", n_hubs = 10 )6.3 细胞通讯预测
# 将模块特征基因作为信号通路输入CellChat cellchat <- createCellChat( sc_data@assays$RNA@data, meta = sc_data@meta.data ) cellchat@netP$pathways <- GetModules(sc_data)$modules
