主题
10 · bulk 转录 / WGCNA / GSEA / 组织定位
代码目录同 08：
/home/research/dr-amd-targets/analysis/。
对应脚本：01_bulk_analysis.R（DEG）、02_wgcna_gsea.R（共表达模块）、03_gsea_fgsea.R（单基因 GSEA）、31_gct_nsr_vs_rpe.R（组织定位）
一、bulk 差异表达（01_bulk_analysis.R）
数据
| 数据集 | 平台 | 对比 | 样本 | |
|---|---|---|---|---|
| AMD | GSE29801 | GPL4133 | 黄斑 RPE-脉络膜 AMD vs 正常 | 多（双色 log 比值） |
| DR | GSE60436 | GPL6884 | PDR 纤维血管膜 vs 正常视网膜 | 仅 9 例 |
全程本地解析 series matrix + GPL 注释，不依赖 GEOquery（离线可复现）。
处理链
r
maybe_log() # 99分位 > 50 判为线性尺度 → log2(x+1);已是 log 的不动
collapse_genes() # probe→基因:同基因取行均值最大的 probe;丢弃 NA/空/含 "///" 的多映射探针
design <- model.matrix(~ group) # 参照自动设为 normal/control(grepl 匹配)
# limma moderated-t + BH FDR为什么用 limma 不用逐基因 t 检验：脚本注释写明——微阵列小样本用经验贝叶斯收缩方差，比逐基因 t 更稳更有功效（2026-07-25 由 base t 升级）。
实际结果（*_DEG_candidates.csv）
DR（GSE60436）——16 个候选里 11 个显著：
| 基因 | logFC | adj.P |
|---|---|---|
| CFH | +0.920 | 1.3e−4 |
| SPRY2 | −0.763 | 7.7e−4 |
| TNFRSF10A | +0.869 | 2.7e−3 |
| MERTK | −0.618 | 0.0115 |
| TNXB | +0.382 | 0.0138 |
| TGFB1 | +0.465 | 0.0142 |
| CSF2 | +0.350 | 0.0204 |
| NFKB1 | +0.309 | 0.0257 |
| F13B | +0.371 | 0.0312 |
| CFI | +1.125 | 0.0359 |
| CASP10 | +0.284 | 0.0494 |
AMD（GSE29801）——全部不显著，最小 adj.P = 0.255（TNFRSF10A）；MERTK adj.P = 0.873（彻底阴性）。
★ DR 侧这批 p 值必须打折看
GSE60436 比较的是 PDR 纤维血管膜 vs 正常视网膜——两种不同组织，不是同组织的病/正常。 纤维血管膜富含成纤维/内皮/免疫细胞，正常视网膜富含神经元/胶质， 所以"16 个候选 11 个显著"主要反映组织构成差异，而非疾病效应。
这直接影响 PPT 里"MERTK ↓ FDR 0.012"这条证据的分量：它不能被当作"DR 使 MERTK 下调"的证据， 只能说"纤维血管膜中 MERTK 低于正常视网膜"。加上仅 9 例样本，这一层证据必须降级陈述。 AMD 侧则是干净的同组织对比，结果是全阴性——如实保留。
二、WGCNA 共表达模块（02_wgcna_gsea.R）
r
if (n < 12) 跳过 # ★ 样本 <12 直接不做 → DR(9例) 被跳过
expr <- expr[order(-var)[1:min(5000, nrow(expr))],] # 取方差前 5000 基因
sft <- pickSoftThreshold(datExpr, powerVector = 1:20)
power <- ifelse(is.na(sft$powerEstimate), 6L, sft$powerEstimate) # 估不出就退回 6
net <- blockwiseModules(datExpr, power = power, minModuleSize = 30,
mergeCutHeight = 0.25, maxBlockSize = 6000)只有 AMD 跑了（DR 9 例 < 12，脚本自动跳过并打印功效不足提示——这个自我保护是对的）。
实际结果（AMD_target_modules.csv）：
| 靶点 | 模块 |
|---|---|
| MERTK | blue |
| SPRY2 | blue |
| AGER | （空） |
| TNFRSF10A | （空） |
| APOE | （空） |
5 个靶点里 3 个没进任何模块（未入选方差前 5000，或被划为灰色未分配）。 所以 WGCNA 这一层只对 MERTK/SPRY2 有话可说，且两者同在 blue 模块。证据分量有限，不宜在正文强调。
三、单基因 GSEA（03_gsea_fgsea.R）
绕开 clusterProfiler，用 fgsea + org.Hs.eg.db 自建 GO-BP 基因集：
r
pathways <- pathways[lengths(pathways) >= 10 & lengths(pathways) <= 500] # 基因集大小过滤
cors <- apply(expr, 1, function(x) cor(x, expr[target,], use="pairwise.complete.obs"))
ranks <- sort(cors[is.finite(cors) & names(cors) != target], decreasing = TRUE)
res <- fgsea(pathways, ranks, minSize = 10, maxSize = 500, eps = 0)方法本质要讲清楚：这不是疾病 vs 正常的 GSEA，而是"全基因按与靶点表达的相关性排序"再跑富集—— 得到的是靶点的通路背景/共表达邻居，属探索性，不能解释为"疾病中该通路被激活"。
跑的组合：{AMD, DR} × {MERTK, AGER, SPRY2}，输出 results/gsea/<数据集>_<靶点>_fgsea.csv。 PPT 结果图 10 的措辞"通路背景（探索性）"是准确的，保持。
四、组织定位 NSR vs RPE（31_gct_nsr_vs_rpe.R）
EGA 视网膜数据，359 个逐样本 RNA-SeQC gene_reads.gct 合并成计数矩阵， 按文件名 MANGT_###.{NSR|RPE}.rnaseqc 解析组织与供体，limma-voom 做 RPE vs NSR。
注意：这是组织间对比（genotype-free），与疾病无关——用途是"选细胞系时该用 RPE 系还是视网膜神经系"。
实际结果（candidate_targets_RPE_vs_NSR.csv，logFC = RPE − NSR）：
| 基因 | NSR logCPM | RPE logCPM | logFC | adj.P |
|---|---|---|---|---|
| CFH | 3.67 | 6.13 | +2.45 | 4.2e−52 |
| TGFB1 | 2.54 | 4.74 | +2.17 | 5.5e−61 |
| TNXB | 2.90 | 3.72 | +0.83 | 5.9e−13 |
| CASP10 | 0.45 | 1.16 | +0.70 | 2.6e−15 |
| MERTK | 4.15 | 4.83 | +0.63 | 1.0e−18 |
| SPRY2 | 5.92 | 5.65 | −0.27 | 6.9e−5 |
| GAS6 | 5.78 | 5.22 | −0.54 | 9.0e−16 |
| AGER | 2.78 | 1.98 | −0.74 | 3.4e−21 |
| APOE | 9.66 | 8.77 | −0.90 | 2.6e−26 |
| TNFRSF10A | −1.21 | −1.19 | +0.07 | 0.600 |
★ 这张表里藏着一个对主线的硬约束
TNFRSF10A 在神经视网膜和 RPE 中的 logCPM 都是负值（≈ −1.2），即基本不表达。
它是 AMD 侧共定位最强的候选（PP 0.92），但在靶组织里测不到表达。 这意味着：① 它的因果作用可能通过循环蛋白而非视网膜局部表达实现； ② 无法用视网膜/RPE 细胞系做表达层验证； ③ 结合实验验证方案里查明的"小鼠、斑马鱼均无可用同源"， TNFRSF10A 的功能验证路径目前是全线受阻的，写作时必须正面交代而不是回避。
相比之下 CASP10 在 RPE 有表达且显著高于 NSR（+0.70, p=2.6e−15）——从组织表达角度它比 TNFRSF10A 更可做。
五、这四层的证据分量小结
| 分析 | 状态 | 可用性 |
|---|---|---|
| bulk DEG · AMD | 全阴性 | 如实报阴性 |
| bulk DEG · DR | 11/16 显著，但组织混杂 + n=9 | ⚠ 必须降级 |
| WGCNA | 只有 AMD 跑成，5 靶点仅 2 个进模块 | ⚠ 分量有限 |
| 单基因 GSEA | 探索性通路背景 | ⚠ 不可作机制结论 |
| NSR vs RPE 组织定位 | n=359，统计扎实 | ✅ 最硬的一层，且直接指导选细胞系 |
输出物
results/deg/ AMD_GSE29801_DEG_all.csv / _candidates.csv
DR_GSE60436_DEG_all.csv / _candidates.csv, _bulk_objects.rds
results/wgcna/ AMD_module_trait.csv, AMD_target_modules.csv (DR 未跑)
results/gsea/ <数据集>_<靶点>_fgsea.csv
results/gct_expr/ DE_RPE_vs_NSR.csv, candidate_targets_RPE_vs_NSR.csv, figures/上一页 → 09 · 细胞通讯与活化轨迹