基于跨细胞类型差异分析的复杂样本反卷积-凯发国际一触即发

期刊菜单

基于跨细胞类型差异分析的复杂样本反卷积
deconvolution of complex samples based on cross-cell type differential analysis

doi: , , html, , 科研立项经费支持
作者: 田中禾, 彭凌, 张伟伟^*：绍兴文理学院数理信息学院，浙江绍兴
关键词: ；；；；；；；；；

摘要: 复杂样本的基因组数据反映了多种细胞类型表达的平均水平。然而，细胞组成的差异会导致许多相关分析结果产生偏差。因此，准确估计细胞组成是分析复杂样本的第一步。目前有许多计算方法已经被开发出来估计复杂样本的细胞组成，但由于缺乏参考数据和先验信息，它们的应用大多有限。针对此问题，本文开发了一种基于最优特征选择的无参考反卷积算法，该算法通过整合一种细胞类型和其他细胞类型之间的跨细胞类型差异表达分析，以及两种细胞类型与其他细胞类型之间的跨细胞类型差异表达分析，迭代地搜索细胞类型特异性的最优特征，并进行细胞组成估计。模拟研究和两个真实数据集上的实证分析表明我们的算法具有更好的性能。

abstract: genomic data from complex samples measure average level of multiple cell types, and differences in cell compositions lead to biased results in many relevant analyses. therefore, accurately estimating cell compositions has been recognized as an important first step in analyzing complex samples. many computational methods have been developed to estimate cell compositions, but they mostly have limited applications due to a lack of reference or prior information. in this work, we develop a feature selection method for reference-free deconvolution, which iteratively searches for cell-type specific features by cross-cell type differential analysis between one cell type versus the other cell types, as well as between two cell types versus the other cell types, and performs composition estimation. extensive simulation studies and analyses of two real datasets demonstrate the favorable performance of our proposed method.

文章引用：田中禾, 彭凌, 张伟伟. 基于跨细胞类型差异分析的复杂样本反卷积[j]. 应用数学进展, 2024, 13(11): 4870-4875.

1. 引言

复杂样本的基因组数据测量的是其含有的多种细胞类型信号的加权平均值，而不是混合物中不同细胞类型的独立测量值，因此，包含差异分析在内的复杂样本的多种相关分析工作会由于细胞混合而导致分析结果出现偏差甚至是错误[1]。此外，有研究表明细胞组成与疾病进展和治疗反应密切相关[2]。例如，神经元细胞死亡是阿尔茨海默病(ad)的显著病理特征[3]，而免疫浸润在多种癌症类型中与肿瘤进展显著相关[4]。因此，准确估计细胞组成对分析复杂样本具有重要的意义。

传统上，许多实验方法，如荧光激活细胞分选(facs)、免疫组织化学(ihc)和单细胞rna测序(single cell rna-seq)能够提供细胞组成信息[5]。然而，这些方法由于通量有限、成本过高、操作繁琐等原因，在大规模临床应用中受限。作为替代方案，有许多被称为“反卷积”的计算方法被开发出来。目前，这些方法大致可分为两类：基于参考基矩阵的方法[6]-[9]和不依赖于参考基矩阵的方法[10]-[12]。基于参考基矩阵的方法需要将从纯化组织中获得的参考面板数据作为输入，并通过如约束线性回归或支持向量回归等方法来估计细胞类型的组成[6]。研究表明，基于参考基矩阵的方法通常比不依赖于参考基矩阵的方法得到的结果更准确[6] [7]。然而，目前的参考数据仅针对少数几种组织类型(如血液、乳腺和大脑)，且仅来自少数个体，因此基于参考基矩阵的方法应用范围受到极大限制。此外，当复杂样本与参考数据在临床条件和表型特征上存在显著差异时，基于参考基矩阵的方法可能会导致细胞组成估算不准确[13]。对于没有适当参考信息的组织，不依赖于参考基矩阵的方法提供了更好的凯发国际一触即发的解决方案[14]。该类方法基于某种类型的因子分析(如非负矩阵分解)，纯细胞类型的特征谱和细胞组成被同时进行估算[15]。然而，由于不依赖于参考基矩阵的方法中需要估计的参数众多，因此其面临着准确度较低的问题。因此，开发一种能够改进不依赖于参考基矩阵的反卷积新算法是一个值得研究的重要问题。

