MuSiC软件 | 单细胞和常规转录组联合谱写优美乐章

2023-07-04 10:44:00, 欧易生物上海欧易生物医学科技有限公司

常规转录组（Bulk RNA-seq）是目前最为常用的高通量测序方法，它能够测量每个基因在组织中的平均表达水平，然而在许多研究中，了解组织内部的细胞类型组成及其比例非常重要。在此之前，小欧带大家认识过几款推断RNA-seq数据中细胞类型占比的软件，如cibersortx、TIMER。这里，我们介绍一款新的能够根据单细胞转录组数据推测Bulk RNA-seq细胞类型组成的反卷积方法软件——MuSiC(Multi-subject Single-cell Deconvolution)，那么该方法有哪些特点，与其他方法又有什么区别呢？跟着小欧一起往下看看吧！

MuSiC(Multi-subject Single-cell Deconvolution）软件能够借助单细胞转录组测序（scRNA-seq）中细胞类型特异性基因的表达谱，通过反卷积算法来表征复杂组织中bulk RNA-seq的细胞类型组成。

MuSiC分析流程框架如下图1所示。首先，MuSiC软件能够接收含有多种细胞类型的多样本表达数据，并构建分层聚类树。通过该聚类树，我们可以观察细胞类型之间的相似性。然后，MuSiC会确定一组一致性基因，并计算这些基因在每种细胞类型中的跨个体平均值（由红到蓝）和跨个体方差（由黑到白）。MuSiC会为个体间变化较小的基因增加权重，为个体间变化差异大的基因减少权重，这里引入MuSiC中的一个关键概念，就是“Marker基因一致性”(marker gene consistency)。当使用scRNA-seq数据作为细胞类型反卷积的参考时，需要考虑两个部分的表达一致性：跨个体和跨细胞，前者是为了避免个体偏差，后者是为了防止scRNA-seq中细胞捕获的偏差。因此在反卷积时，将给同种细胞类型间(跨细胞)和不同个体间（跨个体）表达一致的基因赋予较高的权重。最后，获得bulk RNA-seq样品中细胞类型的比例。

图1. MuSiC分析流程框架

在组织中通常会包含非常相近的细胞类型，这些细胞类型之间的基因表达很接近，使得难以在bulk数据中解析它们的相对比例。为了处理该问题，MuSiC采用聚类树引导方法，递归放大密切相关的细胞类型。简而言之，软件首先将相似的细胞类型分组到同一聚类中并估计聚类比例，然后在每个聚类中递归重复此过程（图1）。在每个递归阶段，只使用具有低簇内方差的基因，也就是在同种细胞类型内表达变化较小的基因，即跨细胞一致基因。

先前已有一些能够对bulk RNA数据的细胞类型比例进行反卷积计算的方法，比如CIBERSORT，它基于支持向量回归，允许使用scRNA-seq基因表达作为参考；TIMER是为癌症数据开发的，专注于免疫细胞浸润的定量。这些方法依赖于预先选择的细胞类型特异性标记基因，因此对显著性阈值的选择很敏感。更重要的是，这些方法没有考虑细胞类型特异性基因表达的跨个体异质性以及单细胞基因表达的细胞内随机性，而这正是MuSiC的特点。

作者将MuSiC应用于人胰岛组织中，结果发现MuSiC优于现有方法（图2），特别是对于具有密切相似细胞类型的组织。

图2. 不同方法预测人胰岛细胞类型组成结果比较

1）软件安装

# install devtools if necessaryinstall.packages(''devtools'')
# install the MuSiC packagedevtools::install_github(''xuranw/MuSiC'')
# loadlibrary(MuSiC)

2）数据准备

MuSiC需要以下两种类型的输入数据:

1.Bulk RNA测序的基因表达量；

2.scRNA-seq的基因表达量和细胞类型注释信息。

*Bulk RNA和scRNA-seq的表达量均使用counts数据。

GSE50244.bulk.eset = readRDS(''https://xuranw.github.io/MuSiC/data/GSE50244bulkeset.rds'')GSE50244.bulk.eset#ExpressionSet (storageMode: lockedEnvironment)#assayData: 32581 features, 89 samples #  element names: exprs #protocolData: none#phenoData#  sampleNames: Sub1 Sub2 ... Sub89 (89 total)#  varLabels: sampleID SubjectName ... tissue (7 total)#  varMetadata: labelDescription#featureData: none#experimentData: use ''experimentData(object)''#Annotation:  
bulk.mtx = exprs(GSE50244.bulk.eset)# Download EMTAB single cell dataset from GithubEMTAB.sce = readRDS(''https://xuranw.github.io/MuSiC/data/EMTABsce_healthy.rds'')EMTAB.sce#class: SingleCellExperiment #dim: 25453 1097 #metadata(0):#assays(1): counts#rownames(25453): SGIP1 AZIN2 ... KIR2DL2 KIR2DS3#rowData names(1): gene.name#colnames(1097): AZ_A10 AZ_A11 ... HP1509101_P8 HP1509101_P9#colData names(4): sampleID SubjectName cellTypeID cellType#reducedDimNames(0):#mainExpName: NULL#altExpNames(0):

3）细胞类型占比计算

# Estimate cell type proportionsEst.prop.GSE50244 = music_prop(bulk.mtx = bulk.mtx, sc.sce = EMTAB.sce, clusters = ''cellType'',                               samples = ''sampleID'', select.ct = c(''alpha'', ''beta'', ''delta'', ''gamma'',                                                                   ''acinar'', ''ductal''), verbose = F)head(Est.prop.GSE50244$Est.prop.weighted)         # alpha       beta       delta        gamma      acinar    ductal# Sub1 0.2079168 0.17460297 0.002203054 0.0015887041 0.095976854 0.5177117# Sub2 0.6964173 0.11979283 0.003683287 0.0007471486 0.030780490 0.1485789# Sub3 0.1913485 0.35211175 0.032037770 0.0003117540 0.191436872 0.2327533# Sub4 0.2339358 0.46732822 0.009425930 0.0226075392 0.053852583 0.2128500# Sub5 0.2160384 0.01280304 0.002134228 0.0006780175 0.254075730 0.5142706# Sub6 0.2026675 0.63455360 0.003597401 0.0005941822 0.006650877 0.1519364

Bulk样本中细胞类型占比热图

图片说明：横坐标为细胞类型名称，纵坐标为样本名称，热图由蓝到红表示占比由低到高。

Bulk样本中细胞类型占比柱状图

图片说明：横坐标为样本名称，纵坐标为细胞类型占比。

MuSiC软件的运算速度很快，结果准确性也较高，感兴趣的老师快来试一试吧！

参考资料：

1. Wang X, Park J, Susztak K, et al. Bulk tissue cell type deconvolution with multi-subject single-cell expression reference[J]. Nature communications, 2019, 10(1): 380.

2. https://xuranw.github.io/MuSiC/index.html