文档介绍:该【迭代式特征选择的单细胞分化轨迹推断算法 何鸿坚 】是由【司棋】上传分享,文档一共【14】页,该文档可以免费在线阅读,需要了解更多关于【迭代式特征选择的单细胞分化轨迹推断算法 何鸿坚 】的内容,可以使用淘豆网的站内搜索功能,选择自己适合的文档,以下文字是截取该文章内的部分文字,如需要获得完整电子版,请下载此文档到您的设备,方便您编辑和打印。:.
计算机科学与探索
JournalofFrontiersofComputerScienceandTechnology
ISSN1673-9418,CN11-5602/TP
《计算机科学与探索》网络首发论文
题目:迭代式特征选择的单细胞分化轨迹推断算法
作者:何鸿坚,殷依婷,谢江
网络首发日期:2022-08-15
引用格式:何鸿坚,殷依婷,[J/OL].计
算机科学与探索.
.
网络首发:在编辑部工作流程中,稿件从录用到出版要经历录用定稿、排版定稿、整期汇编定稿等阶
段。录用定稿指内容已经确定,且通过同行评议、主编终审同意刊用的稿件。排版定稿指录用定稿按照期
刊特定版式(包括网络呈现版式)排版后的稿件,可暂不确定出版年、卷、期和页码。整期汇编定稿指出
版年、卷、期、页码均已确定的印刷或数字出版的整期汇编稿件。录用定稿网络首发稿件内容必须符合《出
版管理条例》和《期刊出版管理规定》的有关规定;学术研究成果具有创新性、科学性和先进性,符合编
辑部对刊文的录用要求,不存在学术不端行为及其他侵权行为;稿件内容应基本符合国家有关书刊编辑、
出版的技术标准,正确使用和统一规范语言文字、符号、数字、外文字母、法定计量单位及地图标注等。
为确保录用定稿网络首发的严肃性,录用定稿一经发布,不得修改论文题目、作者、机构名称和学术内容,
只可基于编辑规范进行少量文字的修改。
出版确认:纸质期刊编辑部通过与《中国学术期刊(光盘版)》电子杂志社有限公司签约,在《中国
学术期刊(网络版)》出版传播平台上创办与纸质期刊内容一致的网络版,以单篇或整期出版形式,在印刷
出版之前刊发论文的录用定稿、排版定稿、整期汇编定稿。因为《中国学术期刊(网络版)》是国家新闻出
版广电总局批准的网络连续型出版物(ISSN2096-4188,CN11-6037/Z),所以签约期刊的网络版上网络首
发论文视为正式出版。
:.
췸싧쫗랢쪱볤ꎺ2022-08-1512:29:01
췸싧쫗랢뗘횷ꎺ.
计算机科学与探索
JournalofFrontiersofComputerScienceandTechnologydoi:.1673-
何鸿坚,殷依婷,谢江+
上海大学计算机工程与科学学院,上海200444
+通信作者E-mail:******@
通过单细胞轨迹推断方法从单细胞转录组学数据或蛋白质组学数据构建细胞的分化轨迹,有助于
理解正常组织的发育过程或者提供病理学相关的信息。然而当前的单细胞轨迹推断算法在精确度和鲁棒
性的提升上仍然是一个难题,原因之一是在单细胞测序中检测到大量不相关的基因而产生的噪声。针对
这一问题,迭代式特征选择的轨迹推断方法iterTIPD被提出。其创新点体现在,将广泛用于筛选差异表
达基因的特征选择方法迭代式地用于线性或分支结构的单细胞RNA测序数据上,通过筛选出对构建的
分化轨迹贡献最大的基因子集来提高细胞伪时间排序的精确度和鲁棒性。在4种scRNA-seq数据集上的
实验结果表明,iterTIPD可以有效地提高单细胞轨迹推断算法的精确度和鲁棒性。同样,iterTIPD也使
其他的轨迹推断算法的性能得到的提升,以此证明了iterTIPD具有泛化性。iterTIPD算法成功重构了神
经干细胞的分化轨迹,并通过对比发现,该分化轨迹与已知的神经干细胞分化轨迹高度一致。同时,发
现Top2a和Gja1可能是定义活化的神经干细胞亚群的新的标志物。
单细胞RNA测序技术;基因差异性表达;单细胞分化轨迹推断;迭代式特征选择;生物信息学
ATP399
Single-cellDifferentiationTrajectoryInferenceAlgorithmwithIterativeFeature
Selection
HEHongjian,YINYiting,XIEJiang+
SchoolofComputerEngineeringandScience,ShanghaiUniversity,Shanghai200444,China
Abstract:Theconstructionofcelldifferentiationtrajectoriesfromsingle-celltranscriptomicdataorproteomic
databysingle-celltrajectoryinferencemethodscanhelptounderstandthedevelopmentalprocessofnormal
,theaccuracyandrobustnessofcurrentsingle-
celltrajectoryinferencealgorithmsarestillachallenge,oneofthereasonsisthenoisecausedbythedetectionof
alargenumberofunrelatedgenesinsingle-,atrajectoryinference
usedtoscreendifferentiallyexpressedgenesofiterativefeatureselectionmethodforlinearorbranchofthe
structureofthesinglecellRNAsequencingdata,throughselectedfortheconstructionofthelargestcontribution
experimentalresultsonfourscRNA-seqdatasetsshowthatiterTIPDcaneffectivelyimprovetheaccuracyand
robustnessofthesingle-
trajectoryinferencealgorithms,provingiterTIPD'
———————————————
基金项目:国家自然科学基金(61873156)。
ThisworkwassupportedbytheNationalNaturalScienceFoundationofChina(61873156).
:.
2JournalofFrontiersofComputerScienceandTechnology
cellswasreconstructedbyiterTIPDalgorithm,andthecomparisonshowedthatthedifferentiationtrajectorywas
,Top2aandGja1maybe
novelmarkersdefiningactivatedneuralstemcellsubpopulations.
KeywordsSingle-cellRNASequencing;DifferentialGeneExpression;Single-cellDifferentiationTrajectory
Inference;IterativeFeatureSelection;Bioinformatics
单细胞RNA测序(Single-cellRNASequencing,成功地重建了已知的分化路径,并且还发现了新
scRNA-seq)技术的出现为探索单个细胞水平的转的候选调控因子。
录组提供了前所未有的机会,已经成为研究细胞综上所述,基于scRNA-seq数据的轨迹推断
异质性[1]和早期胚胎发育[2]等关键生物学问题的理可以揭示细胞分化过程中的联系并且为基因表达
想选择。对scRNA-seq数据进行分析的方法有很的异质性和细胞动力学提供深刻的见解。因此,
多,包括质量控制[3],批次效应校正[4],数据标准基于scRNA-seq数据的轨迹推断方法的研究是一
化[5],特征选择[6],降维[7],聚类[8],细胞亚群的项具有重要意义的工作。
鉴定[9]和轨迹推断[10-13]等。其中,单细胞的分化轨然而,由于scRNA-seq数据具有高噪声、高
迹推断方法是研究的热点之一。异质性等特点[19],并且存在辍学效应(由于技术
传统上研究谱系关系是通过在细胞中注入遗限制,低表达基因可能难以检测)[14],因此提出
传标记,然后追踪细胞的后代实现的[14]。然而,一个更准确、更鲁棒的细胞分化轨迹方法仍然是
有限的遗传标记可能会掩盖细胞种群中的异质一个挑战。目前已有超过70多种单细胞轨迹推断
性。单细胞的分化轨迹推断方法则是从单细胞转的方法[20]。为了对scRNA-seq数据中的高维噪声
录组学数据或蛋白质组学数据中自动重建细胞的进行处理,特征降维和特征选择算法被大量使用
分化轨迹(谱系关系),并根据其在分化轨迹中的在单细胞轨迹推断中。
进展对单个细胞进行排序,探索基于单个细胞间在使用了特征降维的单细胞轨迹推断方法
的基因的动态表达。这些计算方法为每个细胞分中,Monocle[10]是第一个也是最经典的单细胞轨迹
配一个所谓的伪时间(Pseudotime),根据伪时间对推断算法。但是,随着最新的scRNA-seq技术的
细胞的排序称为伪时间排序。伪时间是细胞沿着发展,一次可测量成千上万个细胞,连接大量细
生物系统中连续发展的动态过程的顺序,它可以胞的最小生成树(MinimumSpanningTree,MST)则
定义单个细胞在其分化的过程中所经历的不同的变得非常复杂且难以解释,并且Monocle无法预
过渡阶段。因此伪时间排序可以识别轨迹开始和测分支轨迹。随后在2017年,Monocle2[12]被提
结束状态的细胞类型以及中间阶段的细胞[15]。出,它在低维空间中构建谱系树,保留了更多原
基于单细胞测序数据的轨迹推断方法主要应始数据中关于分化轨迹的信息,解决了预测分支
用在以下两个方面:轨迹的问题。Monocle和Monocle2算法已经被用
(1)识别过渡状态:在许多生物系统中,细胞会于推断各种类型的细胞分化轨迹,其中包括一些
表现出连续的状态的转变,以不同的转录、形神经元谱系轨迹的推断[21,22],造血系统中的祖细
态、表观基因组以及表面标志物的变化为特征胞和干细胞的分化[12]以及胎盘的发育[23]等研究。
[16]。使用基于单细胞测序数据的轨迹分析能够提而后SCUBA[24]使用t分布随机邻居嵌入(t-
供一种更直接,更公正的方法来识别和正确排序distributedStochasticNeighborEmbedding,t-
不同的过渡阶段,即基于scRNA-seq数据重建细SNE)[25]降低数据的维数,然后对降维后的数据拟
胞的谱系轨迹,发现新的过渡状态[17]。合平滑的曲线。SCUBA已用于小鼠小肠中的干细
(2)识别关键调控因子:除了揭示细胞间的基因胞分化过程中的关键基因的计算[26]。Slingshot[13]
表达动态,单细胞轨迹推断还可以帮助识别触发以降维后的数据和聚类的结果作为输入构造
状态转变的关键调控因子。例如,在人类定形内MST,然后在MST各个分支上拟合平滑的曲线来
胚层(DefinitiveEndoderm)细胞发育的研究中[18],进一步细化此树。最后将所有细胞投影到其最接
使用单细胞轨迹推断将细胞沿着发育轨迹排列,近的曲线上,从而产生带有分支的有序的谱系轨:.
3
迹。Slingshot最近已被用来预测嗅觉干细胞谱系分化上的实际应用进行描述。
轨迹中的细胞命运和分支点[27]。1iterTIPD
不同于特征降维的方法,使用特征选择方法iterTIPD是一种能够广泛运用在线性或分支结
可以降低基因的维度并识别最相关的特征,提高构的单细胞RNA测序数据上的筛选差异表达基因
信噪比和下游分析(如聚类或伪时间推断)的计的迭代式的特征选择算法。通过筛选出对构建的
算效率。有效的特征选择方法可以使研究者专注分化轨迹贡献最大的基因子集来提高细胞伪时间
于关键基因的分析,从而提高基于基因表达谱相排序的精确性和鲁棒性,从而提高单细胞轨迹推
似性对细胞进行伪时间排序的轨迹推断方法的性断算法的性能。下文基于TIPD[28]算法来阐明本文
能。目前一些用于单细胞轨迹推断中的特征选择提出的iterTIPD算法。
来筛选高度可变基因,从而减少数据的维度。但TIPD算法通过计算所有单细胞的信号熵(SR)
是由于单细胞数据中的辍学效应,直接从数据中来表示单个细胞的分化能力并由此计算出细胞种
筛选高度可变基因可能会遗漏一些重要的携带有群中的异质状态的分布。SR可用于定量描述单个
效信息的低表达基因。大量细胞RNA测序(bulk细胞的分化状态,高SR值意味着高的分化能力,
RNASequencing,scRNA-seq,bulkRNA-seq)的实验而低SR值代表低的分化能力。TIPD结合信号熵
中通常将筛选差异表达基因作为特征选择方法,和基因表达数据的聚类结果,从而描述细胞群中
但是scRNA-seq数据一般没有细胞亚群的信息异质状态的概率分布情况。使用对称的JS散度
(实际上这也是分析scRNA-seq数据的目标之(Jensen-ShannonDivergence,JSD)衡量细胞簇间的
一),因此差异表达基因的筛选很难直接应用于概率分布距离,然后在这些细胞簇中心和概率分
scRNA-seq数据。布距离构成的完全图上计算MST。以上步骤确定
基于以上提到的问题,本文在一种基于概率了细胞分化轨迹的“骨架”,然后将所有细胞向曲
分布的轨迹推断方法——TIPD算法[28]上提出了迭
线拟合后的“骨架”上投影,最终就得到了细胞
代式特征选择的优化方法(IterativeTrajectory的伪时间排序。由此构造的MST是在信号熵的背
InferenceBasedonProbabilityDistribution,景下基于系统级别的信息计算得到的,因此它可
iterTIPD)。前提假设为:如果一个基因对最终构建以在全局上反映真实的生物过程,使得TIPD算法
的单细胞分化轨迹起作用的话,它应该沿着细胞相较于其他轨迹推断的方法展现出了更高的准确
的轨迹差异表达[12]。
性。
本文的工作和贡献可总结如下:
(1)将iterTIPD在四个单细胞数据集上运行,
单细胞轨迹推断方法的输出是按照细胞分化
通过实验确定了合适的迭代参数。并且,为了探
的伪时间进行的排序,可以认为单个细胞之间就
究这种迭代式特征选择方法的有效性,本文将
是处于不同的状态,因此可以在伪时间顺序上计
iterTIPD,TIPD以及另外两种常用的特征选择方
算差异表达基因作为scRNA-seq数据的特征选择
法在精确度、鲁棒性和检测金标准基因的能力上
方法。这里我们使用TIPD算法进行伪时间的计算
进行了比较。
(图1中绿色虚线框标注的部分)。
(2)基于优化后的iterTIPD算法成功重构了神
已有的研究表明,广义加性模型(Generalized
经干细胞的分化轨迹,并通过对比发现,该分化
AdditiveModel,GAM)[10]被用来计算伪时间上差异
轨迹与已知的神经干细胞分化轨迹高度一致。
(3)发现了在活化的神经干细胞种群的三种不表达的基因和检测金标准基因。我们可以通过
同状态下差异性表达的基因Top2a和Gja1,这可GAM拟合基因的表达与伪时间之间的关系,来进
能是活化的神经干细胞种群的新的标志物。行特征选择。首先通过Tobit模型[29]对每个基因在
本文将在第1节中介绍iterTIPD算法。第2节细胞中的表达水平进行建模,然后每个基因的表
在四个单细胞数据集上进行实验,并与其他算法达与伪时间之间的关系都使用GAM模型拟合,这
在对精确度、鲁棒性和检测金标准基因的能力上里每个基因的表达水平𝒀取决于潜在变量𝒀∗,表达
进行比较。第3节中将对iterTIPD在神经干细胞式为::.
4JournalofFrontiersofComputerScienceandTechnology
值;𝑭𝑫𝑹𝒑−𝟏是通过Benjamini-Hochberg对𝒑−𝟏进
𝒀∗𝒀∗≥𝝀
行FDR校正[30]后获得的基因的p值。我们在每
𝒀={(1)
𝝀𝒀∗<𝝀一次迭代中保留𝑭𝑫𝑹𝒑−𝟏
一次迭代,并认为这些基因是差异表达基因。以
其中𝝀是检测的阈值,默认情况下𝛌设置为
上过程中关于GAM模型的实现和相关的测试函数
。潜在变量𝒀∗取决GAM模型中的变量𝒙,𝒙在
𝒊𝒊由“VGAM”包[31]提供。
此处表示每个细胞对应的伪时间值,𝒊∈
以上步骤计算出的伪时间顺序上差异表达的
{𝟏,…,𝒏}。因此,计算伪时间上差异表达基因的
基因将形成新的特征子集进行新一轮的计算。迭
GAM模型的表达式定义为:
代的次数由指定的参数控制,默认为100,可根据
不同数据集在迭代过程中精确度的变化趋势来决
定()。算法参照迭代参数进行多次迭
代后,在精确率呈现明显下降趋势的时候停止迭
代,输出所有计算中精确度最高的细胞伪时间排
序的结果以及对应的基因子集。
如果轨迹推断方法构建的轨迹是线性轨迹,
如图2A,则直接计算伪时间顺序的差异表达基因
作为筛选的特征子集;如果是分支轨迹,如图2B所
示是含有一个分支结构的分化轨迹,则分别计算
每一条谱系(例如a-b-c和a-b-d)伪时间顺序上差
异表达的基因,最后求这些差异表达基因的并集
作为筛选的特征子集。因为每一轮使用的特征子
集都是基于上一轮的基因集合计算出来的基因子
集,随着迭代次数的增加,筛选的基因的个数逐
渐减少,可能会丢失一些包含有效信息的基因,
因此迭代的次数不宜过多。
图1iterTIPD算法流程图
𝐸(𝑌)=((𝑥,))+𝜀(2)
图2线性轨迹和分支轨迹的细胞伪时间排序示意图
其中,𝒈𝒊表示单个细胞,𝒊∈{𝟏,…,𝒏},-temporal
𝝍(𝒙,𝒈𝒊)代表细胞𝒈𝒊被分配的伪时间值,表示有orderingoflinearandbranchingtrajectories
效自由度为3的三次平滑函数。误差项𝜺是正态分
布的,平均值为零。最后通过式(3)迭代筛选差异2
表达基因。
𝐺𝑡={|𝑭𝑫𝑹𝒑−𝟏<,∈𝐺𝑡−1}(3)
在本文中,使用了四个scRNA-seq数据集,
𝑡−𝟏其规模大小不同并且包括线性结构和分支结构两
其中𝐺表示迭代筛选次后的基因集合;𝒑为在
𝑡−12种类型,物种涵盖人类和小鼠。它们分别是:
𝐺上利用近似似然比检验计算出基因的p
(1)小鼠肺泡2型细胞(LungAlveolarType:.
5
2,AT2)2
AT2数据集包含101个小鼠肺泡2型细胞,这𝑆𝜋1,𝜋2=∑ℎ(𝜋1,𝜋2,,)(4)
|𝐴|(|𝐴|−1),∈𝐴;≠
些细胞是在胚胎小鼠的肺中的四个时间点(、
[32]其中𝜋1和𝜋2表示两个伪时间的排序,𝐴是𝜋1和𝜋2中
、)收集的。该分化轨迹为线
所有细胞的并集,|𝐴|表示集合𝐴的基数。对于一对
性结构。
细胞和,如果它们在𝜋1中的前后顺序与在𝜋2中相
(2)人类骨骼肌成肌细胞(HumanSkeletal
MuscleMyoblasts,HSMM)同,那么ℎ(𝜋1,𝜋2,,)=1,否则ℎ(𝜋1,𝜋2,,)=0。
鲁棒性评分越高表示两个伪时间排序越相似。
HSMM数据集包含271个人类骨骼肌成肌细
本文以金标准基因的平均排名[11]作为衡量指
胞,包括三种细胞类型:增殖细胞(Proliferating
标来测试各特征选择方法能否提高检测金标准基
Cells)、分化的成肌细胞(DifferentiatingMyoblasts)
因的能力。所谓金标准基因是对于给定的基因表
和间质间充质细胞(InterstitialMesenchymal
Cells)[10]。该分化轨迹为线性结构。达数据,通过查阅文献收集已知在该生物过程中
(例如细胞分化或细胞的发育过程)差异表达的
(3)小鼠胚胎干细胞(EmbryonicStemCells,
基因。计算得到的所有差异表达基因按照FDR校
ESC)[31]
正后的p值由小到大的顺序排列,差异表达基
ESC数据集是白血病抑制因子(Leukemia
因中金标准基因的排名则可以反映该算法检测已
InhibitoryFactor,LIF)撤离后在四个时间点(d0、
d2、d4和d7)收集的2717个小鼠胚胎干细胞[33]。知差异表达基因的能力,排名越靠前说明算法检
测已知差异表达基因的能力越高。
该分化轨迹为线性结构。
(4)小鼠常规树突状细胞(
Cells,cDC)iterTIPD算法的特点之一在于迭代式地筛选特
cDC数据集包含251个小鼠常规树突状细征基因,迭代次数过少可能使筛选的基因子集包
胞,由树突状细胞祖细胞(macrophagedendritic含大量无效基因从而使得算法的效果不明显,迭
cellsprogenitor,MDP),普通DC祖细胞(common代次数过多可能导致重要的特征基因被过滤掉。
dendriticcellsprogenitor,CDP)和前树突状细胞(pre-
DC)组成[34]。小鼠的cDC由两个主要的谱系cDC1据集上,初步设定迭代次数为100次,观察每一
和cDC2组成[35,36],每条谱系都由不同的转录程序个数据集上迭代次数与精确度之间的变化情况,
控制,这些程序从巨噬细胞和MDP开始分化,接结果如图3所示。迭代次数为0时对应的精确度是
着是正常的CDP,然后是preDC,最后分别分化初始的精确度,即TIPD得到的精确度。黄色虚线
成cDC1谱系和cDC2谱系[34]。该分化轨迹为分支标注的是迭代过程中取得的最高精确度以及对应
结构。的迭代次数,红色虚线是使用多项式回归拟合实
。由图3可见在迭代过程
在单细胞分化轨迹推断的研究中,精确度、中,四个单细胞数据集上的精确度上下波动,并
且随着迭代次数的增加呈现先升高后降低的总体
鲁棒性和检测金标准基因的能力三个方面常用来
[11]趋势。精确度上下波动是由于iterTIPD算法中的
全面衡量算法的性能,本文也使用这三个指标
降维步骤也会对精度造成影响。精确度的趋势呈
来评价轨迹推断算法的性能。
现先升高后降低,这是由于当数据包含较多冗余
在精确度方面,本文使用细胞的伪时间排序
和数据集中记录的实验阶段的信息(即外部信特征时,iterTIPD迭代次数的增加会缩小特征子集
从而提高算法的精确度,当精确度达到最高后,
息)之间的肯德尔等级(KendallRank)相关系数来
[37,38]特征的数量继续减少使得携带有效信息的基因被
进行衡量。
排除在外,导致精确度呈现下降趋势。迭代次数
鲁棒性可通过计算原始数据和被扰动数据的
[11]接近100次时,四个单细胞数据集上精确度的下
伪时间排序的鲁棒性评分来评估。扰动的数据
集是通过从原始数据集中随机二次抽样90%,降趋势都较为明显。因此,默认情况下,iterTIPD
的迭代次数设置为100。该参数为后续研究中使用
80%和70%的细胞生成的,每次重复50次。鲁棒
iterTIPD在其他数据集上进行细胞分化轨迹推断实
性评分的公式为::.
6JournalofFrontiersofComputerScienceandTechnology
验提供参考依据。迭代完成后,iterTIPD将输出迭最高精确度的迭代次数以及对应的最高精确度)。
代过程中精确度最高的细胞伪时间排序的结果以总体来看,iterTIPD在四个单细胞数据集上均有较
及对应的基因子集,图3显示了iterTIPD在AT2、好的优化效果,尤其是在ESC数据集上,iterTIPD
HSMM、ESC和cDC数据集上分别在第40、56、的精确度提升了大约9%。
43和62次达到最高的精确度(括号里标注了达到
图3iterTIPD迭代次数与精确度的关系
%,%,%和4%,说明
本文选择了两个流行的软件包“Monocle”包iterTIPD迭代式特征选择方法有效地提高了TIPD
[10]和“Seurat”包[39]提供的筛选高度可变基因的精确度,并且在除了ESC数据集的其他三个单
(HighlyVariableGenes,HVGs)的特征选择方法与细胞数据集上的精确度均为最高。M_HVGs+TIPD
TIPD算法结合(分别命名为M_HVGs+TIPD和方法只在ESC数据集上对TIPD的精确度有提
S_HVGs+TIPD),来比较iterTIPD的迭代特征选择升,S_HVGs+TIPD方法只在ESC和cDC数据集
方法的优化效果。其中,Monocle和Seurat中用于上对TIPD的精确度有提升。
筛选HVGs的函数参数均使用默认参数,其筛选表1使用不同的特征选择方法的精确度比较
Table1Comparisonofaccuracyusingdifferentfeature
HVGs