本文中，我们开发了一种基于跨细胞类型差异分析来筛选细胞类型特异性的特征并进行细胞组成估计。该算法通过整合一种细胞类型和其他细胞类型之间的跨细胞类型差异分析，以及两种细胞类型与其他细胞类型之间的跨细胞类型差异分析，迭代地搜索细胞类型特异性的最优特征，并进行细胞组成估计。广泛的模拟研究和两个真实数据集上的实证分析证明了我们方法的优越性。

2. 模型构建

给定一个 $m \times n$ 的矩阵y，其表示m个特征在n个复杂样本上高通量实验生成的数据。我们假设这些复杂样本是k种细胞类型的混合物，w为 $m \times k$ 的矩阵，表示m个细胞类型特异性特征在k种细胞类型上的表达谱，h为 $k \times n$ 的矩阵，表示细胞组成比例矩阵，h的各元素要求非负，并且每列元素之和为1。假设观测数据y为细胞类型特异性特征谱的线性组合，组合系数为细胞组成比例，即

$y = w h ε$

其中 $ε$ 是 $m \times n$ 的误差矩阵。将我们的算法命名为decd，我们算法的目标是仅利用观测数据y来估算细胞组成比例矩阵h。

figure 1. workflow of the proposed algorithm decd

图1. 算法decd的流程图

图1展示了decd的工作流程，该流程由两个主要步骤组成：步骤(a)从原始数据y中选取变异系数最大的前1000个特征(m₀)，在矩阵y中取m₀个特征，得到矩阵 $y_{m_{0}}$ ，对 $y_{m_{0}}$ 利用deconf进行无参考反卷积，得到比例矩阵h。步骤(b)利用步骤(a)估算的细胞比例矩阵h和原始数据矩阵y进行跨细胞类型差异分析，将一种细胞类型和其他细胞类型之间的跨细胞类型差异分析，以及两种细胞类型与其他细胞类型之间的跨细胞类型差异分析所得到的差异基因整合起来，将其作为新的m₀，并在新一轮迭代中用于步骤(a)，该算法迭代多次后停止。根据我们的经验，对于包含四种细胞类型和100个样本的基因表达数据，30次迭代就得到满意的结果。对于样本量较小(例如少于50个)或细胞类型较多(例如6种或更多)的研究，则需要更多的迭代次数。在我们的软件中，用户可以指定总迭代次数。在每次迭代中，decd会计算重构观测值 $\hat{y}$ 与真实观测值y之间的均方根误差(rmse)，并选择对应于最小rmse的迭代结果所对应的细胞组成比例矩阵h作为最终估计结果。

decd中跨细胞类型差异分析的统计模型如下：

设所有观测数据为 $y = {[y_{1}, y_{2}, \dots, y_{n}]}^{t}$ ， $μ_{k}$ 表示细胞类型k的表达谱。假设所有样本的细胞类型比例已知，并将样本s的细胞组成比例表示为 $θ_{s} = (θ_{s 1}, θ_{s 2}, \dots, θ_{s k})$ 。因此，观测数据y可描述为一个线性模型： $e [y] = w β$ ，其中

$w = [\begin{matrix} θ_{11} \\ θ_{21} \\ ⋮ \\ θ_{n 1} \end{matrix} \begin{matrix} θ_{12} \\ θ_{22} \\ ⋮ \\ θ_{n 2} \end{matrix} \begin{matrix} \dots \\ \dots \\ ⋮ \\ \dots \end{matrix} \begin{matrix} θ_{1 k} \\ θ_{2 k} \\ ⋮ \\ θ_{n k} \end{matrix}]$ , $β = [μ_{1}, μ_{2}, \dots, μ_{k}] .$

我们关注两个假设检验：

(1) 细胞类型k与其他细胞类型之间的差异表达分析： $h_{0} : μ_{k} - \frac{1}{k - 1} \sum_{i \neq k} μ_{i} = 0$ ；

(2) 检验细胞类型 $k, l$ 与其他细胞类型之间的差异表达分析： $h_{0} : (μ_{k} μ_{l}) - \sum_{i \neq k, l} μ_{i} = 0$ 。

第一个假设检验执行的是任一种细胞类型与其余细胞类型之间的差异分析，第二个假设检验执行的是任两个细胞类型与其余细胞类型之间的差异分析。decd分别从这两个假设检验中选择差异最显著的500个特征，并且确保这两部分特征之间的重叠为0，然后将这两部分特征结合在一起，在下一轮的反卷积中使用。

