Bulk RNA-seq反卷积¶
一句话概述¶
利用BisqueRNA、MuSiC和CIBERSORTx等反卷积算法，从bulk RNA-seq混合表达谱中估计各细胞类型的比例，是在缺乏单细胞数据的大队列研究中推断组织细胞组成的关键计算方法。
核心知识点表格¶
| 知识点 | 关键内容 | 常用工具 |
|---|---|---|
| 反卷积原理 | 线性混合模型: bulk = Σ(比例×细胞表达) | - |
| 参考基方法 | 需要单细胞参考的细胞类型特征 | CIBERSORTx, MuSiC, BisqueRNA |
| 无参考方法 | 不需要外部参考数据 | CDSeq, Linseed |
| 特征选择 | marker基因/特征矩阵构建 | CIBERSORTx createSignature |
| 统计框架 | ν-SVR/加权NNLS/约束回归 | CIBERSORTx(ν-SVR), MuSiC(加权NNLS), BisqueRNA(NNLS) |
| 验证方法 | 模拟数据/FACS验证/IHC比较 | - |
| 批次校正 | 参考与bulk数据的平台差异 | BisqueRNA(RNA-decomposition) |
各步骤详解¶
第一步：理解反卷积原理¶
白话解释： 一个组织的bulk RNA-seq数据就像把不同颜色的颜料混合在一起——最终看到的是混合色。反卷积就是从混合色（bulk表达）推断各种原色（细胞类型）的比例。数学上，假设混合信号是各细胞类型表达的线性加权和，权重就是我们要估计的细胞比例。
技术细节： - 线性模型: Y_bulk = Σ(p_k × S_k) + ε - Y_bulk: bulk表达向量（每个基因的表达量） - p_k: 第k种细胞类型的比例（要估计的） - S_k: 第k种细胞类型的特征表达谱（从单细胞获取） - ε: 噪声 - 约束条件: 所有比例之和=1，每个比例≥0 - 挑战: 细胞类型特征的稳定性、参考数据与目标数据的平台差异
第二步：CIBERSORTx分析¶
白话解释： CIBERSORTx是最知名的反卷积工具。首先从单细胞数据构建一个"特征矩阵"（每种细胞类型最有代表性的基因表达模式），然后用支持向量回归(ν-SVR)解出每个bulk样本中各细胞类型的比例。还能进行样品级别的基因表达校正（GEP imputation）。
代码示例：
# CIBERSORTx主要通过web界面使用
# https://cibersortx.stanford.edu/
# 步骤1: 准备单细胞参考数据（特征矩阵构建）
library(Seurat)
sc_data <- readRDS("scRNA_reference.rds")
# 导出CIBERSORTx格式的单细胞数据
# 格式: genes × cells, 第一行是细胞类型标签
sc_counts <- GetAssayData(sc_data, slot="counts")
cell_labels <- sc_data$cell_type
# 创建CIBERSORTx输入文件
# 方法A: 直接提交单细胞矩阵（在线创建signature）
sc_for_cibersort <- as.data.frame(as.matrix(sc_counts))
sc_for_cibersort <- rbind(cell_labels, sc_for_cibersort)
rownames(sc_for_cibersort)[1] <- "CellType"
write.table(sc_for_cibersort, "sc_reference_for_cibersort.txt",
sep="\t", quote=FALSE)
# 方法B: 自己构建signature matrix（pseudo-bulk）
library(Biobase)
# 计算每种细胞类型的平均表达
avg_expr <- AverageExpression(sc_data, group.by="cell_type")$RNA
write.table(avg_expr, "signature_matrix.txt", sep="\t", quote=FALSE)
# 步骤2: 准备bulk数据
bulk_counts <- read.csv("bulk_rnaseq_tpm.csv", row.names=1)
# CIBERSORTx需要TPM或Microarray标准化值
write.table(bulk_counts, "bulk_mixture.txt", sep="\t", quote=FALSE)
# 步骤3: 在CIBERSORTx网站上运行
# 上传signature matrix和mixture file
# 参数: permutations=100, disable quantile normalization (RNA-seq)
# 输出: 每个样本各细胞类型比例
# 步骤4: 读取结果
cibersort_results <- read.csv("CIBERSORTx_Results.csv", row.names=1)
# 列: 各细胞类型比例 + P-value + RMSE + Correlation
第三步：MuSiC反卷积¶
白话解释： MuSiC (Multi-Subject Single Cell deconvolution) 的特色是：1）利用单细胞数据中多个受试者的信息来估计基因的跨个体变异性；2）对变异小的基因（更稳定的marker）给予更高权重。这样能更robust地处理个体差异和批次效应。
代码示例：
# 安装MuSiC
# devtools::install_github("xuranw/MuSiC")
library(MuSiC)
library(Biobase)
library(SingleCellExperiment)
# 准备单细胞参考（MuSiC v1.0.0+使用SingleCellExperiment）
# 需要包含subject信息（MuSiC的核心优势）
# 从Seurat对象转换为SingleCellExperiment
sc_sce <- as.SingleCellExperiment(sc_data)
# 确保colData中包含细胞类型和受试者列
# sc_sce的colData应包含 cell_type 和 patient_id 列
# 准备bulk数据（提取为矩阵）
bulk_mtx <- as.matrix(bulk_counts) # genes × samples矩阵
# 运行MuSiC反卷积
music_result <- music_prop(
bulk.mtx = bulk_mtx,
sc.sce = sc_sce, # SingleCellExperiment对象
clusters = "cell_type", # colData中细胞类型列名
samples = "patient_id", # colData中受试者列名（关键）
select.ct = c("T_cell", "B_cell", "Macrophage",
"Epithelial", "Fibroblast", "Endothelial"),
verbose = TRUE)
# 提取估计比例
estimated_props <- music_result$Est.prop.weighted
# 结果: samples × cell_types 的比例矩阵
head(estimated_props)
# 验证: 每行和应为1
rowSums(estimated_props)
# MuSiC还提供了不使用权重的NNLS版本（作对比）
props_nnls <- music_result$Est.prop.allgene
# 可视化
library(ggplot2)
library(reshape2)
props_long <- melt(estimated_props)
colnames(props_long) <- c("Sample", "CellType", "Proportion")
ggplot(props_long, aes(x=Sample, y=Proportion, fill=CellType)) +
geom_bar(stat="identity") +
theme(axis.text.x = element_text(angle=90, hjust=1)) +
labs(title="MuSiC Cell Type Deconvolution")
第四步：BisqueRNA反卷积¶
白话解释： BisqueRNA的独特优势是提供两种模式：1）有参考模式——当有配对的bulk+scRNA-seq数据时，能更好地校正平台差异；2）无参考模式——利用已知marker基因进行反卷积。它特别适合处理不同平台/批次间的系统偏差。
代码示例：
# 安装BisqueRNA
# install.packages("BisqueRNA")
library(BisqueRNA)
library(Biobase)
# === 有参考模式(Reference-based) ===
# 准备数据
sc_eset <- ExpressionSet(
assayData = as.matrix(sc_counts),
phenoData = AnnotatedDataFrame(data.frame(
cellType = sc_data$cell_type,
SubjectName = sc_data$patient_id,
row.names = colnames(sc_data))))
bulk_eset <- ExpressionSet(
assayData = as.matrix(bulk_counts),
phenoData = AnnotatedDataFrame(bulk_metadata))
# 运行BisqueRNA (Reference-based)
bisque_result <- ReferenceBasedDecomposition(
bulk.eset = bulk_eset,
sc.eset = sc_eset,
cell.types = "cellType",
subject.names = "SubjectName",
use.overlap = TRUE) # 利用配对信息校正
# 提取结果
estimated_props_bisque <- bisque_result$bulk.props
# 行: 细胞类型, 列: bulk样本
# === 无参考模式(Marker-based) ===
# 需要marker基因的data.frame（列: gene, cluster, 可选weight）
markers <- data.frame(
gene = c("CD3D","CD3E","CD2","IL7R",
"CD19","MS4A1","CD79A",
"CD68","CD163","CSF1R",
"EPCAM","KRT18","KRT19",
"COL1A1","DCN","LUM"),
cluster = c(rep("T_cell",4), rep("B_cell",3),
rep("Macrophage",3), rep("Epithelial",3),
rep("Fibroblast",3)))
# 运行marker-based反卷积
bisque_marker <- MarkerBasedDecomposition(
bulk.eset = bulk_eset,
markers = markers,
weighted = TRUE)
bisque_marker_props <- bisque_marker$bulk.props
第五步：方法比较与验证¶
白话解释： 不同反卷积方法可能给出不同结果，需要：1）在模拟数据上评估准确性；2）用多种方法取一致性结果；3）与实验验证数据（FACS/IHC）比较。
代码示例：
# 1. 模拟数据验证
# 用单细胞数据人工混合生成pseudo-bulk
create_pseudobulk <- function(sc_eset, cell_types, n_samples=50) {
proportions_matrix <- matrix(0, n_samples, length(unique(cell_types)))
colnames(proportions_matrix) <- unique(cell_types)
for (i in 1:n_samples) {
# 随机生成比例
props <- runif(ncol(proportions_matrix))
props <- props / sum(props)
proportions_matrix[i, ] <- props
}
# 按比例混合
expr_data <- exprs(sc_eset)
pseudo_bulk <- matrix(0, nrow(expr_data), n_samples)
for (i in 1:n_samples) {
for (ct in colnames(proportions_matrix)) {
cells_of_type <- which(cell_types == ct)
if (length(cells_of_type) > 0) {
sampled <- sample(cells_of_type,
size = round(proportions_matrix[i, ct] * 100),
replace = TRUE)
pseudo_bulk[, i] <- pseudo_bulk[, i] + rowSums(expr_data[, sampled])
}
}
}
list(bulk = pseudo_bulk, true_props = proportions_matrix)
}
sim <- create_pseudobulk(sc_eset, sc_data$cell_type)
# 对模拟数据运行反卷积，比较估计比例与真实比例
# 计算RMSE和相关系数
rmse <- sqrt(mean((estimated - sim$true_props)^2))
cor_per_ct <- diag(cor(estimated, sim$true_props))
# 2. 多方法比较
library(ggplot2)
comparison_df <- data.frame(
CIBERSORTx = cibersort_results[, "T_cell"],
MuSiC = estimated_props[, "T_cell"],
BisqueRNA = t(estimated_props_bisque)["T_cell", ],
Sample = rownames(cibersort_results))
ggplot(comparison_df, aes(x=CIBERSORTx, y=MuSiC)) +
geom_point() + geom_smooth(method="lm") +
geom_abline(slope=1, intercept=0, lty=2) +
annotate("text", x=0.1, y=0.4,
label=paste("r =", round(cor(comparison_df$CIBERSORTx,
comparison_df$MuSiC), 3)))
第六步：下游分析应用¶
代码示例：
# 1. 细胞比例与临床特征关联
library(survival)
library(survminer)
# 合并反卷积结果与临床数据
clinical_props <- merge(clinical_data, estimated_props, by="sample_id")
# 免疫细胞比例与生存分析
clinical_props$immune_high <- ifelse(
clinical_props$CD8_T > median(clinical_props$CD8_T), "High", "Low")
fit <- survfit(Surv(OS_time, OS_status) ~ immune_high, data=clinical_props)
ggsurvplot(fit, pval=TRUE, risk.table=TRUE)
# 2. 细胞比例与治疗反应
library(pROC)
roc_result <- roc(clinical_props$response, clinical_props$CD8_T)
cat("AUC for CD8 T cell proportion predicting response:", auc(roc_result), "\n")
# 3. 多变量Cox模型（包含细胞比例）
cox_model <- coxph(Surv(OS_time, OS_status) ~
CD8_T + Macrophage + Fibroblast + stage + age,
data = clinical_props)
summary(cox_model)
面试常问点¶
Q1: CIBERSORTx、MuSiC和BisqueRNA的主要区别？ A: CIBERSORTx用ν-SVR求解，最成熟但需要注册使用。MuSiC利用多受试者单细胞数据估计基因权重（变异小的基因权重高），更robust。BisqueRNA提供配对数据校正模式，处理平台差异更好。
Q2: 反卷积的线性假设合理吗？ A: 基本合理但有局限。大部分基因表达是细胞比例的线性函数（已被实验验证）。但存在非线性情况：细胞间互作改变表达（如T细胞激活取决于周围细胞）、微环境效应等。对于大多数应用，线性近似足够。
Q3: 反卷积结果能精确到什么程度？ A: 主要细胞类型（>5%比例）估计较准确（r>0.8），稀有细胞类型（<1%）估计不可靠。相对比例的排序比绝对数值更可信。建议用模拟数据评估每种细胞类型的估计精度。
Q4: 参考单细胞数据和目标bulk数据来自不同平台怎么办？ A: 1）BisqueRNA有跨平台校正模式；2）CIBERSORTx的batch correction功能；3）选择在两个平台上都稳定表达的基因作为特征；4）使用比例而非绝对值比较。最佳实践是使用与bulk同组织、同处理条件的单细胞参考。
易错点¶
参考数据细胞类型不全：如果目标组织含有参考中未包含的细胞类型，其信号会被错误分配到其他类型。
细胞类型定义粒度不当：粒度太粗（如只分T/B/Myeloid）丢失信息；粒度太细（如CD8 naive/memory/effector/exhausted）则标记重叠导致不可区分。
未做基因过滤：参考和bulk之间的基因交集可能不够。低表达基因的噪声会干扰估计。
TPM vs Counts混用：不同方法对输入数据标准化的要求不同。CIBERSORTx需要TPM，MuSiC可用counts。
单细胞参考的批次效应：参考数据本身有批次效应时，构建的特征矩阵会不准确。
补充知识¶
反卷积方法选择指南¶
| 场景 | 推荐方法 |
|---|---|
| 有配对单细胞+bulk | BisqueRNA (reference-based) |
| 有多个受试者的单细胞 | MuSiC |
| 无单细胞数据 | EPIC / quanTIseq (预建signature) |
| 免疫细胞专门分析 | CIBERSORTx (LM22) / TIMER |
| 大规模队列 | CIBERSORTx (batch correction) |
| 需要比例+表达校正 | CIBERSORTx (Hi-Res mode) |