[0001] 本公开涉及用于获得用于预测来自基因转录物的时变特性的预测器的计算机实现方法以及对应的计算机程序产品、计算机可读介质和系统。具体地但非排他地，本公开涉及用于疾病、病况或病症的进展、生物学年龄或实足年龄的预测器。本公开的另一个非排他性方面涉及估计对不同基因或基因转录物的预测的贡献。

[0002] 许多疾病具有衰老成分，例如帕金森病、阿尔茨海默病和骨关节炎。越来越有兴趣寻找通过新疗法诱导细胞和组织再生的方法，这些新疗法可以释放非常年轻的细胞中存在的潜在再生能力。在过去五年中，在部分表观遗传重编程的领域中的科学已经取得了许多
进展，该领域前景广阔。

[0003] 此前已知的诱导细胞变换为更年轻状态的唯一综合方法是创建iPSC(例如使用山中(Yamanake)因子)。不幸的是，经历这种向多能性转变的细胞也会改变它们的身份，因此该技术不能用于创建抗衰老疗法，也不能用于延长健康寿命。

[0004] 现在从部分表观遗传重编程的研究中得知，年龄逆转成分可以与细胞身份成分解耦合，并且目前正在努力将该过程变换为临床。

[0005] 衰老时钟是一种了解如何驱动细胞再生过程的优雅方式。第一衰老时钟是由Horvath等人开发的(例如参见US20160222448A1和US20190185938A1)并且是基于甲基化水
平，因此被描述为表观遗传时钟。尽管表观遗传时钟能够高度准确地预测年龄，但它也有一些局限性，包括难以做出生物学推论，以及目前无法验证或针对个体位点以获得潜在的治
疗益处。因此，注意力转向了转录组时钟，它根据RNA表达水平预测年龄。例如在
US10325673B2中和由Holzscheck等人(npj Aging Mech Dis 7,15(2021))已经描述转录组
时钟。然而，这些转录组时钟的一个显著特征是它们在对应的基因通路的总和的转录水平
上运行，因此需要预先了解基因通路才能制作此类时钟。发明人已经认识到这具有许多缺
点，如下所述。因此，本领域需要一种克服这些限制的时钟(衰老的预测器)。

[0031] 参考图1，在步骤110针对一批次的细胞样本接收基因转录数据集。细胞样本和数据集可以以任何合适的方式获得，例如如上所述。该数据集是根据原始基因转录计数(个体转录物的计数或按基因求和)、每个转录物或基因的一个计数、细胞样本和测量时间点生成的，以给出每个细胞样本和时间点的表达水平的一个表达载体。数据样本可以从单个细胞
(单个细胞样本)的计数获得，或者可以从几个细胞的样本的池化计数获得。数据集被处理
或已经被处理以使用计数数据的常规数值条件来得到表达水平，包括标准化数据、对数变
换数据以及将经对数变换的数据标准化以具有零均值和单位标准差，例如，注意到总体比
例因子当然可以随意改变。至关重要的是，每个转录水平因此是针对个体基因转录物或针
对该基因池化(例如，求和)的个体基因的基因转录物。因此，在生成表达载体的过程中不需要基因通路或生物活性的现有知识。此外，后续描述的处理可以在不使用基因通路或生物
活性的现有知识的情况下完成。当然，应当理解，附加的处理步骤，诸如回归后对预测系数的进一步调整可以使用此类现有知识或可以使用从回归本身获得的生物学知识，例如关于
某些基因如何预测时变特性，这可以例如从预测系数得到，如下所述。