本文通过两个指标评估decd的性能：(1) 估计比例矩阵 $\hat{h}$ 与真实比例矩阵h在四种细胞类型上的平均皮尔逊相关系数mc；(2) 估计比例矩阵 $\hat{h}$ 与真实比例矩阵h在四种细胞类型上的平均绝对误差mae。更高的mc值和更低的mae值表明算法具有更优的性能。

3. 模拟分析及真实数据分析

模拟数据y基于两个随机生成的矩阵来构造：细胞类型特异性参考矩阵w和细胞组成比例矩阵h。比例矩阵h根据dirichlet分布进行模拟，y由四种细胞类型组成时比例矩阵h的参数设置为(0.968, 4.706, 0.496, 0.347)。参考矩阵w根据基因表达综合数据库(geo)中登录号为gse11058的免疫数据集来生成。该数据集包含四种类型免疫细胞(jurkat, lm-9, raji, thp-1)的基因表达谱，并且每个数据集具有来自三个样本的测量值。我们对这三个样本的基因表达值做对数变换，计算对数变换后的每一个基因在三个样本上的均值和方差，并用估计的均值和方差根据正态分布生成了w。最后，矩阵y为这两个矩阵w和h的乘积，且加上小高斯噪声。图2的展示结果来自于100个蒙特卡罗实验结果的总结。

图2(a)展示的是模拟结果。由于在模拟时真实的比例矩阵是生成的，因此是已知的。我们可以计算真实比例矩阵和估计比例矩阵之间的mc和mae。更高的mc值和更低的mae值表明算法具有更优的性能。图2(a)展示了在迭代次数为0时估计比例与真实比例之间的皮尔逊相关系数(此时的估算比例是使用观测数据中变异系数最高的前1000个特征作为反卷积的输入)和使用decd方法进行6次、13次和20次迭代后估计比例与真实比例在四种细胞类型上的皮尔逊相关系数。图2(b)、图2(c)分别显示了四种细胞类型中估计比例与真实比例之间的平均皮尔逊相关系数和平均绝对误差。“before”时的估计细胞类型比例矩阵实际上就是迭代次数为0时的估计比例矩阵，而decd的估计比例矩阵是基于30次迭代中rmse最小的迭代对应的估计比例矩阵。显然，在所有四种细胞类型上，随着迭代次数的增加，估计比例与真实比例之间的相关性总体上不断提高。与“before”(迭代次数为0)相比，细胞类型比例估算精度的提升非常显著。四种细胞类型上的平均相关性显著增加(p值为3.16e−34)，平均绝对误差显著减少(p值为2.89e−34)。

figure 2. performance of decd on synthetic mixtures. (a) boxplots of pearson correlations between estimated and true proportions for four cell types by number of iterations; (b) boxplots of mean pearson correlations between estimated and true proportions across four cell types from “before” and decd; (c) boxplots of mean absolute errors between estimated and true proportions across four cell types from “before” and decd

图2. decd在合成混合物上的性能表现。(a) 不同迭代次数下，四种细胞类型上估计比例与真实比例之间的皮尔逊相关系数的箱线图；(b) decd和“before”方法在四种细胞类型上的平均皮尔逊相关系数箱线图。(c) decd和“before”方法在四种细胞类型上的平均绝对误差箱线图

我们从基因表达综合数据库(geo)上下载了2个基因表达数据集：mouse-mix数据和immune数据，这两个数据具有实验算出的真实比例矩阵h。mouse-mix数据由3种细胞类型组成，immune数据有4种细胞类型组成，每一个细胞类型在每个样本上所占比例都是已知的，即真实比例矩阵h是已知的。

首先，我们利用deconf函数分别计算出mouse-mix和immune数据的估计比例矩阵，即图3中“before”对应的每一个细胞类型在每一个样本上的估计比例。然后，利用decd计算出这两个数据集上的估计比例矩阵。图3(a)和图3(b)显示了mouse-mix和immune数据集上每种细胞类型在初始点(before，即迭代次数为0)和应用decd后的估算比例与真实比例在每一个样本上的散点图，并计算了各细胞类型的平均皮尔逊相关系数。显然，应用decd后细胞比例矩阵估计的准确性显著提升，平均相关系数(mc)在mouse-mix数据集上从0.75增加到0.995，在immune数据集上从0.451增加到0.963。

figure 3. scatter plots of estimated and true proportions of mouse-mix and immune data for each cell type at initial point and after applying decd

