🤩 BLEND | 新一代反卷积利器！~（多参考集估算更准更快）

生信漫卷

发布于 2025-11-13 18:30:31

960

文章被收录于专栏：R语言及实用科研软件R语言及实用科研软件

写在前面

随着单细胞转录组和bulk RNA-seq技术的快速发展，我们获得了前所未有的组织水平基因表达数据。😂

然而，这些数据往往代表的是多种细胞类型的混合信号, 每一个样本的转录组，实质上是多个细胞类型的“加权平均”。🫠

这就带来了一个核心问题：👇

我们如何准确地“反卷积”（deconvolution），即从混合信号中推断出每种细胞类型的比例？

传统的反卷积方法（如 CIBERSORTx、MuSiC、SCDC、BisqueRNA 等）通常依赖单一的参考单细胞转录组数据集。

然而，不同实验条件、组织来源甚至测序平台之间的系统性偏差（batch effect）常常会导致参考数据与目标样本的不匹配，从而降低反卷积精度。

BLEND包，可在无需标记基因筛选、无需数据标准化、无需参考质量评估的情况下，直接利用多个参考单细胞数据集，对bulk RNA-seq 计数数据进行高精度的细胞类型比例估计。

BLEND 的独特之处在于：

🔹 Individualized references — 为每个批量样本构建个性化参考；

🔹 Bag-of-words representation — 对表达谱进行表征，避免复杂的预处理步骤；

🔹 Probabilistic modeling — 通过概率建模框架实现稳健、灵活的反卷积。

用到的包

rm(list = ls())
if (!"BLEND" %in% rownames(installed.packages())) {
  devtools::install_github('Penghuihuang2000/BLEND')
}
library(BLEND)
library(tidyverse)

示例数据

load("BLEND_example.RData")
names(BLEND_example)

Bulk数据必须是counts, 或者做了线性转换也是可以的。😘

dim(BLEND_example$bulk)

BLEND_example$bulk[1:6,1:6]

参考数据集

参考数据集需要时一个list，其中每个元素代表来自多个参考数据集的细胞类型的特异性基因表达矩阵。🥸

需要注意的是，不同的细胞类型可能有不同数量的参考数据。🫢

names(BLEND_example$reference)

dim(BLEND_example$reference$Astrocytes)

colnames(BLEND_example$reference$Astrocytes)

colnames(BLEND_example$reference$Excitatory)

BLEND是一种基于RNA的模型。😅

因此，估计的分数是细胞类型的RNA分数。😂

如果需要估算细胞组分，则需要调整细胞大小。😀

细胞大小可量化不同细胞类型的相对RNA丰度，可以从scRNA-seq（sn）数据中估计或由其他来源提供。🧐

BLEND_example$cell_size

开始反卷积

这里提供三种方法：mixSQP、Gibbs和EM-MAP算法。😅

mixSQP算法是最新的算法，比EM-MAP算法和Gibbs更快。🥸

速度：mixSQP > EM-MAP >> Gibbs。🙊

time.mixSQP <- system.time(res.mixSQP <- BLEND(bulk = BLEND_example$bulk,
             phi = BLEND_example$reference,
             method = "mixSQP",
             ncore = 4))[3]
cat("Average running time per sample using one core: ",round((time.mixSQP*5)/30,1), "sec")

time.EMMAP <- system.time(res.EMMAP <- BLEND(bulk = BLEND_example$bulk,
             phi = BLEND_example$reference,
             method = "EMMAP",
             ncore = 4))[3]
cat("Average running time per sample using one core: ",round((time.EMMAP*5)/(30*60),1), "min")

time.GIBBS <- system.time(res.GIBBS <- BLEND(bulk = BLEND_example$bulk,
            phi = BLEND_example$reference,
            method = "GIBBS",
            ncore = 5))[3]
cat("Average running time per sample using one core: ",round((time.GIBBS*30)/(30*60),1), "min")

估算结果如下：👇

round(res.mixSQP$cellular_fraction[1:6,], 3)

round(res.mixSQP$reference_mix_prop$Excitatory[1:6,],3)

可视化

df <- as.data.frame(res.mixSQP$cellular_fraction)

df_long <- df %>%
  rownames_to_column("Sample") %>%
  pivot_longer(-Sample, names_to = "CellType", values_to = "Fraction")

ggplot(df_long, aes(x = Sample, y = Fraction, fill = CellType)) +
  geom_bar(stat = "identity", width = 0.8) +
  labs(y = "Cellular Fraction", x = "Sample") +
  theme_minimal(base_size = 13) +
  theme(axis.text.x = element_text(angle = 60, hjust = 1),
        legend.position = "right") +
  scale_fill_brewer(palette = "Set2")

library(pheatmap)

pheatmap(
  t(df), 
  cluster_rows = TRUE, 
  cluster_cols = TRUE,
  color = colorRampPalette(c("navy", "white", "firebrick3"))(50),
  main = "Cellular Fraction Heatmap"
)