[0032] 术语“基因通路”是指共同发挥作用以执行具体生物过程的基因网络。此类生物过程也可称为“生物通路”，即细胞中分子之间的一系列相互作用，其导致生物效应，例如细胞的变异或产物的产生。这些分子是由基因编码的，因此可以看出，基因通路中基因网络的结果将是生物通路。可以获得有关基因通路和生物通路的知识，例如根据“Hallmark”通路汇集(Liberzon,A.et al.The molecular signatures database Hallmark Gene Set Collection.Cell Syst.1,417–425(2015))或公开地可用的数据库(诸如KEGG通路数据库
(https：///www.kegg.jp/))。

[0033] 所得到的基因转录数据集被组织(或接收)为矩阵E，该矩阵E具有作为行向量的转录向量，其中每个细胞样本和时间点一行。使用诸如特征分解或更典型的奇异值分解之类
的任何合适的技术来找到特征向量矩阵W和特征向量的对角矩阵(用基向量W变换的E的方
差)。

[0034] ETEW＝WΛ              等式1

[0035] 在步骤120，使用与最大特征值(或所解释的方差)Λi，i相关联的k列特征向量(或主成分)Wi的矩阵υ和对角缩放矩阵S形成嵌入数据集X，

[0036] υ＝[w1 w2 w3 … wk]              等式2

[0037] 该对角缩放矩阵S通过细胞样本之间的逆标准差来缩放主成分，以便在主成分的较高方差成分与较低方差成分之间的回归的初始贡献中保持公平的竞争环境，如上所述。

[0038]

[0039] k可以被选择为合适的，具有需要更多的计算、但可能包括更多的生物学相关信息的较高的值。已找到k＝50在大多数设置中是合适的最大值，并且在一些实施方式中，k可以例如在20与30之间。还可以通过比较下面描述的针对k的不同值的系数拟合的性能并且选
择实现最佳或至少令人满意的性能的值来以迭代方式选择k。在一些实施方式中，替代选择具有k个最大特征值的成分，可以根据不同的准则来选择成分，例如在一些实施方式中基于如上所述的性能，在特征值范围的中间或在特征值的特定序数处。

[0040] 数据集的每个数据样本还包括细胞样本(或获得细胞样本的生物体)的时变特性的实际值，注意数据集包含多个时间点的多个细胞样本并且每个细胞样本和时间点有一个
此类的值。实际值可以在每个时间点测量，例如通过测量诸如指示生物学年龄或与疾病轨
迹有关的生物标志物的量，实际值可以对于生物体单独已知，诸如疾病进展或阶段，或者可能只是时间点本身，就像实足年龄一样。诸如表观遗传时钟测量之类的衰老时钟测量可以
用作指示生物学年龄的生物标志物。除了生物学年龄或实足年龄或疾病病症(例如，诸如阿尔茨海默病或帕金森病等神经退行性病症)的阶段或进展之外，可以使用细胞样本或细胞
样本来源的生物体的任何其他时变特性。

[0041] 实际值被组织或接收在列向量y中，该列向量y具有与E相同的行数，对于每个数据样本一个。通过在步骤130将嵌入数据集应用到线性预测器来针对嵌入数据集X训练线性预
测器，以通过在步骤140调整包含针对回归中的主成分的线性权重的向量β中的预测系数和*
偏移β0来预测时变特性y的值y。

[0042] y＊＝Xβ+β0          等式4

[0043] 调整系数以最小化y与y*之间的误差的测量，例如平方误差的平均值或绝对误差的中值 可以使用各种最小化方法，包括简单的最小二乘回归。在一些
实施方式中，已经找到使用弹性网络线性回归是有利的(参见Zou，H.，&Hastie，T.(2005)Regularization and Variable Selection via the elastic net；Journal of the 
Royal Statistical Society：Series B(Statistical Methodology)，67(2)，301‑320；
https：//doi.org/10.1111/j.1467‑9868.2005.00503.x，通过引用并入本文，其还讨论了可能在某些实施方式中使用的几种替代回归方法)。有利的是，弹性网络回归促进了预测系数的稀疏性，即，对于大多数系数而言，系数往往较小，其中大量系数集中在更具预测性的回归变量(在这种情况下，更具预测性的主成分)。这有助于解释主成分(以及对应的转录水平)在时变特性捕获的过程中的生物相关性。

[0044] 系数的训练和调整可以以任何合适的方式来实现。为了减少过度拟合问题，使用n折交叉验证来训练参数可能是有利的。此外，某些数据可能会保留为纯测试数据，以评估模型在未见过的数据上的性能。任何线性预测器都可以根据特定的实施方式来使用，并且可
以组合嵌入和回归步骤。可以使用的一种线性预测器是偏最小二乘法或其变体，其包括E和y两者的嵌入。然而，本公开不限于线性预测器，并且可以使用诸如前馈神经网络或循环神经网络之类的其他预测器来提供时变特性的值的预测器。值得注意的是，线性预测器的优
势不仅在于其算法简单性和效率，而且还在于预测系数的可解释性，如下所述。

[0045] 为了评估每个表达水平对时变特性的预测的贡献，即确定哪些表达水平比其他表*
达水平更具预测性，预测系数可以通过以下方式投影回表达水平的空间中的贡献系数β：

[0046]

[0047] 其中R“取消缩放”系数以补偿回归期间S的缩放。因此，β*的元素提供了与相应转录水平对应的基因或转录对于时变特性的预测程度的测量。

[0048] 在可选步骤150，可以接收新转录样本并且可以使用经训练的预测器来预测新转录样本的时变特性的值。新转录样本可以是从用于获得用于训练的样本的同一实验/事件
或实验/事件的集合获得的样本，其时变特性的值不可用，或者新转录样本可以是新获得
的。为了获得良好的预测，必须仔细控制获得新获得的样本的条件，以匹配获得训练样本的条件，以避免由于技术噪声的差异而降低预测性能而导致显著的批次效应。在许多情况下，这可能是一个挑战，下面讨论纠正此类批次效应的方法，将新的训练数据添加到现有的训
练数据中，或者将未标记的新数据与一个或多个训练数据集相结合以提高预测性能。

[0049] 在进一步可选的步骤160，可以生成报告，提供用于每个基因/转录物的β*的一个或两个元素，以允许评估它们的预测性以及一个或多个新数据样本(如果适用)的时变特性
的预测值。报告的其他元素可以是回归系数或其他拟合优度指标、残差和/或可以有助于数据和预测器的解释的任何其他量。

[0050] 用于训练预测器并且使用组合数据集进行预测的过程包括：接收第一基因转录数据集E的步骤210和接收第二(另外的)基因转录数据集 的步骤212，每个步骤均如上文步
骤110所述。在步骤214，将两个数据集组合成一个组合数据集，

[0051]

[0052] 其中{·|·}是数据集组合操作，在最简单的实施方式中是两个数据集的级联。在一些实施方式中，组合操作包括对个体数据集进行适当的标准化，例如用针对每个细胞样
本计算的余弦范数来替换表达水平 在一些实施方式中，数据集组合操作
包括对不同批次的数据集之间的差异(通常由技术噪声引起)的校正。在一些实施方式中，
从第二批次中的每个数据样本中减去批次校正向量，或者按照批次校正行向量的批次校正
矩阵B来减去批次校正向量，

[0053]

[0054] 然后可以在步骤220以与步骤120类似的方式形成嵌入数据集X，如下：

[0055]

[0056] 其中 知 分别是在步骤120和等式3中用C代替E解释的特征向量/主成分和特征值/方差。训练预测器的步骤230、调整预测系数的步骤240以及提供报告的步骤250则类
似于上述步骤130、140和160，并且据此应用相应的公开内容。有利的是，通过组合来自不同批次的数据集，例如来自不同实验、同一实验随时间的不同实例、具体生物体的不同个体
等，可以创建并且使用更丰富的数据集来获得改进的预测器。

[0057] 在一些实施方式中，组合步骤214、等式6和7的数据集包括将数据集变换到不同的坐标系。在特定的实施方式中，找到组合数据集的主成分，并且使用与k个最大特征值相关联的组合数据集的主成分的矩阵 对组合数据集进行变换，

[0058]

[0059] 然后如上所述使用变换后的数据集。

[0060] 计算组合数据集的主成分包括：以要合并的每个数据集的均值的平均值为中心(而不仅仅是组合数据集的均值)，并且通过要合并的相应数据集中的细胞样本的数目的倒
数来加权每个细胞样本对协方差矩阵的贡献(或者等效地，通过使用要合并的数据集的协
方差矩阵的平均值作为用于主成分分析的协方差矩阵)。然后以常规方式，例如使用特征分解或奇异值分解，针对组合数据集计算主成分。

[0061] 然后，如上文参考等式6和7所述对组合数据集的选定主成分进行批次校正。在这些实施方式中，组合数据集的维度通过批次校正保持正交，并且虽然可以使用等式8来形成组合的嵌入数据集X，但却没有必要这样做，并且可以使用相同选定的k个维度的C来形成组合的嵌入数据集，

[0062]

[0063] 其中Vi，i是C的协方差矩阵V的非零对角线项。当然，可以选择较少数目的维度的C。

[0064] 用于计算批次校正向量B的各种方法是已知的并且可以在实施方式中使用。在一些实施方式中，使用相互最近邻(MNN)方法，参见Haghverdi，L.，Lun，A.T.L.，Morgan，M.D.，&Maroni，J.C.(2018)Batch effects in single‑cell RNA‑sequencing data are corrected by matching mutual nearest neighbors；Nature Biotechnology，36(5)，
421‑427；https：//doi.org/10.1038/nbt.4091和https：//marionilab.github.io/
FurtherMNN2018/theory/description.html，每一个都通过引用并入本文。通过首先针对
中的每个En创建K个最近邻的列表，并且然后针对E中的每个 创建K个最近邻的列表，定
义MNN。如果在 的列表中找到n并且在 的列表中找到n，则相应数据集中的两个细胞样本n
和 是MNN。K是基于经验或以经验为主地针对具有大量最近邻的每个数据集选择的，增加
了对噪声的鲁棒性，并且将最近邻采样到每个细胞样本云的更深处，但增加了计算成本。实际上K＝20是一个合适的选择。

[0065] 针对MNN的MNN批次校正向量是差向量 在直接基于表达水平找到MNN的实施方式中，无需正交化和/或降维，诸如PCA，如上所述。在这些实施方式中，可以使用高度可变的基因(HVG)来找到MNN，正如本领域公共的那样。虽然在一些实施方式中可以使用HVG找到MNN，但是在计算批次校正向量的这个阶段可以包括所有感兴趣的基因或所有可用的
基因，或者可以针对每个感兴趣的基因集计算单独的批次校正向量。

[0066] 上述导致产生了针对MNN的批次校正向量集，或MNN批次校正向量。然后从MNN批次校正向量中找到针对非MNN的其他数据样本的批次校正向量，例如，通过将它们与高斯
(Gaussian)核组合，使用另一种形式的加权平均，仅使用每个细胞样本的最近邻的MNN批次校正向量，等等。这为所有数据样本提供了局部变异的批次校正向量，然后将其用于如上所述的等式7中。

[0067] 在一些实施方式中，在应用如上所述的MNN批次向量之前，将每个批次中的细胞样本投影到垂直于每个数据集中的MNN批次向量的平均向量的相应平分平面上(但针对MNN细
胞样本本身的投影进行调整)。这确保了合并的细胞样本混合在一起，而不仅仅是像接触的云一样聚集在一起，即使K不够大，无法对批次的理论上面对的表面之外的最近邻进行采
样。替代地，批次校正之后的结合数据集中的细胞样本可以被投影到垂直于平均MNN批次校正向量的公共平分平面中，或者可以省略此步骤，尤其是对于足够大的K值。

[0068] 上文引用的Haghverdi等人(2018)阐述了与上述一致的批次校正方法的全部细节。，其中Haghverdi等人(2018)也描述了可用作batchelor R包的一部分的补充信息和软
件包。请参阅https://marionilab.github.io/FurtherMNN2018/theory/
description.html了解进一步的相关实施方式，该实施方式在应用批次校正向量之前压缩
每个数据集中平均MNN批次校正向量的变异，如上所述。同样可以使用输出校正数据集的降维嵌入的批次校正的替代方法，例如Seurat v3，其以与上述类似的方式在识别“锚”之前实现规范相关性分析。

[0069] 参考图3，接收第一基因转录数据集310和第二基因转录数据集312的步骤、生成组合数据集314、生成组合嵌入数据集320、训练预测器330以及调整预测系数340类似于上述
步骤210、212、214、220、230和240，并且据此应用对应的公开内容，例外的是，仅第一基因转录数据集包括时变特性的实际值y，并且在步骤312该信息不与第二基因转录数据一起接收
(或被忽略)。因此，在步骤330和340，仅使用来自特性值可用的第一数据集的数据样本来训练预测器并且调整预测系数，然后所得到的预测器被用来预测第二数据集的数据样本的特
性的相应值。通过这种方式，可以预测特性的未知值，例如，针对从生物体的新个体获得的样本，其时变特性是未知的。准备报告的步骤360类似于上述步骤160，包括第二数据集中的样本的预测值。

[0070] 所描述的实施方式使用主成分分析(例如使用SVD来实现)来计算嵌入，并且选择多个主成分以用于降维。获得嵌入的其他方法同样适用于各种实施方式并且可以用来代替
PCA来进行嵌入。例如，可以使用非线性方法(诸如核方法，例如核PCA(kPCA))或非线性方法(诸如训练自动编码器(AE))来找到嵌入。kPCA将特征分解或SVD应用于使用核函数从数据
得到的核矩阵，其方式与PCA应用于协方差矩阵的方式类似。基因的预测系数可以用与上面针对PCA描述的类似的方式被恢复，使用逆映射找到基因空间中的权重。可以通过数值优化来找到逆映射，并且可以至少近似地恢复所得到的基因预测系数。AE是经过训练以将其输
入与其输出相匹配的神经网络，并且包括隐藏嵌入层，该隐藏嵌入层的单元数少于提供嵌
入的输入层和输出层的单元数。可以使用网络的隐藏嵌入层与输出层之间的经训练的解码
网络从嵌入的预测系数至少近似地恢复基因预测系数。一般而言，通过将嵌入变换的逆映
射应用到嵌入的预测系数，可以从嵌入的预测系数中找到至少近似的基因预测系数。逆映
射可以对应于数学逆运算，或者可以是从嵌入到基因空间的任何其他操作映射，即从嵌入
数据样本到对应的数据样本的映射。因此，对(非嵌入)数据集的维度的这种投影可以是近
似的(例如通过数值方法或神经网络训练找到)或数学上精确的(例如通过矩阵求逆或转置
发现，如上面详细描述的PCA作为嵌入的情况)。

[0071] 硬件实施方式示例

[0072] 图4示出了计算设备400的一种实施方式的框图，该计算设备中具有可以执行用于使计算设备执行本文所讨论的任何一种或多种的方法中的指令集。在一些替代的实现方式
中，计算设备可以连接(例如，联网的)到局域网(LAN)、内联网、外联网或互联网中的其他机器。计算设备可以以客户端‑服务器的网络环境中的服务器或客户端机器的能力进行操作，或者作为对等(或分布式)网络环境中的对等机器进行操作。计算设备可以是个人计算机
(PC)、平板计算机、机顶盒(STB)、个人数字助理(PDA)、蜂窝电话、网络设备、服务器、网络路由器、交换机或网桥，或者任何能够执行指令集(顺序或其他)的机器，这些指令指定该机器要采取的操作。此外，虽然仅示出了单个计算设备，但术语“计算设备”还应当被理解为包括单独或联合执行一个(或多个)指令集以执行本文讨论的任何一种或多种的方法的机器(例
如，计算机)的任何汇集。

[0073] 示例计算设备400包括处理设备402、主存储器404(例如，只读存储器(ROM)、闪存、动态随机存取存储器(DRAM)(诸如同步DRAM(SDRAM)或Rambus DRAM(RDRAM)等)、静态存储器406(例如闪存、静态随机存取存储器(SRAM)等)和辅助存储器(例如，数据存储设备418)，它们经由总线430彼此通信。

[0074] 处理设备402代表一个或多个通用处理器，诸如微处理器、中央处理单元等。更特别地，处理设备402可以是复杂指令集计算(CISC)微处理器、精简指令集计算(RISC)微处理器、超长指令字(VLIW)微处理器、实现其他指令集的处理器或者实现指令集的组合的处理
器。处理设备402还可以是一个或多个专用处理设备，诸如专用集成电路(ASIC)、现场可编程门阵列(FPGA)、数字信号处理器(DSP)、网络处理器等。处理设备402被配置为执行处理逻辑(指令422)以用于执行本文讨论的操作和步骤。

[0075] 计算设备400还可以包括网络接口设备408。计算设备400还可以包括视频显示单元410(例如，液晶显示器(LCD)或阴极射线管(CRT))、字母数字输入设备412(例如，键盘或触摸屏)、光标控制设备414(例如，鼠标或触摸屏)和音频设备416(例如，扬声器)。

[0076] 数据存储设备418可以包括一个或多个机器可读存储介质(或更具体地，一个或多个非暂态计算机可读存储介质)428，其上存储了体现本文描述的任何一种或多种方法或功
能的一个或多个指令集422。指令422还可以在由计算机系统400执行期间完全地或至少部
分地驻留在主存储器404内和/或处理设备402内，主存储器404和处理器402也构成计算机
可读存储介质。

[0077] 上述各种方法可以通过计算机程序来实现。计算机程序可以包括被布置为指示计算机执行上述各种方法中的一种或多种方法的功能的计算机代码。用于执行此类方法的计
算机程序和/或代码可以在一个或多个计算机可读介质上或者更一般地在计算机程序产品
上提供给诸如计算机之类的装置。计算机可读介质可以是暂态的或非暂态的。一个或多个
计算机可读介质可以是例如电子、磁、光、电磁、红外或半导体系统，或者用于数据传输(例如用于通过互联网下载代码)的传播介质。替代地，一个或多个计算机可读介质可以采用一个或多个物理计算机可读介质的形式，诸如半导体存储器或固态存储器、磁带、可移动计算机软盘、随机存取存储器(RAM)、只读存储器(ROM)、硬磁盘和光盘，诸如CD‑ROM、CD‑R/W或DVD。

[0078] 在一种实施方式中，本文描述的模块、组件和其他特征可以被实现为分立组件或集成在诸如ASICS、FPGA、DSP或类似设备的硬件成分的功能中。

[0079] “硬件组件”是能够执行某些操作并且可以以某种物理方式配置或布置的有形(例如，非暂态)物理组件(例如，一个或多个处理器的集合)。硬件组件可以包括被永久配置为执行某些操作的专用电路或逻辑。硬件组件可以是或包括专用处理器，诸如现场可编程门
阵列(FPGA)或ASIC。硬件组件还可以包括由软件临时配置以执行某些操作的可编程逻辑或
电路。

[0080] 因此，短语“硬件组件”应当被理解为涵盖可以被物理构造、永久配置(例如，硬连线)或临时配置(例如，编程)为以某种方式操作或执行本文描述的某些操作的有形实体。

[0081] 另外，模块和组件可以被实现为硬件设备内的固件或功能电路。此外，模块和组件可以以硬件设备和软件组件的任何组合来实现，或者仅以软件(例如，存储或以其他方式体现在机器可读介质或传输介质中的代码)来实现。

[0082] 除非另有明确说明，从以下讨论中可以明显看出，应当理解，在整个描述中，讨论使用诸如“接收”、“确定”、“比较”、“使能”、“维持”、“标识”、“获得”、“接收”、“生成”、“应用”、“调整”、“产生”、“缩放”、“得到”等术语，是指计算机系统或类似电子计算设备的动作和过程，其操纵计算机系统的寄存器和存储器内表示为物理(电子)量的数据并且将其转换为类似地表示为计算机系统存储器或寄存器或其他此类信息存储、传输或显示装置内的物
理量的其他数据。

[0083] 示例

[0084] 为了证明衰老时钟方法，发明人分析了来自Tabula Muris Senis的单细胞基因表达数据(Tabula Muris  Consortium.A single‑cell transcriptomic  atlas 
characterizes ageing tissues in the mouse.Nature583,590–595(2020)。https://
doi.org/10.1038/s41586‑020‑2496‑1，通过引用并入)，该数据包含来自已知实足年龄的小鼠的多个组织的细胞的转录组。从微流体(“液滴”)方法获得的数据包含四个具有足够时间点来尝试时钟训练的组织：心脏、肺、四肢肌肉和脾脏。在这些组织内，选择了最普遍的细胞类型(由原作者注释)，从而得出表1中概述的选择。分析也仅限于雄性小鼠，以防止任何性别影响。

[0085] 表1‑表格概述了Tabula Muris Senis中的组织及其贡献细胞类型，其中有足够的细胞来可靠地训练单细胞衰老时钟。

[0086]

[0087] 中值绝对误差(MAE)被发现对于本文描述的训练过程来说是一个很好的损失函数。图5是显示中值绝对误差(MAE)在衰老时钟的预测年龄中的分布的示例，衰老时钟使用
来自心脏内皮细胞的不同数目的主成分作为输入进行训练。在经较少的主成分训练的时钟
中，最后的已被丢弃。对于每个数目的主成分，训练过程重复10次。虚线是接近最佳模型的阈值，阴影框是可以产生具有这种性能的时钟的主成分的最低数目。通过首先标识产生具
有最低平均值MAE的时钟的主成分的数目来计算阈值。然后将由使用该数目的主成分训练
的10个时钟产生的MAE的标准误差添加到该组时钟的平均值MAE中。

[0088] 时钟测试

[0089] 图6A至G示出了在来自小鼠心脏(A和B)、肺(C和D)、四肢肌肉(E和F)和脾脏(G)的测试集的单细胞中通过本文描述的方法直接在基因表达上进行训练的时钟的性能(在图中
称为“RD时钟”)；个体细胞类型被标记如下。

[0090] 字母 组织:细胞类型A 心脏：内皮细胞
B 心脏：成纤维细胞
C 肺：经典单核细胞
D 肺：支气管平滑肌
E 四肢肌肉：间充质干细胞
F 四肢肌肉：骨骼肌卫星细胞
G 脾脏：B细胞

[0091] 每个箱线图代表来自单个供体小鼠的细胞的预测年龄的分布：分组箱线图对应于相同年龄的小鼠，在x方向上抖动以辅助可视化。上枢纽和下枢纽分别对应于第75个百分位数和第25个百分位数，中间枢纽表示中值。胡须从外部枢纽扩展1.5*四分位数范围，落在该范围之外的点由黑点表示。每个图的每个细胞的中值绝对误差(MAE)以月为单位显示，皮尔逊(Pearson)相关系数(Cor.)也是如此，并且y＝x由黑色虚线表示。

[0092] 当在单个数据集中进行训练和测试时，本文描述的方法的误差(由MAE测量)与直接在前2000个高度可变基因上训练的时钟的误差相似，如图6A至G所示。然而，由于后来的主成分已从本文描述的时钟中被丢弃，因此模型将受到技术噪声的影响较小，并且可能较
少因“过度拟合”而产生偏差。因此，本文描述的时钟的准确性比直接基因表达时钟的准确性更少膨胀。

[0093] 从图7中可以看出，另一个好处是减少了时钟训练所需的时间，并且这种时间减少随着用于训练的细胞的数目的增加而增加。图7显示了执行(直接在基因表达上的
[“Expr.”，方形]或通过本文描述的方法的[“RD”，圆形])时钟训练的单次迭代所需的平均时间与用于训练过程的细胞的数目的关系。这些点根据用于训练的组织和细胞类型进行字
母标注，并且已经使用线性回归拟合直线。插图：来自时钟方法的值使用与主图相同的x轴来被显示，但使用截断的y轴以辅助可视化。使用AMD Ryzen 7 5800X 8核处理器(3.80GHz)和32GB RAM进行训练。

[0094] 鉴于在训练和优化期间训练过程往往需要重复数千次，这种时间减少是显著的。作为参考，“现实”训练集大小约为5000个细胞(脾B细胞)，本文描述的方法大约快60倍。

[0095] 在数据集之间转移时钟

[0096] 衰老时钟方法可以用于预测数据集中单个细胞的供体的年龄，这些数据集很少或没有先前的年龄注释。发明人再次使用Tabula Muris Senis来证明这一点，这是因为它还
包含针对之前使用的四种组织的单细胞表达数据，这些数据是通过基于荧光激活细胞分选
(FACS)的不同测序方法收集的。在3、18和24个月时收集雄性小鼠的细胞。然而，由于在通过液滴法概况的所有组织中缺乏1个月与18个月之间的雄性小鼠样本，因此3个月大的细胞被
排除在进一步分析之外。

[0097] 图8显示了针对本文所述的衰老时钟(“RD”)和针对直接在基因表达上训练的时钟(“Expr.”)的性能度量的比较。在每个面板中，度量被归一化到本文描述的时钟的度量。A：
单个数据集中每次ELN训练迭代所花费的时间；B：在单个数据集上训练和测试时钟时每个
细胞的MAE；C：在时钟在一个数据集中进行训练并且用于预测单独数据集中的细胞的年龄
时每个细胞的MAE。在C中，本文描述的时钟在由针对液滴细胞的MNN方法产生的校正PCA矩
阵上进行训练，并且在FACS细胞上进行测试；将图6A至G中训练的直接表达时钟直接应用于FACS细胞；时钟还在根据MNN校正PCA矩阵重建的表达矩阵(“Expr.recon.”)上被训练。

[0098] 如上所述，将之前在液滴数据上训练的直接表达时钟应用于FACS数据；平均而言，在批次校正后，它们的表现比本文描述的在液滴数据上训练并且在FACS数据上测试的时钟差(图8C)。本文所述时钟的MAE平均降低了37％；在批次效应更显著的数据集中，这种改进也可能更大。

[0099] 当直接表达时钟被批次转移时，对技术噪声的过度拟合也会导致直接表达时钟相对于本文描述的时钟的误差增加。一般来说，过度拟合会降低直接表达时钟的普遍性以及
由此得到的任何生物学结论。这意味着在数据集中训练的直接表达时钟在该数据集的生物
复制中表现不佳，即使在绝对零批次效应(极不可能)的情况下也是如此。由于后者的条件
很难满足，可以研究这些时钟的通用性的一种方法是使用批次校正的输出。由于MNN方法的输出是类似于校正PCA矩阵的校正矩阵，因此可以通过与本文描述的方法类似的方法根据
该矩阵重建“校正”基因表达矩阵。需要注意的是，由于细胞在PCA空间中的强制运动，所得到的“表达”矩阵将高度扭曲，并且通常不应用作针对真实基因表达的数学替代。然而，该矩阵代表唯一实用的方法，通过该方法可以将本文所述的时钟的通用性与不存在批次效应的
直接表达时钟进行比较。在这些条件下，本文描述的时钟方法产生了减少的误差(平均MAE
减少＝30％，图8C)，这表明根据本公开的时钟中没有过度拟合而产生了显著的益处。

[0100] 应当理解，以上描述旨在是说明性的，而不是限制性的。在阅读并理解以上描述后，许多其他实施方式对于本领域技术人员来说将是显而易见的。尽管已经参考特定示例
实施方式描述了本公开，但是应当认识到，本公开不限于所描述的实施方式，而是可以在所附权利要求的精神和范围内通过修改和改变来实践。因此，说明书和附图应当被认为是说
明性的意义而不是限制性的意义。因此，本公开的范围应当参考所附权利要求以及这些权
利要求所享有的等同物的完整范围来确定。

[0101] 公开的方面和实施例包括以下编号的条款：

[0102] 1.一种获得用于基于基因转录数据预测时变特性的预测器的计算机实现方法，该方法包括：

[0103] 接收数据集，数据集包括从具有不同值的时变特性的相应的细胞样本获得的数据样本，每个数据样本包括一定数目的转录水平以及针对每个数据样本的细胞样本的时变特
性的相应的实际值，其中每个转录水平是个体基因转录物的转录水平或个体基因的多个基
因转录物的池化转录水平；

[0104] 生成嵌入数据集，该嵌入数据集包括针对每个数据样本的嵌入样本，其中嵌入样本的维度的数目小于转录水平的数目；

[0105] 将嵌入数据集作为输入应用到预测器，以针对每个嵌入样本产生时变特性的预测值；

[0106] 通过调整预测器的预测系数以减少时变特性的相应的预测值与实际值之间的预测误差的误差测量，获得预测器。

[0107] 2.根据条款1所述的方法，包括：对数据集应用变换以生成嵌入数据集，所述方法还包括：通过对数据集进行操作来获得变换。

[0108] 3.根据条款2所述的方法，包括：不使用基因通路的知识来获得变换。

[0109] 4.根据条款2或3所述的方法，包括：通过对数据集的协方差矩阵进行操作来获得变换。

[0110] 5.根据任一前述条款所述的方法，包括：缩放嵌入数据集以具有跨维度基本恒定的方差。

[0111] 6.根据任一前述条款所述的方法，包括：对转录数据集应用线性变换以生成嵌入数据集。

[0112] 7.根据条款6所述的方法，其中嵌入数据集包括转录数据集的主成分的子集。

[0113] 8.根据任一前述条款的方法，包括：对预测系数应用从嵌入数据样本映射到数据样本的逆映射，以将预测系数投影到数据集的维度上，从而得到对预测针对每个基因或转
录物的时变特性的值的贡献的测量。

[0114] 9.根据条款8所述的方法，当从属于权利要求6和7时，其中逆映射包括线性变换的矩阵求逆。

[0115] 10.根据任一前述条款所述的方法，其中预测器是线性预测器。

[0116] 11.根据条款10所述的方法，其中线性预测器包括正则化方法以促进预测系数的稀疏性。

[0117] 12.根据任一前述条款所述的方法，其中调整预测系数包括弹性网络回归。

[0118] 13.根据任一前述条款所述的方法，其中预测误差是中值绝对预测误差。

[0119] 14.根据任一前述条款所述的方法，还包括：

[0120] 接收另外的数据集，该另外的数据集包括从具有时变特性的不同值的相应的另外的细胞样本获得的另外的数据样本，每个另外的数据样本包括一定数目的另外的转录水平
以及针对每个另外的数据样本的另外的细胞样本的时变特性的相应的另外的实际值，其中
每个另外的转录水平是个体基因转录物的转录水平或个体基因的多个基因转录物的池化
转录水平；

[0121] 将数据集和另外的数据集变换成包括数据样本和另外的数据样本的公共数据集，从而减少数据样本和另外的数据样本对于数据集和另外的数据集非公共的变异性，以及

[0122] 其中生成嵌入数据集包括：针对公共数据集中的每个数据样本生成嵌入样本。

[0123] 15.根据条款1至14中任一项所述的方法，还包括：

[0124] 接收另外的数据集，该另外的数据集包括从具有不同值的时变特性的相应的另外的细胞样本获得的另外的数据样本，每个另外的数据样本包括一定数目的另外的转录水
平，其中每个另外的转录水平是个体基因转录物的转录水平或个体基因的多个基因转录物
的池化转录水平；

[0125] 将数据集和另外的数据集变换成包括数据样本和另外的数据样本的公共数据集，从而减少数据样本和另外的数据样本对于数据集和另外的数据集非公共的变异性，

[0126] 其中生成嵌入数据集包括：针对公共数据集中的每个数据样本生成嵌入样本，并且

[0127] 其中将嵌入数据集作为输入应用到预测器包括：仅将与基因转录数据样本相对应的嵌入样本作为输入进行应用，以产生针对与数据样本相对应的嵌入数据样本的时变特性
的相应的预测值；以及

[0128] 在获得预测器之后，将与另外的数据样本相对应的嵌入样本应用到预测器，以预测针对另外的细胞样本的时变特性的相应的值。

[0129] 16.根据任一前述条款所述的方法，其中嵌入样本的维度的数目是基于具有不同的相应数目的维度的嵌入数据集相应的预测性能来选择的。

[0130] 17.根据任一前述条款所述的方法，其中时变特性是关于已经从其获得细胞样本的一个或多个生物体或受试方的。

[0131] 18.根据条款16所述的方法，还包括：生成标识一个或多个生物体或受试方的时变特性的报告。

[0132] 19.根据条款1至16中任一项所述的方法，其中细胞是已经从细胞培养物中获得的。

[0133] 20.根据任一前述条款所述的方法，其中细胞样本均为单个细胞的单细胞样本。

[0134] 21.根据任一前述条款所述的方法，其中时变特性是生物学年龄。

[0135] 22.根据条款1至20中任一项所述的方法，其中时变特性是实足年龄。

[0136] 23.根据条款1至20中任一项所述的方法，其中时变特性是病症或疾病的进展状态。

[0137] 24.根据条款23所述的方法，其中病症或疾病是神经退行性疾病。

[0138] 25.根据条款24所述的方法，其中神经退行性是阿尔茨海默病。

[0139] 26.根据条款24所述的方法，其中神经退行性是帕金森病。

[0140] 27.根据条款1至20中任一项所述的方法，其中时变特性是癌症的进展状态。

[0141] 28.根据任一前述条款所述的方法，其中转录水平以及在存在情况下的另外的转录水平，在不使用基因通路的知识的情况下已经由细胞样本中基因转录物的转录计数得
到。

[0142] 29.根据条款28所述的方法，包括：不使用基因通路的知识的情况下，生成嵌入数据集。

[0143] 30.根据条款28所述的方法，包括：在不使用基因通路的知识的情况下，应用嵌入数据集并且获得预测器。

[0144] 31.根据任一前述条款所述的方法，还包括：在获得预测器之后，使用基因通路或生物活性的现有知识、任何其他现有生物学知识或从预测系数得到的知识来优化预测器。

[0145] 32.一种计算机程序产品，包括计算机代码指令，该计算机代码指令当在处理器上执行时，实现根据任一前述条款所述的方法。

[0146] 33.一种计算机可读介质，包括计算机代码指令，该计算机代码指令当在处理器上执行时，实现根据条款1至31中任一项所述的方法。

[0147] 34.一种系统，包括处理器和如条款33中所定义的计算机可读介质，其中处理器被配置为执行计算机代码指令。

[0148] 35.一种系统，包括用于实现如在条款1至31中任一项所定义的方法的部件。