图3. mouse-mix和immune数据集中上每种细胞类型在初始点和应用decd后的估计比例与真实比例的散点图

4. 结论

本文开发了一种基于最优特征选择的无参考反卷积算法decd，该算法通过跨细胞类型的差异分析，迭代搜索细胞类型特异性的特征，将一种细胞类型与其他细胞类型之间的差异分析和两种细胞类型与其他细胞类型之间的差异分析相结合筛选特征，进而估算细胞组成。decd完全依赖数据驱动，不需要任何关于细胞类型或其比例的先验知识，并且不局限于特定的无参考反卷积方法。目前，大多数现有的不依赖参考基矩阵的反卷积算法都可以与该过程结合使用。模拟研究和两个真实数据集上的实证分析证明了decd的优越性。

基金项目

研究得到了江西省自然科学基金(20212bab202001)的支持。

notes

^*通讯作者。

参考文献

[1]	zhang, w., li, z., wei, n., wu, h. and zheng, x. (2019) detection of differentially methylated cpg sites between tumor samples with uneven tumor purities. bioinformatics, 36, 2017-2024.
[2]	giannakopoulos, p., herrmann, f.r., bussière, t., bouras, c., kövari, e., perl, d.p., et al. (2003) tangle and neuron numbers, but not amyloid load, predict cognitive status in alzheimer’s disease. neurology, 60, 1495-1500.
[3]	yang, y., mufson, e.j. and herrup, k. (2003) neuronal cell death is preceded by cell cycle events at all stages of alzheimer's disease. the journal of neuroscience, 23, 2557-2563.
[4]	li, b., severson, e., pignon, j., zhao, h., li, t., novak, j., et al. (2016) comprehensive analyses of tumor immunity: implications for cancer immunotherapy. genome biology, 17, article no. 174.
[5]	sturm, g., finotello, f., petitprez, f., zhang, j.d., baumbach, j., fridman, w.h., et al. (2019) comprehensive evaluation of transcriptome-based cell-type quantification methods for immuno-oncology. bioinformatics, 35, i436-i445.
[6]	newman, a.m., liu, c.l., green, m.r., gentles, a.j., feng, w., xu, y., et al. (2015) robust enumeration of cell subsets from tissue expression profiles. nature methods, 12, 453-457.
[7]	teschendorff, a.e., breeze, c.e., zheng, s.c. and beck, s. (2017) a comparison of reference-based algorithms for correcting cell-type heterogeneity in epigenome-wide association studies. bmc bioinformatics, 18, article no. 105.
[8]	hattab, m.w., shabalin, a.a., clark, s.l., zhao, m., kumar, g., chan, r.f., et al. (2017) correcting for cell-type effects in dna methylation studies: reference-based method outperforms latent variable approaches in empirical studies. genome biology, 18, article no. 24.
[9]	gong, t., hartmann, n., kohane, i.s., brinkmann, v., staedtler, f., letzkus, m., et al. (2011) optimal deconvolution of transcriptional profiling data using quadratic programming with application to complex clinical blood samples. plos one, 6, e27156.
[10]	houseman, e.a., molitor, j. and marsit, c.j. (2014) reference-free cell mixture adjustments in analysis of dna methylation data. bioinformatics, 30, 1431-1439.
[11]	rahmani, e., schweiger, r., shenhav, l., wingert, t., hofer, i., gabel, e., et al. (2018) bayescce: a bayesian framework for estimating cell-type composition from dna methylation without the need for methylation reference. genome biology, 19, article no. 141.
[12]	kang, k., meng, q., shats, i., umbach, d.m., li, m., li, y., et al. (2019) cdseq: a novel complete deconvolution method for dissecting heterogeneous samples using gene expression data. plos computational biology, 15, e1007510.
[13]	li, z. and wu, h. (2019) toast: improving reference-free cell composition estimation by cross-cell type differential analysis. genome biology, 20, article no. 190.
[14]	rahmani, e., zaitlen, n., baran, y., eng, c., hu, d., galanter, j., et al. (2017) correcting for cell-type heterogeneity in dna methylation: a comprehensive evaluation. nature methods, 14, 218-219.
[15]	houseman, e.a., kile, m.l., christiani, d.c., ince, t.a., kelsey, k.t. and marsit, c.j. (2016) reference-free deconvolution of dna methylation data and mediation by cell composition effects. bmc bioinformatics, 17, article no. 259.

为你推荐

凯发国际一触即发的友情链接