多变量表型分析能够对哺乳动物基因功能进行全基因组推断
乔治·尼科尔森 ，休·摩根，哈比卜·甘吉加希，史蒂夫?布朗，安-玛丽?马龙，克里斯·福尔摩斯
出版日期： 2022年08月09日
抽象
人类和小鼠基因组中大多数基因的功能尚不清楚。研究和阐明这个黑暗的基因组是生物医学科学面临的一个重大挑战。国际小鼠表型分析联盟（IMPC）正在通过为每个蛋白质编码基因生成和广泛基调（KO）小鼠品系来解决这个问题，从而产生一个多维数据集，该数据集是从基因到表型的全基因组注释图谱的基础。在这里，我们开发了一种多变量（MV）统计方法，并将其应用于包含4，548 KO系中测量的148种表型的IMPC数据。
单变量（UV）模型分别分析每个表型的命中率为4，256（302，997个观测数据测量值的1.4%），而MV模型的观测数据结果中有31，843个（10.5%）命中，对应于MV模型的功率相对于UV模型估计增加了7.5倍。数据集的一个关键属性是其55.0%的缺失率，这是由于质量控制过滤器和一些KO线的不完全测量造成的。这就提出了一个问题，即是否有可能在没有数据可用的表型基因对上推断扰动，即使用统计分析而不是实验来推断一些体内效应。我们证明，即使在缺失表型的情况下，MV模型也可以检测扰动，其功率可与单表型分析相媲美，从而以良好的灵敏度填充完整的基因 - 表型图。
对 MV 模型的拟合协方差结构进行因子分析可识别 20 个表型聚类，每个聚类倾向于集体扰动。这些因素累积解释了KO诱导的数据变异的75%，并有助于扰动的生物学解释。我们还证明，MV方法加强了IMPC表型与现有基因注释数据库之间的对应关系。对跨多个实验室重复测量的KO线子集的分析证实，中压模型以高可复制性提高了功率。
引文： Nicholson G，Morgan H，Ganjgahi H，Brown SDM，Mallon A-M，Holmes C（2022）多变量表型分析能够对哺乳动物基因功能进行全基因组推断。PLoS Biol 20（8）：e3001723。https://doi.org/10.1371/journal.pbio.3001723
学术编辑： Nicole Soranzo，Wellcome Trust Sanger Institute，英国
收到： 2020年1月13日;接受： 六月 22， 2022;发表： 八月 9， 2022
版权所有： ? 2022 Nicholson et al.这是一篇根据知识共享署名许可条款分发的开放获取文章，允许在任何媒体上不受限制地使用，分发和复制，前提是注明原作者和来源。
数据可用性： 作者证实，研究结果背后的所有数据都是完全可用的，不受限制。用于生成论文结果的数据和代码可在 https://github.com/georgenicholson/multivariate_phenotype_data_and_code 和 https://zenodo.org/record/6787112 获得。
资金： 这项工作得到了医学研究委员会（https://mrc.ukri.org/）计划赠款MC_UP_A390_1107（G.N.，H.G.和C.H.）和美国国立卫生研究院（https://www.nih.gov/）赠款U54 HG006370（H.M.，S.D.M.B.，A.-M.M.）的支持。资助者在研究设计，数据收集和分析，出版决定或手稿准备方面没有任何作用。
竞争利益： 作者宣布不存在相互竞争的利益。
缩写： 血压， 生物过程;EM，期望最大化;eQTL，表达数量性状位点;Fdr，错误发现率;Fsr，错误符号率;GO，基因本体论;IMPC，国际小鼠表型分析联盟;KO，淘汰赛;LOO-MV，一个程序出MV;MAP，最大后验;MAR，随机失踪;麦芽浆，多元自适应收缩;MCMC，马尔可夫连锁店蒙特卡洛;MV，多变量;或者，比值比;紫外线，单变量;WT，野生型;XD，极端反卷积
介绍
人类和小鼠基因组中大多数基因的功能尚不清楚。研究和阐明这个黑暗的基因组是生物医学科学的主要挑战[1]。制定哺乳动物基因功能的综合目录将是罕见和常见疾病研究和精准医学进展的重要基础。国际小鼠表型分析联盟（IMPC）是全球21个研究机构之间的合作，旨在通过为每个蛋白质编码基因（www.mousephenotype.org）生成和广泛基调敲（KO）小鼠系来应对暗基因组的挑战。
对每只动物进行了300多次测量，范围从临床血液化学、量热法和身体成分分析，到行为表型[2]。通过从生成的多维数据集中推断，IMPC正在编制从基因到表型的全基因组注释图，该图谱已经为哺乳动物基因功能和各种疾病的基因组景观提供了独特的见解[3，4，5，6，7，8]。
到2022年3月，到目前为止，已经产生了大约10，000个KO小鼠品系，其中许多是针对鲜为人知的基因，其中8，623个品系已经使用标准化程序对各种疾病系统进行了表型分析。在本文中，我们分析了部分IMPC数据集，该数据集包含截至2018年3月26日的148个定量表型中的一些表型数据，其中包含4，548个KO品系。
在IMPC成人表型分析管道中，对单基因KO和年龄在9至16周之间的对照小鼠进行一系列标准化测量。我们将测量称为表型，这些表型以称为程序的组进行测量;在给定的一周内测量给定程序中的所有表型（S1图）。科学目的，实验设计和每个程序的详细描述都在IMPC网站上提供[9]。主要科学目标是识别具有统计学意义的KO诱导的表型扰动，也称为表型命中或阳性注释。IMPC的实验设计平均测量每个KO线的14只动物（每个性别7只），与对照动物的滚动基线同时进行。这在图1中其中一个表型中心MRC Harwell的一对表型中可视化。统计目标是估计和测试每个KO线和共享的对照动物集之间的表型平均值的差异。从概念上讲，KO和控件之间的不成对t检验是基本的统计思想，但在实践中，由于复杂的实验结构，多级建模是必要的。例如，凋落物和其他实验地层偶尔会与感兴趣的基因-KO效应混淆，因此需要使用分层模型来识别感兴趣的效应[10，11，12]。
缩略图 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 1. IMPC的实验设计。
每个点对应于 1 只动物，来自 2 条 KO 线的数据（标记为 g 和）与来自大量对照（野生型或 WT）动物的同期数据一起显示为灰色（见图例）。图谱（a）和（b）分别显示表型p（甘油三酯）和（体脂百分比）的数据。我们的目标是在存在结构化实验噪声的情况下，量化滚动WT基线（用平滑的黑色曲线表示）的红色/蓝色点的潜在预期扰动。在图上注释的每个基因表型对的红色/蓝色测量数据的右侧是UV和MV模型的后平均估计值，带有空正方形和填充正方形（见图例），以及表示±2个后部SD的误差线。在目前的论文中，我们将多个相关表型中的信息组合在一起，例如p和，从而生成改进的估计器 。图中所示的UV和MV估计器的相对均值和SD说明了它们的一般特性 - 在90.2%的情况下（表型 - 基因对）相对于UV后均值，MV后均值向零（此处为黑色曲线）缩小;而>99.9%的病例中，MV后部SDs小于UV后部SDs。用于生成此数字的数据和代码可在 [13，14] 中找到。IMPC，国际小鼠表型分析联盟;KO，淘汰赛;MV，多变量;紫外线，单变量;WT，野生型。
https://doi.org/10.1371/journal.pbio.3001723.g001
到目前为止，关于高通量表型的文献完全集中在使用所谓的单变量（UV）模型一次调用一种表型的命中（测试扰动）[11，12]。然而，IMPC的初步结果显示，不同表型的扰动之间存在很强的相关性。多变量（MV）关联方法已在许多遗传应用中被证明是成功的，例如全基因组关联研究[15，16，17]和多组织eQTL研究[18，19，20，21，22]。这指出了通过使用MV方法跨表型共享信息来改善IMPC中推理的机会。特别是，当样本量因道德和财务原因而受到严重限制时，希望MV方法可以在计算上增加从收集的数据中提取的信息。此外，在我们的IMPC数据集中，并非所有表型在每个KO系上都可用。这就提出了一个问题，即是否有可能在没有数据可用的（表型，KO线）对上推断扰动，即使用统计分析而不是实验来推断一些体内效应。我们着手实现一个MV模型，当某些数据丢失时，它可以有效地执行这种类型的推理。
我们采用可组合的MV建模方法，该方法在计算上具有吸引力，同时有效地捕获了IMPC数据集中的重要变化。首先，我们为每个表型分别拟合一个UV多级模型[11]。其次，我们采用UV模型的输出效应估计值和标准误差，并将MV模型拟合到这些，在方法论上建立在[22，23]的工作之上。我们在方法 - 与现有方法的比较中，在现有工作的背景下，将我们方法的性能进行情境化和比较。
IMPC的一个主要目标是创建一个全面的基因 - 表型注释图。从统计学的角度来看，这涉及检验没有表型扰动的原假设。除了MV模型之外，我们还设计了一种基于排列的假设检验方法，旨在在仔细控制和监测假阳性率的情况下进行强大的推理。我们的方法基于通过对对照动物进行结构化随机重采样来生成合成零KO线（详细信息请参阅方法 - 错误率控制）。通过分析合成零线和真KO线，我们能够选择有效控制错误率的显著性阈值。
我们通过多种方式验证我们的 MV 方法。我们通过人为地屏蔽数据并将屏蔽的数据结果与完全观察到的数据结果进行比较来评估在存在缺失数据的情况下推理的有效性。我们通过检查在多个实验室测量的相同KO线上调用的命中的可复制性，独立评估MV命中调用方法。我们还围绕结果的生物学合理性进行了一些额外的检查，并评估了模型鲁棒性和拟合性的定量测量。我们的检查表明，MV方法可以大大提高IMPC中的命中率，同时保持错误率控制和可复制性，即使在缺少表型数据的情况下调用命中时也是如此。开发灵敏、可复制和全面的基因表型图谱将确保将IMPC后续实验中使用的动物数量降至最低，这与置换、还原和细化的3R相一致[24]。
结果
我们之前设计了一个UV贝叶斯线性多级模型，针对基因KO动物相对于野生型（WT）动物的表型扰动[11]。我们将该模型分别拟合到每个（表型，中心）组合，从而产生表型扰动的估计值（和SE），对于每个（表型p = 1,...,P，基因g = 1,...,G）对，其中有可用的测量值。示例数据和估计值如图1所示。
在本文中，我们开发了一个MV建模框架，基于[22，23]的方法学工作。该方法以所有（表型p，基因g）组合的UV结果（）作为输入，并输出MV估计值（），包括那些无法获得数据的对。MV模型基于协方差结构Σ，允许扰动在不同表型之间相关，如图1所示，甘油三酯和体脂百分比之间所示。该方法还包含一个相关矩阵R，以解释跨表型的实验噪声中的结构。这种2阶段模型的一个实际有用的特性是它的可组合性，因此结果可以在2种不同的分析或计算工具之间有效地转移 - 在这里从任意复杂的UV模型到高度结构化的MV模型。
我们将结果分为3个概念阶段。首先，我们提供UV和MV模型的高级技术描述。其次，我们表征了 IMPC 命中调用结果，对比了 UV 和 MV 模型，并侧重于展示统计功效和可复制性。最后，我们研究应用来演示MV方法如何阐明表型扰动与潜在生物学机制之间的关系，并且与UV相对有效。这些示例通过说明其结果如何具有直观意义并与现有科学知识保持一致，为MV方法的有效性和可复制性提供了额外的证据。
单变量模型
整个感兴趣的参数由θ表示页码并表示相对于WT动物而言，第g个基因KO中第p个表型的预期扰动（图1）。该UV模型仅拟合来自KO线g的数据以及来自WT动物的整个滚动基线的数据，采用线性多级模型（或混合效应模型）的形式：
(1)
其中 y我是第 i个小鼠上第 p个表型的 Box–Cox 变换 [25] 测量值。β中的参数根据性别，性别 - 基因型相互作用，菌株，研究者和其他实验元数据进行加性调整，而日和凋落物效应则通过α日和α垃圾使用方差分量和 .在本文中，我们专注于θ的估计。页码，基因型g对表型p的主要影响。在基因型效应因性别而异的情况下[5]，θ页码可解释为这些性别特异性效应的平均值。测量基线的纵向变化使用惩罚性样条进行建模，该样条具有固定和随机分量[26]。为 θ 指定非信息先验页码、β和σr，模型通过马尔可夫链蒙特卡洛（MCMC）拟合，并从边际后验分布p（θ）输出样本页码|y）（有关更多详细信息，请参阅 S1 注释和 [11]）。
UV 推断输出每个 θ 的估计值和标准误差页码，即后均值和后部 SD 分别为。我们对紫外线结果进行仔细的质量控制，保守地过滤掉（从下游MV分析）任何在紫外线结果中表现出异常纵向模式的中心程序组合;这种模式可以表明未建模的实验文物，而不是生物效应（S2图）。接下来，为了确保每种表型都有足够的数据，我们应用了QC后启发式过滤器，从而仅保留至少500 KO系的具有UV效应估计值的表型。在QC和过滤之后，对UV估计值（和SE）进行缩放，以便每个表型中心内每种表型的单位SD，然后作为MV模型的输入向前推进。
多元模型
在收集UV多级模型的结果时，我们获得了θ的无偏估计值（和SE）·g受中压实验噪声影响，具有协方差结构。此外，潜在的P维MV扰动θ·g趋向于表现出较强的P×P协方差结构。数据的这些方面表明模型遵循[22，23]的形式：
(2)
(3)
其中参数 Σs表示 θ 的协方差·g，即KO线的预期表型扰动，超参数R模拟实验噪声的相关结构。这些是 UV 模型输出的标准误差的已知对角矩阵。潜伏扰动的密度 p（θ·g|Σ， π1：米），是混合概率π的中压高斯混合模型1：M;1：S在由常数 ω 给出的指定阶梯上1：米 [22，27]和 S≥1 协方差矩阵 Σ1：S待学习[23]。我们将我们的方法更详细地与[22，23]联系起来，并在方法中比较性能 - 与现有方法进行比较。
我们约束 Σs到因子模型形式（例如，参见 [28]）：
(4)
其中 Ws是一个 P×K 矩阵，而 Ψs具有正对角线元素的对角线 P×P 矩阵。我们对固定 K∈{15，20，30，40} 进行了全面分析。手稿中提出的结果是针对K = 20和S = 1的选择，参考错误发现率控制（Fdr控制）命中率[29，30]和交叉验证的模型拟合似然度量（参见方法-与现有方法的比较）进行选择。
我们采用经验贝叶斯方法在（2）–（4）指定的MV模型中进行推理。实验相关超参数R是根据合成零数据估计的，并预先固定在[22]。期望-最大化 （EM） 算法用于获得超参数 Σ 的最大后验 （MAP） 估计值1：S和π在平坦先验下（EM算法的推导和更多细节在S2 Note中）。以 MAP 估计为条件，θ 的后验·g以封闭形式提供（当数据缺失时，请参阅方法–MV 模型）。
结果的可视化概述
为了在 UV 和 MV 模型之间进行全局比较，我们通过定义为 的 z 统计量可视化 UV 和 MV 分析的输出。为了增强解释，z 统计量按其相应的显著性阈值进行缩放，即我们绘制 ，以便对应于显著扰动（图 2）。MV 模型中显著 z 统计量的比例更大，在同一过程中，显著性通常跨表型同时发生，并且方向在程序内趋于相关。在 UV 模型热图中，由白色区域表示的缺失数据的实例。
缩略图 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 2. MV 模型灵敏度提高的全局表示。
每行对应于一个表型，多个表型按程序分组，左标记。每列对应于一个KO线，多个线按表型实验室分组，底部标记。为了有效地可视化，仅显示500 KO行的随机子集。热图显示缩放的z统计量，因此和对应于导致表型中显着增加/减少的基因KO（方法 - 错误率控制）。（a） 紫外线模型，其中白色方块表示缺少表型数据;（b） 中压模型。用于生成此数字的数据和代码可在 [13，14] 中找到。KO，淘汰赛;MV，多变量。
https://doi.org/10.1371/journal.pbio.3001723.g002
我们继续展示对应于 和 的相关性估计值的热图（图 3）。在几个程序中存在明显的相关性块，这是预期的，因为类似的表型倾向于聚集在程序中。虽然几乎所有的实验相关性都发生在同一程序中的表型之间，但其中有大量的块状对角线相关结构，表明不同程序之间存在相关的表型扰动。例如，在图3A中，KO基因扰动与开放场和明暗测试，临床化学和身体组成，听觉脑干反应和声学惊吓，身体组成和回声以及血液学和临床化学之间相关。
缩略图 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 3. 相关矩阵的热图，这些矩阵是 IMPC 中表型的系统性共扰动的基础。
（a） 基因KO诱导的生物协变的估计相关矩阵， .（b） 估计的实验相关矩阵，可归因于测量过程而不是目标生物学。用于生成此数字的数据和代码可在 [13，14] 中找到。IMPC，国际小鼠表型分析联盟;KO，淘汰赛。
https://doi.org/10.1371/journal.pbio.3001723.g003
能够检测 KO 扰动
我们首先比较MV和UV模型的统计功效，以检测扰动，即在观察到数据的基因 - 表型对上，以及UV和MV结果都可用的扰动。我们使用基于专门创建的合成零数据的Westfall-Young排列过程[31，32]将Fdr控制在5%以下;详细信息请参阅方法 – 错误率控制。
图4直观地表示了各种方法的相对功率和重叠，包括与现有IMPC数据库的比较，该数据库包含来自不同UV方法的结果[12]-我们将在结果-与IMPC数据库的比较中更详细地讨论这种比较。UV模型调用了4，256次（302，997次观测数据测量值的1.4%）命中，而MV模型的观测数据结果中有31，843次（10.5%）命中，对应于相对于UV模型，MV模型的功率估计增加了7.5倍。当我们检查UV和MV结果之间的一致性时，只有UV模型调用了95（0.0%）命中，而MV模型调用的命中率为27，682（9.1%）。在 4，256 个 UV 模型表型命中中，MV 模型在 4，161 个 （97.8%） 案例中共同调用具有相同方向性的命中，并且从不调用与 UV 模型相反方向的命中。另请参阅 S3 图，了解比较 UV 和 MV 模型输出的缩放 z 统计量的散点图。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 4. 直观比较这些方法的命中率和重叠，以观察和缺失的测量值。
黑色和红色的大轮廓圆圈表示按面积划分的观测和缺失的测量值的数量。对应于一个方法的每个圆圈都有表示调用的命中数（在观察到的或缺失的数据上）的区域。圆圈之间的重叠区域表示两种方法调用的命中数。用于生成此数字的数据和代码可在 [13，14] 中找到。IMPC，国际小鼠表型分析联盟;MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.g004
我们继续更详细地检查MV和UV模型的相对灵敏度，通过比较每个表型（图5A）和KO基因（图5B）的2个模型调用的命中次数。MV模型在所有148种表型和2，750种（占60.5%）KO系中发现了比UV模型更多的扰动;UV模型在33条（0.7%）KO线中识别出更多的扰动。平均而言，与UV模型相比，MV模型每个表型的点击次数增加了186.4次，每KO线的点击次数增加了6.1次。图5C检查了MV方法的程序性功率增强，显示了每个程序中至少有一个显着扰动表型的KO线的比例（另见S4图的表型比较，以及每种表型缺失数据比例的详细信息）。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 5. 功率增强：中压方法提高了检测基因KO诱导的扰动的灵敏度。
（a） 中压模型与紫外线模型确定的每种表型的扰动次数。（b） 中压与紫外线型号确定的每条KO线的扰动数。（c） 用于比较UV方法和MV方法（在测量和缺失数据上）在手术中至少击中1次的线的比例（即在该过程中至少有1个表型受到干扰）。程序按UV模型的命中率排序。用于生成此数字的数据和代码可在 [13，14] 中找到。KO，淘汰赛;MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.g005
数据缺失时的推理
即使对于没有测量数据的基因 - 表型对（此处称为缺失数据），MV模型也可用于通过未测量和测量的表型之间存在的相关结构来推断基因KO效应。MV模型识别出4，819个扰动（370，107个缺失数据案例中的1.3%），这与UV模型在观测数据上的命中率为1.4%相比是有利的。当缺失数据结果与观测到的数据结果相结合时，MV模型总共检测到36，662个扰动，比UV方法检测到的4，256个扰动增加了8.6倍。
重要的是要注意，θ的估计·g当仅部分观测到时，只要统计模型相对于基础数据生成机制进行了很好的指定，并且未观测到的数据是随机缺失的（MAR）[23，33，34]，则可以连贯地执行。虽然有很大比例的缺失数据，但从图2中可以清楚地看出，大部分数据在明显的块中缺失，并且是某些中心未执行某些测量/程序的结果。在某些中心系统地不进行一部分测量的背景下，MAR假设是合理的，因为缺失数据机制，“给定缺失数据和观测数据的值，对于缺失数据的所有可能值都是相同的。[33].
尽管有这种令人放心的观察结果，但在如此庞大而复杂的数据集中，自然仍然有一些相对较小比例的数据违反了MAR假设。因此，我们对MAR假设的实际合理性进行额外检查。这些在结果 - 验证可复制性（参考图6C和6D）和方法 - 预测屏蔽数据中进行了描述。我们对从业者的建议是，根据Rubin及其同事的工作，仔细检查MAR假设在其特定背景下的适当性[33，34]。如果对MAR假设的有效性有任何疑问，我们建议进一步进行实证检查。特别是，方法 - 预测掩码数据中描述的交叉验证掩码和预测方法可以在具有缺失数据的各种MV数据集中实现，我们建议将其作为在缺失率高时事后检查准确性的工具。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 6. 可复制性验证散点图，用于比较表型实验室的结果。
每个面板检查一对复制结果的不同类型的比较：（a）UV模型与UV模型，（b）MV模型（测量）与MV模型（测量），（c）MV模型（缺失）与MV模型（缺失），（d）MV模型（测量）与MV模型（缺失）。我们通过散射拼分尺度z统计来检查KO参考线的实验室间一致性，对于相同的KO线，但在不同的实验室中测量。显著扰动对应于 ，如在图上用红色虚线分隔的那样。图中的每个点对应于2个不同的实验室，在同一KO线上测量相同的表型。估计错误符号率 （Fsr） 的最有用情况发生在两个实验室都检测到显著扰动时，这些扰动对应于散点图上蓝色/红色阴影区域中的点。实验室在蓝色阴影区域中意见一致，但在红色阴影区域中意见一致。 估计值（95% CI）显示在每个面板的顶部，并基于在阴影区域中观察到的一致/不同意水平（方法 - 可复制性和错误符号率）。叠加每个显著性组合的计数 （%）;当轴扩展到 [–3， 3] 时，计数适用于所有数据，包括超出图比例的数据。用于生成此数字的数据和代码可在 [13，14] 中找到。Fsr，错误符号率;KO，淘汰赛;MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.g006
验证可复制性
我们通过利用实验数据的多实验室性质来验证UV和MV结果。作为IMPC的一部分，在多个实验室中多次测量了少量KO系，对其特殊状态视而不见，即同一基因KO，在多个实验室中表型;我们将这些称为参考线。我们在UV和MV模型下分析它们，同时确保模型对它们作为复制样品之间的对应关系视而不见。分析后，我们揭示了参考线，并检查了同一参考线上的多个表型中心的结果的可复制性。
S5 图绘制了 UV 和 MV 模型下参考线的注释结果。这在视觉上强化了MV模型的影响：它强烈地提高了命中率（用更高的交叉密度表示），并以可复制的方式做到这一点。MV模型的成对参考线命中的方向性在295种情况下是一致的，在7种情况下是不一致的。观察到的可复制性水平可以用相应的错误符号率（Fsr）进行有用的解释，方法 - 可复制性和错误符号率中进行了描述，并使用IMPC参考线进行估计。对于 MV 模型，我们获得了 1.2% （95% 置信区间：0.6% 至 2.4%） 的低全球估计值，这令人放心。
图6进一步了解了参考品系重复中实验室之间的可复制性程度。每个面板中的蓝色/红色阴影区域包含结果在实验室之间分别一致/不一致的实例。MV模型（图6B，6C和6D）比UV模型（图6A）识别出更多的扰动，并且具有高水平的可复制性，这是通过红色阴影区域中的少量点来测量的，并通过每个面板顶部显示的Fsr估计（）进行量化，对应于该面板中显示的结果。重要的是，无论数据是缺失（图6C）还是测量（图6B），或者是在一个实验室测量的，而不是在另一个实验室（图6D）中，实验室之间的一致性水平都是好的。
与 IMPC 数据库的比较
我们将UV和MV模型的有符号表型调用与IMPC数据库中基于不同UV方法的现有调用进行比较[12]。IMPC数据库相关子集中的命中率为1.9%，而我们的UV模型命中率为1.4%，测量数据的中压模型命中率为10.5%。由于错误率控制（标称p<10的差异，与IMPC数据库命中率进行直接比较并不简单?4在 IMPC 数据库中，与 Fdr <我们的 UV 和 MV 型号为 5%）。但是，当我们检查方法与现有数据库的一致性时，我们会看到良好的一致性（表1），指出在所有情况下都有有效的错误率控制。我们的UV模型与IMPC数据库一致，在所有情况下，两者都称为显着的表型命中。我们的 MV 模型仅在 3 种情况下与 IMPC 数据库不一致（在 0.1% 的实例中，它们都调用了命中）。我们在方法中更详细地研究了这些分歧，我们得出结论，在这3种情况下，任何一个模型的表现都优于另一个模型， 几乎没有证据。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 1. 签名命中与现有 IMPC 数据库的比较。
（a） 紫外线模型;（b） 中压模型。将每个模型与现有 IMPC 数据库中调用的相应命中数进行比较（上图）。我们用 {?1，0，1} 中的数字表示呼叫，其中 1 和 ?1 分别表示显著的正表型和负表型扰动，0 表示缺乏统计显著性。
https://doi.org/10.1371/journal.pbio.3001723.t001
杂合子与纯合子
对于IMPC中的一些基因，同时测量杂合子和纯合子KO系。从生物学上讲，杂合子和纯合子表型扰动（如果存在的话）可能朝同一方向起作用，这在生物学上是合理的。因此，我们可以比较杂合子/纯合子对在结果中的一致性（S6图）。正如生物学上所预期的那样，MV模型（7.6%）比相应的杂合子系（2.3%）更频繁地将纯合子系称为命中。在基因的杂合子和纯合子系都被称为命中的情况下，我们在594种情况下观察到方向一致性，在46种情况下观察到不一致。假设所有杂合子/纯合子对确实在同一方向上扰乱表型，那么这种不一致程度与3.7%（95%CI：2.8%至5.0%）一致。这种低 Fsr 估计值进一步证明了我们对命中呼叫中假阳性率的控制是有效的，增加了参考线可复制性分析提供的证据。实际上，可能存在例外，杂合子和纯合子实际上在不同方向上扰动表型，在这种情况下，这种基于合子的估计仍然可以被有效地解释为实际Fsr的上限。复制.
基因本体共富集
基因本体论（GO）使用有向图来注释和相互关联生物学上有意义的GO术语[35，36]，例如声音的感觉感知（GO：0007605）和运动行为（GO：0007626）。每个GO术语都有自己的基因集，通过手动管理已发表的实验文献或通过自动计算方法分配的小鼠基因列表。类似地，我们对IMPC数据库的分析为每个IMPC表型生成一组引起显着表型扰动的基因;我们说每个IMPC表型都有自己的基因集。
通过鉴定具有重叠基因集的GO术语和IMPC表型，我们的目标是增加对影响表型的基因的一般生物学特征的理解。此外，观察GO项和IMPC表型之间的共富集增加了统计方法表现良好的证据（假设检测共富集的假阳性率得到适当控制）。因此，在本节中，我们提出一个问题：哪些对GO术语基因集和IMPC基因集共享的基因集比偶然预期的要多？图7中提供了这种类型的共富集分析的一个例子，其中我们量化了GO：运动行为和IMPC：Locomotor活性的基因集之间的重叠。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 7. 用于共富集测试的 2 x 2 列联表说明性表。
（a） 紫外线模型;（b） 中压模式。每个列联表根据每个基因是否在GO术语基因集（左）和/或IMPC表型基因集（上）中，将每个基因分配到4个类别之一。Fisher 精确的检验 p 值显示在每个表的上方。GO，基因本体论;IMPC，国际小鼠表型分析联盟;MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.g007
我们测试了148个IMPC基因集和生物过程（BP）子本体中的5，368个GO项之间的共富集，这些项被注释为一个或多个IMPC KO基因。我们专注于表现出大扰动的基因（≥2个群体SD），并将每个表型的家庭错误率控制在所有BP GO术语的测试中低于5%;有关更多详细信息，请参阅方法-基因本体分析。在所有基因集中，MV模型识别出1，359对IMPC和GO基因集之间的共富集，而UV模型显示342对的共富集。在80个（54%）IMPC基因集上，MV模型比UV模型识别出更多的共富集GO基因集，而UV模型在仅17个（11%）IMPC基因集上识别出更多的共富集GO基因集（表2提供了更详细的比较）。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 2. 在UV和MV模型中比较了共富集计数。
（a） 与每种 IMPC 表型共同富集的 GO 术语数;例如，有 23 种表型在 UV 模型中有 0 个 GO 项，但在 MV 模型中有 1 到 5 个 GO 项。（b） 与每个GO术语共同富集的IMPC表型的数目;例如，有 180 个 GO 项在 UV 模型中富集了 0 个表型，但在 MV 模型中富集了 1 到 5 个表型。
https://doi.org/10.1371/journal.pbio.3001723.t002
图8提供了IMPC表型和生物学途径之间系统关系的总体表征。与相应的UV模型图（S7图）的比较表明，MV模型产生的GO注释数量大大增加，同时也证明了在两个模型都显示基因富集的情况下，UV和MV结果之间的定性一致性。在图 8 中，行和列在共富集模式的基础上聚类。这种聚类是在不参考按程序对表型进行分组的情况下进行的，因此值得注意的是，来自同一过程的表型倾向于在水平方向上聚类在一起（表型标签根据程序着色 - 参见图8左下角的图例）。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 8. GO项（左）与MPC表型（底部）的共富集，用于MV模型调用的命中。
GO项和IMPC表型之间具有统计学意义的共富集用粗体轮廓正方形表示（控制每种表型的家庭错误率<5%）。正方形的颜色表示GO项处明显扰动的KO基因的百分比，这些基因在正方向上改变了表型（参见顶部的比例尺）。IMPC表型由沿水平轴的GO项模式聚类，而BP GO项按表型模式垂直聚类。表型标签根据左下角图例中的程序进行着色。GO项的子集，由右侧的行（a-h）标记，在图9中更详细地检查。为了便于阅读，我们只在图中包括那些至少具有3个显着共富集实例的IMPC表型和GO项。用于生成此数字的数据和代码可在 [13，14] 中找到。血压， 生物过程;GO，基因本体论;IMPC，国际小鼠表型分析联盟;KO，淘汰赛;MV，多变量。
https://doi.org/10.1371/journal.pbio.3001723.g008
图8中的全球图景是共富集分析与现有科学知识之间的一致性之一。为了说明这一点，我们现在详细检查图8中带有代表性的几行，标记为（a-h），在图9的子表中更详细地呈现。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 9. GO和IMPC基因集共富集的例子。
每个表都列出了GO项（标记为顶部）和IMPC表型（左栏）之间显着共富集的实例，以及Fisher精确检验p值量化共富集证据（右栏）。GO，基因本体论;IMPC，国际小鼠表型分析联盟。
https://doi.org/10.1371/journal.pbio.3001723.g009
GO：脂质生物合成过程的调节与IMPC：总胆固醇，IMPC：HDL-胆固醇和IMPC：甘油三酯表型共同富集，来自IMPC：临床化学程序（图9A）。
GO：运动行为与来自IMPC：Open Field程序的表型共同丰富，用于评估焦虑和探索行为;IMPC：骨区也显示与GO：运动行为重叠，与导致运动受损的异常骨骼结构相容（图9B）。
GO：大脑发育与IMPC：身体成分（DEXA瘦肉/脂肪）（DEXA）程序的代谢表型以及IMPC：胰岛素（图9C）显示出有趣的重叠。
GO：循环系统开发与IMPC：心电图（ECG）程序中的心脏功能表型重叠（图9D）。
GO：生长和GO：解剖结构发育与广泛的IMPC表型共同丰富，代表影响体型，强度和新陈代谢的全身扰动（图9E和9H。
GO：化学突触传递与来自IMPC：Open Field程序的表型共富集，从而指出突触功能障碍与运动受损，焦虑和探索性行为表型之间的联系（图9F）。
GO：声音的感觉感知与IMPC表型共同丰富，主要在IMPC：听觉脑干反应（ABR）和IMPC：声学惊吓和脉冲前抑制（PPI）程序中。这是有道理的，因为ABR程序[37]直接针对听力敏感性，而PPI程序主要用于评估感觉运动门控（感觉事件抑制运动反应的能力）[38]（图9G）。
为了结束GO共富集部分，我们注意到这种类型的分析将具有最大的功率，并在IMPC中的所有基因KO都经过表型分析并分析数据后提供最佳见解。我们在这里的讨论只是为了让我们了解对完整数据集的最终分析将提供哪些见解。
中压扰动的因子分析
MV模型拟合协方差结构的特征分解（方法-交叉验证和模型平均中的（26））表明，75%的相关结构由前20个特征向量解释;S8 图绘制了所解释的累积方差。我们将这些特征向量旋转到一组稀疏的、可解释的载荷或因子，如图10A所示。在这种情况下，稀疏性的重要概念（如图10A所示）是，在任何特定因子处的绝大多数表型的载荷接近于零（即，它们呈绿色）。因此，每个因子都定义了一小簇表型，这些表型具有大的正负载荷或小的负载荷。从生物学的角度来看，每个表型簇往往受到基因KO的集体干扰。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
图 10. 潜在稀疏因素的表征 。
（a） 20个系数的稀疏载荷;对每个荷载向量进行符号和缩放，使得最大荷载的大小为 1，>50% 的有效因子评分扰动为正，如图 （c） 所示。（b） 优势比作为对因子之间扰动依赖性的度量（结果-MV扰动的因子分析）。（c） 每个系数的撞击百分比（即受到严重扰动的百分比），分别用红色/蓝色表示在与载荷相同/相反方向上受到扰动的线路的百分比。例如，标有“Body size （?）”的条形图中红色的大比例表示大多数扰动与因子载荷的方向相同，即它们倾向于导致Body尺寸减小。用于生成此数字的数据和代码可在 [13，14] 中找到。
https://doi.org/10.1371/journal.pbio.3001723.g010
通过检查每个表型簇，并考虑其负载的迹象，我们手动策划描述每个因子的生物学解释的标签。例如，第一个因子根据骨矿物质含量，骨面积，瘦体重，体长和心脏重量的负负荷来定义一个簇;因此，这个因素被标记为“身体大小（?）”。因子标签显示在图10B的轴上。后缀（+）或（?）表示图10A中载荷符号所隐含的效应的方向性，例如，“体型（?）”表示“体型减小”，“耳聋（+）”表示“耳聋增加”。我们测试对应于每个因子的分数中的扰动（方法 - 因子模型）;图10C绘制了具有显着扰动因子分数的KO线的比例，根据扰动是正还是负（沿着图10A中相应的载荷矢量定义的轴）进行分离。这些因子根据具有显著扰动的线的比例进行排序;最常见的扰动在27.5%的线路中被识别出来，而最不常见的扰动在3.7%的线路中被识别出来。荷载向量的符号被定义，使得大多数扰动是正的，导致正方向上的平均率为69.1%。
我们通过分析每对因子的扰动显著性计数的 2×2 列联表来描述因子的统计共扰，其中因子扰动以二进制形式分层 {0 ≡因子没有显著扰动，1 ≡因子被显著扰动}。在几乎所有情况下都有共同扰动的重要证据;具体而言，在190个案例中，有187个案例拒绝了跨成对因子的扰动独立性的原假设（Fisher对2×2个表的精确检验，Fdr控制在5%）。图10B显示了比值比（OR），量化了每对因子的统计共扰动。在这里，对特定因子对之间的OR的解释是，观察一个因子中的扰动会将观察另一个因子中的扰动的几率乘以OR。在图10B中，有一些因素组倾向于相对强烈地共扰，例如，（体型，心功能不全）;（活动/探索因素，协调/平衡，睡眠时间长短，每日睡眠百分比）。一些因素，如活动/探索3，每日睡眠百分比，中性粒细胞：淋巴细胞比率与许多其他因素强烈共扰。其他因素，如耳聋，与其他因素的共扰性较弱（由图10B的耳聋（+）行中相对接近1的OR表示）。
讨论
IMPC揭示了KO诱导的表型扰动中明显的依赖结构。在这里，我们已经证明，一些相关性可归因于对单个潜在表型的多次测量（主要是图3A中的程序内相关性），其中一些可归因于多效性基因效应（特别是图3A中的一些程序间相关性）。鉴于这种结构，可以预期跨表型共享信息可以极大地帮助注释。在这种情况下进行中压分析具有挑战性，尤其是因为数据集的大小和实验设计引起的复杂的采样间相关结构。我们开发了一种可组合的2阶段MV建模方法，可以解决这些问题。
在测量数据设置中，命中率从UV模型的1.4%增加到MV的10.5%是值得注意的。在缺少测量数据的情况下，MV模型的命中率为1.3%，而UV模型在观测数据上的命中率为1.4%，这实际上是有用的。为了验证结果的有效性和一致性，我们实施了几项单独的措施和检查。我们方法的一个基本要素是通过同时期测量的WT动物的计算机重新标记来生成逼真的合成零线。这些合成零值支持多项分析，包括增强实验相关性R的估计，但它们最重要的作用是在表型注释的假设检验中校准假阳性率。这一特殊应用强调了对照动物在表型分析管道中盲法表型的重要性，以及WT动物遵循与KO动物相同的实验设计的效用，例如，在共享窝，天和其他实验协变量方面。
历史上一直存在关于动物表型注释在不同实验室中的可重复性的担忧，尤其是行为表型[39，40]。在这里，我们通过{?1，0，1}中的有符号注释的镜头专注于可复制性，±1对应于特定方向上具有统计显着性的扰动，0表示没有显着影响。有符号注释对于IMPC的有效科学影响至关重要，因此其可复制性是任何统计方法的基本下游要求。在模型物种的表型筛选中，样本量和功率受到严格限制。因此，我们希望在我们的可复制性研究中观察到许多在一个实验室中称为参考系表型命中，但在另一个实验室中则不然。然而，当 2 个实验室都调用命中时，一致的命中对的比例是可复制性的有用度量。在这里，我们已经表明，IMPC参考线中的一致性水平和可复制性很高，特别是显示出与= 1.2%（95%置信区间：0.6%至2.4%）的低Fsr的兼容性。值得注意的是，在图6和S5 Fig的所有程序中观察到实验室的成对一致命中，包括行为测试，如开放场和明暗测试。这种对相干性的检查之所以可行，是因为包含了参考线，并证明了实验设计的价值，该设计在潜在的异构测量环境中结合了技术复制。
作为进一步的相干性检查，我们量化了同一KO系的杂合子和纯合子之间表型命中的方向一致性。在这里，在MV模型的结果中观察到的一致与Fsr估计的3.7%（95%CI：2.8%至5.0%）兼容，即使这里观察到的不一致水平可能被任何具有真正相反效应的杂合子/纯合子KO夸大。
我们已经在MV模型下证明，对于未进行测量的基因 - 表型对，可以调用具有相对较高功率的命中。这有可能增强IMPC数据库的科学影响力，因为与具有缺失数据块的矩阵相比，完整的基因逐表型注释矩阵为最终用户提供了更百科全书和通用的工具（图2）。我们已经通过多种方式评估了在存在缺失数据的情况下推理的准确性。对于参考线，缺失数据结果的可复制性与测量数据结果的可重复性相当（图6）。另外，我们已经通过一个过程输出（LOO-MV）交叉验证技术演示了整个缺失过程的高功率（图5和S4图）和准确性（S9图）。这些结果表明数据中存在一定程度的冗余，包括跨程序，从某种意义上说，某些表型中的大多数信息都是在其他表型中捕获的。这里开发的缺失数据方法有可能用统计分析取代一些动物实验，符合NC3Rs[24]。如果成本较低的实验可以使更昂贵的实验变得多余，那么这将是特别有效的，其中成本包括动物痛苦的道德成本以及对资金和其他资源的考虑。
我们专注于148种定量表型，但IMPC还包括许多分类表型。有用的是，我们的MV模型可以直接扩展到混合反应类型的多种表型。这是因为它接受UV估计的效应大小和标准误差作为输入，并且这些可以像从普通线性模型估计和标准误差一样轻松地采用由逻辑回归输出的估计对数OR的形式。我们预计该方法将朝着这个方向扩展和应用，从使用广义线性多级模型的紫外线分析开始。
[22，23]中介绍的关于估计效应大小和标准误差的MV模型具有许多优点，可能使其在其他应用领域有用。初始紫外分析可以任意复杂，允许对样品之间的相关结构进行仔细的紫外建模;在UV阶段没有必要同时考虑不同响应变量之间的相关结构，这在分析和计算上都是困难的;UV模型拟合可以并行进行;输入到MV分析中的数据集的大小大大减少，可能减少一个数量级或更多。对于此方法的富有成效的应用，某些数据集属性是首选属性。我们需要足够数量的独立中压观测值（在本例中为 KO 线）来有效地估计 Σ 和 R 中的协方差结构。如果没有足够的数据来估计完整的 P×P 协方差矩阵，则可以使用简化的秩因子模型更简洁地表示 Σ 和 R，就像我们在这里所做的那样。UV模型相关结构（在我们的例子中跨动物）或MV模型中的独立性假设（在我们的情况下是跨KO线的i.i.d.效应向量）的错误指定可能导致输出校准错误。我们使用排列生成的合成零线校准了模型输出，我们希望这在其他应用程序中也是一个有用的（如果不是必需的元素）。
经验贝叶斯推理方法具有主要优点，因为它允许预先完成估计Σ和R并定义适当显著性阈值的计算机密集型工作。此后，随着新KO生产线上的数据或现有生产线上的进一步测量值的可用，可以快速地进行计算处理并快速更新估算值。另一种方法是针对后验执行完整的贝叶斯推断。我们确实通过MCMC实现了完整的贝叶斯，但发现它不如MAP估计和经验贝叶斯推理那么可行;这主要是由于担心MCMC采样器的缓慢混合。
一个更雄心勃勃的目标是将完整的多级因子模型直接拟合到原始数据[41，42，43，44，45]，即以后验p（Σ，R，Θ|Y），其中 Y 是原始的动物级数据。原则上，这是非常有吸引力的，因为它可能允许更多的信息将原始数据流向感兴趣的参数，并且可以提供更多的功率。当然，有效的推理将依赖于（更复杂的）模型是数据的良好表示。参考（1）处的多级UV模型，联合多级因子模型可能需要基于每个不同随机效应的跨测量（P×P）协方差结构;这将是不平凡的实现，特别是对于不可识别性的考虑。一个更基本的挑战是数据集的大小增加了一个数量级，这可能会对计算复杂性产生相当大的影响，具体取决于实现。虽然本文的范围是构建和扩展[22，23]的模块化框架，但我们确实认为联合多级因子建模是未来探索的一个有前途的领域，特别是随着复杂模型的可扩展优化方法的持续开发[46]。
稀疏因子模型的开发将包含75%MV扰动的空间的维数从148降低到20。即使在缩小的20维空间中，我们也观察到注释中存在很强的因子间相关性，这表明有效维数仍然较小。在这种潜在的因素空间中，更容易将特定的KO线置于更广泛的上下文中。例如，我们可以确定哪些因素受到扰动并检查它们的特定性质，例如它们的生物学解释，扰动在IMPC中更普遍的罕见程度，以及扰动的方向性是常见还是罕见（图10C中的正效应或负效应）。
总之，我们开发了一种可组合的MV方法，用于分析来自IMPC的高维数据集，展示了与现有UV方法相比的4个主要改进。首先，通过纯分析手段可以大大提高检测KO扰动的能力，在观察到的数据上产生7.5倍的基因表型命中率。其次，即使KO系缺少一些测量值，我们也可以以良好的功率调用缺失测量值的命中数，并输出完整的基因 - 表型图。第三，MV方法的更大功能增强了IMPC表型与现有GO数据库之间的对应关系，随着IMPC数据库的完成，有望获得更强大的生物学见解。最后，高维表型扰动可以在更小的20维亚空间中进行信息查看，从而促进基因KO效应的解释并阐明小鼠基因组表型景观中的丰富结构。
方法
缺少数据时的 MV 模型
在这里，我们将（2）–（4）中引入的MV模型推广到未观察到某些测量子集的情况下。下面 MV 缺失数据模型的有效性依赖于满足 MAR 假设 [23，33，34] 的数据集，该假设在结果 – 数据缺失时的推断中讨论。如果，在基因g处，只有Pg然后使用 * 下标 θ 观察≤P 测量值*g，R 表示对 P 的限制**g测量数据的索引，写入缺少数据的模型
(5)
(6)
with Eq (6) unchanged from the fully observed model, i.e., still with θ·g denoting the full P-vector of latent perturbations for gene g. When data are MAR, the posterior for θ·g under model (5)–(6) is a Gaussian mixture available in closed form:
(7)
EM algorithm
A detailed derivation of the EM algorithm used to fit the MV model introduced at (2)–(4) is given in S2 Note. Algorithm 1 outlines the computations required where, for notational brevity, we use superscripted assignment notation in place of for loops; for example, performs an assignment for each (g, m, s)∈{g = 1…,G}×{m = 1…,M}×{s = 1…,S} .
Algorithm 1 EM algorithm targeting Σ1:S, π.
Inputs:
UV model outputs
Estimated noise correlation
Known scale parameters ω1:M
Initialize: repeat where
until convergence
Outputs:
Initialisation.
The key parameters to initialise are the S covariance matrices Σ1:S. We initialise Σs在使用指定的样本子集计算的样本协方差矩阵处，记作：
(8)
其中，UV结果的P×G矩阵;中的任何缺失数据都是零填充的，仅用于上述计算目的。如果这产生了一个正的半定性，我们添加εI以确保初始化时的正确定性;对于此处显示的结果，我们使用 ε = 0.05。
我们的主要结果基于单个协方差矩阵（S = 1）的情况，在这种情况下，包含训练集中的所有（非合成零）样本。当模型被指定为具有 1 个以上的协方差矩阵 （S>1） 时，我们选择子集，通过基于模型的聚类分析对训练集进行分区，该聚类分析使用 R 包 mclust 中的函数 Mclust（） 和默认参数设置。
MV模型（2）–（4）的可能性是多模态的，因此，EM算法的收敛对初始化很敏感。作为方法 - 模型检查和灵敏度分析中的模型检查部分的一部分，以多种方式研究此灵敏度。通过重复我们对数据子集的整个分析，例如，大小为500，我们捕获了经验初始化的变化（因为仅基于训练数据）以及数据子采样的可能性面的变化;我们证明我们的结果对这两种类型的变异的组合是稳健的。
经验协方差矩阵 （8） 处的初始化有助于使 EM 算法能够以全局最优为目标。为了证明这一点，我们研究了Σ的随机，香草初始化，设置
(9)
我们在单协方差矩阵情况 S = 1 中对大小为 2，000 的 10 个数据子集执行此随机初始化（这些数据子集与用于主要交叉验证分析的数据子集相同）。然后，在每种情况下，我们检查交叉验证似然拟合的值，并将其与样本协方差初始化的交叉验证似然拟合进行比较（S10图，方法 - 交叉验证和模型平均和方法 - IMPC数据的交叉验证似然）。在10倍中，随机初始化的拟合在所有情况下在CV可能性方面的表现都表现得更差，说明了在这种情况下使用监督初始化来减轻优化的非凸性的好处。
EM 算法有潜在的增强功能，以增加收敛到全局最大值的概率，例如 [23] 的拆分和合并算法。虽然我们的基本EM实现似乎为这里考虑的数据集提供了良好的性能，特别是在合理的初始化下，但它可以有效地扩展以在将来包含此类增强功能。
收敛。
当连续迭代之间目标函数的变化低于公差阈值时，EM 算法被视为已收敛。我们自适应地选择公差阈值，参考样本之间对数似然贡献的变化。具体而言，表示第 g个样本在第 t次迭代时对数似然的贡献，公差设置为 ，其中 MAD（） 表示中位数绝对偏差，N特拉训练样本的数量，以及ε托用户指定的常量（我们使用了ε托= 10?4).
错误率控制
当真正的数据生成机制位于假设的统计模型空间内时，模型拟合的统计度量（例如 p 值和贝叶斯因子）是特别有用的工具。当模型空间排除真实机制时，具有统计学意义的度量可能会被错误校准[47]。这在高度结构化的科学数据中尤其重要，其中错误地假设模型中的条件独立性可能导致人为地紧密置信区间和夸大的测试假阳性率。
我们的解决方案是使用一种非参数方法来控制错误率，称为Westfall-Young排列方法[32，48]。这种方法的本质非常简单 - 我们生成合成的零数据，这些数据尽可能精确地模仿实际数据的结构，但通过设计不会系统地偏离WT动物;在我们的方法中，合成零线是从WT样本中随机绘制的。这些空线作为我们实现Westfall-Young排列方法的错误率控制[32，48]中的真零假设集。
我们参考了许多不同的错误率，因此在表3中，我们提供了一个术语表，该术语表将符号与简要描述相关联，并且定义了每个错误率。对于我们论文中的主要分析，我们使用Westfall-Young排列方法进行控制。我们监控和.最后，我们还在方法中控制lfsr（使用Westfall-Young，并且名义上也使用），以便与现有方法进行比较以进行基准测试。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 3. 本文中提到的错误率词汇表。
https://doi.org/10.1371/journal.pbio.3001723.t003
合成空数据。
我们将合成零线定义为通常随机选择的10至20 WT动物的子样本，以反映实际KO线的实验设计特性。合成零值在我们的推理中起着重要的方法论作用：最重要的是在调用表型命中时对Fdr的基于排列的控制，以及在（2）处估计MV模型中的实验相关矩阵R。合成零线是通过从单个中心随机选择WT动物组来生成的，以便与该中心特定真KO线的实验设计特征相匹配。具体来说，对于真正的KO线的每一窝，我们从同一中心的计算匹配的WT窝中取样。对于在d天首次表型化的含有l动物的KO窝，我们从同一表型中心的所有可能的WT窝中抽取WT窝，其中至少有l只动物，并从该窝中随机选择l动物。在时间接近d天的凋落物被选择的概率更高[11]。
假设检验和 Fdr.
IMPC的一个重要输出是数据驱动的（表型，KO基因）对列表的汇编，其中有证据表明表型受到KO基因的干扰。这就引出了检验原假设 H 的分析目标。0： θ页码= 0，在受控的误报率下具有高统计功效。IMPC 数据具有许多级别的复杂结构，导致模型错误指定和参数化测试的假阳性率升高。此外，IMPC的大量通常密切相关的测试需要有效的功率保持多重测试校正。通过使用Westfall-Young排列方法控制Fdr来应对这些挑战，该方法为模型错误指定提供了鲁棒性，并在测试相关时结合高统计功效[32，48];合成空线用作真空假设的集合。
为了检验 KO 线 g 中表型 p 无扰动的原假设，我们使用定义为后均值与后部 SD 之比的 z 统计量，即
(10)
具有 UV 模型输出的相应定义。我们选择一个显著性阈值，记作τ，这样如果
(11)
则线g在表型p处称为显着扰动，方向性由z的符号决定页码.
我们选择τ来控制Fdr。我们使用“贝叶斯”Fdr定义[30]：
(12)
其中 H0是原假设，T 是检验统计量，是临界（拒绝）区域，选择该区域以控制相应的经验形式 [30]
(13)
（13）中Fdr的定义是保守的，因为我们的控制意味着对Benjamini-Hochberg FDR的类似控制[30，50]。我们选择（13）主要是出于方便：合成零数据允许估计和控制（13）分子中的项。（13）中的其他项可以直接处理：分母是已知的，并且先验的P（H）0true） 可以指定，保守地指定为 1，就像我们在这里所做的那样，或者在先前信息可用时提供信息。我们在2个粒度级别上估计Fdr：表型 - 基因对和基因。在表型-基因对水平（即，对于单个测试），Fdr单是因否定原假设而产生的 Fdr，在每个表型 p 处，基因 g 对|z页码|≥τ：
(14)
我们估计Fdr单(τ） 通过
(15)
其中 z页码从合成零线包含在分子中以估计分子中的第二项（13）。
在基因水平上，Fdr完成是在每个基因上否定所有表型的完全原假设得出的 Fdr，其最大值p{|z页码|}≥τ，
(16)
其经验形式定义类似于（13）。我们估计如下
(17)
其中再次 z页码从合成零行包含在分子中。
我们监控Fdr的相关但不同的估计。然而，正是我们习惯于选择τ来控制，同时通过优化来最大化功率。
(18)
基于排列（即基于合成零）的控制是Westfall-Young排列过程的实现[31，32，48]。
可复制性和错误符号率。
当KO线在多个实验室中进行表型分析时，在两个实验室中以相同的方向调用命中（识别显着扰动）支持方法的可复制性。在图 6 中，一致命中对应于蓝色阴影区域中的点。相比之下，在一个实验室中表型增加而在另一个实验室中表现为减少表型的命中意味着2个命中中至少有一个是假阳性（由图6中的红色区域表示）。在我们的分析中，我们检查了跨上下文的这种一致性（刚刚介绍的实验室之间，并且还比较杂合子与纯合子）。能够将观察到的可复制性程度与基础错误率相关联是很有用的，因为这提供了对有效错误率控制的额外验证。因此，我们开发了一种量化（不）一致性的方法：基于可复制性的Fsr估计，我们在下面表示并推导。该方法将有符号注释的列联表映射到兼容的 Fsr。
感兴趣的一般情况是在 2 个条件独立的上下文中输出成对的有符号显著性调用。在本节中，我们用 {?1，0，1} 表示调用，其中 1 和 ?1 表示显著的正或负表型扰动，分别表示 z>τ 和 z
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 4. 方法 A 和 B 之间的可复制性表。
https://doi.org/10.1371/journal.pbio.3001723.t004
我们将证明，虽然表 4 提供的关于 Fdr 的信息很少，但表 4 背后的概率可能与 Fsr 相关，Fsr 定义为在否定零效应的原假设的情况下，错误估计效应符号（产生“S 型错误”）的概率 [49]：
(19)
我们通过考虑以下比率来激励我们对 Fsr 估计器的推导，该比率随着表 4 中的不一致程度而增加：
(20)
我们注意到，其中 A，B∈{?1，0，1}，并且我们表示
(21)
(22)
其中我们定义了ψ?P（θ = 0|AB≠0），其中ψ可解释为“双重错误发现率”，即，在2次有条件的独立测试中，如果零（θ = 0）为真，则为真的概率为真，例如，在2个不同实验室收集的数据集上（请注意，ψ = O（Fdr）2）在Fdr的合理控制下很小）。此外，在从（21）到（22）的步骤中，我们假设假阳性同样可能处于正或负方向。我们还在从 （21） 到 （22） 的步骤中使用了以下内容：
求解 Fsr 的 （22） 得到：
(23)
（23） 的右侧是 ψ 的递减函数，因此我们通过设置 ψ = 0 来定义 Fsr 的保守（略微向上偏倚）估计器：
(24)
其中 q 已被 （20） 处定义的估计器替换。我们通过替换（代替 in （24））精确的二项式置信区间来获取 Fsr 的近似置信区间，该限值在一个模型下推导，其中分歧数遵循二项式分布，成功概率为 q：
交叉验证和模型平均。
高度参数化的统计模型可能会过度拟合数据，从而导致样本输出性能不佳。这种过度拟合问题适用于这里的MV模型，因为它具有灵活的高维协方差矩阵参数化，这是通过经验学习的，尽管通过（4）处的因子模型表示的结构正则化在一定程度上缓解了它。为了防止过度拟合，所有MV结果都是在交叉验证框架内推断的，通过该框架，我们将数据集拆分为“训练”和“测试”集C次，然后使用贝叶斯模型平均在拆分之间合并测试集结果。
我们将每行g表示为在真KO行的集合中，或者在合成零行的集合中;这两组中的每组都包括 N托特= 4，548 行（因为每条合成空线都与真正的 KO 线的设计相匹配）。
对于每个 1,...,C 的交叉验证拆分 c，我们随机划分为大小为 N 的训练集特拉和大小为 N 的测试集托特?N特拉.我们以类似且独立的方式随机划分为 和 。我们继续使用训练基因进行估计，即.然后，我们估计使用测试基因的条件，即.使用贝叶斯模型平均在交叉验证拆分中合并的测试集估计值，即
(25)
它表示组合的后验是分裂特异性后验的混合物，每个后验都是高斯混合体。我们将后兆 MV 估计值定义为 （25） 中组合后部的平均值和标准差，并将这些估计值转发到受控 Fdr 下的表型调用。我们用于交叉验证的经验贝叶斯推理的框架如表5所示。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 5. 交叉验证 ed 经验贝叶斯推理。
https://doi.org/10.1371/journal.pbio.3001723.t005
此外，我们使用测试集数据，通过贝叶斯模型在C分割特定模型中求平均值来计算合并协方差估计，即：
(26)
从（26）到因子分析（方法-因子模型）。
模型检查和灵敏度分析
敏感性分析。
我们通过比较不同交叉验证折叠的结果来验证我们的下游因子分析是否可靠。具体来说，我们将最终估计值（26）中的varimax旋转因子载荷与每个折叠c中的可变最大旋转因子载荷进行比较，并根据对称的KL背离（），
(27)
载荷在S11图中绘制并进行比较，并且在质量上看起来相似，仅在少数因子上存在差异。这些有限的差异是由于2个分解中的因子的合并或分裂。我们的结论是，因子分析对数据子采样相对不敏感：在最坏的情况下，因子分解之间的差异应该只发生在少数因子上。
数据子采样。
在这里，我们检查了结果对KO线的潜在MV异方差性的稳定性，这些KO线可能无法被我们的MV混合模型捕获。我们执行灵敏度分析，其中我们随机对总共4，548条线中的500行进行子采样，并将它们用作训练集以重新拟合MV模型。我们执行此子采样 c = 1,...,10 次，每次从训练集进行估计，并保留来自测试集的 MV 模型表型调用。我们发现褶皱 c′在 （27） 处与 as 之间具有最大的对称 KL 散度。
在表6中，我们将fold c′的有符号表型调用与完整分析中的相应调用进行比较（即（26）处的交叉验证和模型组合分析）。不一致程度很低，在8，316个实例中总共有37个分歧，其中两个模型都称之为命中。（请注意，我们无法从表 6 中有效地估计 Fsr，因为完整和子采样分析的测试结果之间存在条件依赖性。一个令人放心的质量控制检查是，选择大小为500的最不一致的子样本只会导致这种小程度的不一致，这表明减少样本量，同时降低功率，不应该导致与全数据分析结论的任何定性分歧。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 6. MV模型的有符号表型命中数比较，应用于训练大小N = 500的最大KL发散子采样数据集（左图）与训练大小N = 2000的完整数据集（顶部）。
我们用 {?1，0，1} 中的数字表示呼叫，其中 1 和 ?1 分别表示显著的正表型和负表型扰动，0 表示缺乏统计显著性。
https://doi.org/10.1371/journal.pbio.3001723.t006
预测屏蔽数据。
从图2可以看出，KO线的整个过程中经常缺少数据。为了检查缺失的数据推断，我们执行以下“掩码-预测-比较”算法：（i）对于测试集中的每一行KO，依次人为地屏蔽来自其每个测量过程的数据;（ii） 预测被屏蔽数据背后的扰动;（iii）将预测的扰动与UV模型在未屏蔽数据上估计的扰动进行比较。我们将屏蔽数据的推断称为一个程序输出MV（LOO-MV）。LOO-MV和UV结果之间的不一致程度较低（S9图），并且兼容= 0.4%（95%置信区间：0.2%至0.7%）。这与即使在存在缺失数据的情况下执行推理时，也要很好地校准误报率。
检查MV模型和IMPC数据库之间的不一致。
回到结果 -与 IMPC 数据库和表 1（B） 的比较，在这里，我们通过检查哪种方向性（我们的 MV 模型或 IMPC 数据库）在生物学上更明智，更仔细地检查了这 3 个不一致的案例。我们使用经验贝叶斯来量化关于任何特定表型p的命中方向性的先验信念作为概率：P事先(θ页码>0|θ页码≠0）。这涉及汇总来自该表型在所有基因中的命中方向性的信息，这可以通过简单的平均值来完成：
(28)
由于我们的UV模型的调用与IMPC数据库中的调用之间没有分歧，因此我们在（28）中包括来自这两种方法的命中，但我们不包括来自MV模型的调用。我们使用（28）中预先定义的方法来分析3个有符号表型命中，显示我们的MV模型和IMPC数据库之间存在分歧，产生1.45的贝叶斯因子支持MV模型，但这是一个弱贝叶斯因子，没有提供实质性证据，无论哪种方式，模型的输出方向性在这些分歧实例上最明智;这个负结果是有意义的，因为3的小样本量导致低功耗。
生物应用方法
基因本体分析。
我们使用R包GOfuncR来测试GO项和IMPC表型之间的共富集。该软件包的一个重要特征是，它纠正了多个测试和测试的相互依赖性，使用基因相关变量的随机排列来控制家庭错误率。我们创建IMPC-表型基因集，这些基因不仅受到显着扰动，而且在所有基因中表现出至少2倍于效应大小的SD的效应大小;在我们引入的数学符号中，表型p的IMPC基因集定义为：
(29)
(30)
我们应用额外的过滤器来仅关注纯合子KO。我们使用5%的家庭错误率阈值（当测试单个IMPC基因集对所有BP GO项的共富集时，一个或多个假阳性的概率被限制为小于5%）。
我们对每个 IMPC 表型执行 1，000 个 GO 图排列。背景基因集（也称为基因宇宙）被定义为MV模型的所有纯合子 - KO基因，其中有一些表型测量值可用（总共2，628个基因）;对于UV模型，背景基因集包括所有纯合子-KO基因，其中这种特定表型可用。基本的推理工具是 Fisher 精确测试，用于在 2 x 2 列表中对行和列分类的独立性进行精确测试，如图 7 所示。在GOfuncR中实施此基本测试可确保正确控制错误率（此处为家庭错误率）。
因子模型。
为了便于解释表型扰动，我们计算了相关矩阵的特征分解，即
(31)
在 varimax（） 函数中，默认参数在 R [51] 中。
通过 λ 表示旋转的稀疏载荷 P 向量l， l = 1,...,20，第g个KO基因的第l个因子评分为。
执行假设检验以识别因子评分中的显著扰动。表示 和 ，我们形成检验统计
(32)
控制Fdr类似于方法 - 错误率的控制，其中表型p = 1,...,P，我们现在的因子评分l = 1,...,20。
与现有方法的比较
极端反卷积 （XD）。
[23]中使用的模型和随附的软件包Extreme Deconvolution（XD）的基础模型类似于我们使用（2）–（3）的可能性，但有一些区别。XD 具有约束 M = 1，即没有多个缩放参数 ω1：米在 [22] 中引入。XD 推广到具有非零均值μ的基础混合物组分1：S这些值本身是估计的，即 N（θ·g|μs，Σs).出于方法比较的目的，XD 与 μs≡0，因为零均值模型适用于我们在这里分析的数据集。
XD 使用与我们类似的 EM 算法来最大化 Σ 联合的可能性1：S和π（和μ1：S更一般地）。在整个 XD EM 算法优化过程中，每个 Σ 的秩s与其初始化值的秩相同 [22]。我们将 XD 初始化为相同的值 Σ1：S并在初始化自己的模型时π（除非我们运行XD为mash生成数据驱动的矩阵，在这种情况下，我们遵循[22]中的指示）。总之，我们模型的拟合和XD的拟合之间的任何差异主要是由没有缩放参数ω驱动的。1：米在 XD 和 Σ 的因子模型正则化中1：S.
多变量自适应收缩（麦芽浆）。
我们的方法的模型似然（2）–（3）与[22]中引入的相同，这是软件包混搭的基础。[22]中一个特别重要的见解是在用ω表示的阶梯上引入多尺度混合物1：米.作者指出了这种方法在多组织eQTL背景下的效用，我们发现它对MV小鼠表型数据也很有用。我们相信[22]的多尺度协方差模型形式有可能增强广泛科学学科的MV推理。
我们的方法和混搭之间的关键区别在于协方差矩阵Σ1：S被定义和估计。我们使用少量正则化协方差矩阵（我们在这里考虑 S = 1，2）参数化模型并优化 Σ1：S并作为模型拟合的一部分集体π。相反，麦芽浆生成并固定一个较大数量的协方差矩阵Σ（S = P + 10）1：S，在优化 （2）–（3） π之前。
更详细地说，mash 生成 2 种不同类型的协方差矩阵：8 个数据驱动矩阵和 P+2 规范矩阵。输入到mash中的8个数据驱动协方差矩阵是经验表型协方差在表现出最大表型效应的样本中的低秩表示。其中三个数据驱动的协方差矩阵是使用 XD 软件生成的 [23]。此外，生成 P+2 规范 P×P 协方差矩阵，包括单位矩阵、1 矩阵，对于 p = 1,...,P，其中 ep是一个 P 向量，除设置为 1 的第 p个元素外，所有位置都有零。
虽然麦芽浆的一个优雅方面是，相对于π给定Σ的优化1：S是凸的，其协方差矩阵的生成涉及XD软件中的非凸优化，因此对初始化可能有一定的敏感性[22]。我们的 EM 算法相对于 Σ 的 MAP 优化1：Sπ是非凸的，我们在方法初始化中研究了初始化的敏感性。
命中率、错误率和模型拟合度。
我们在 IMPC 数据上比较了不同方法的功率（命中率）和错误率（估计的 Fdr 或 Fsr）。我们考虑了许多假设检验框架，由它们的检验统计数据和关键区域（拒绝标准）决定。
用于测试基因 g 处表型 p 扰动的检验统计量：
z 统计量
本地错误符号率 [27]：
(33)
直观地“如果我们使用我们对符号的最佳猜测（正数或负数）”，我们会错误地预测效应迹象的概率“[22]。
关键区域：
控制 Fdr完成通过基于排列的测试统计阈值 τ <5%，即，临界区域|z页码|>τ 或 lfsr页码<τ（方法–错误率控制）和 [32，48]);
名义上控制本地 Fsr lfsr页码<5% [22].
表7显示了3种不同错误率控制方法下的命中率和错误率，子表对应于上面定义的测试统计数据和关键区域：表7A表示A1（测试统计量A与关键区域1）;表7B为B1;表7C表示B2。命中率根据原始数据是测量还是缺失（95% 的非参数自举 CI）分层显示。我们显示错误率估计值 、和方法-假设检验和 Fdr 中定义的误差率估计值。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 7. 不同型号的命中率和错误率比较。
显示论文中主要分析中使用的模型和错误率控制的行以粗体突出显示。最高命中率带有下划线。
https://doi.org/10.1371/journal.pbio.3001723.t007
交叉验证 IMPC 数据的可能性。
在方法 - 交叉验证和模型平均中描述的推理框架中，我们计算仅使用训练数据拟合的模型下测试集数据的可能性。参照 （2）–（3），每样本对数交叉验证的折叠 c 似然为
(34)
其中是折叠c（方法 - 交叉验证和模型平均）测试集中的KO基因集。表 8 显示了基准测试结果，其中我们显示了跨 10 个交叉验证折叠获取的 （±2×SEM） 的平均值。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 8. 在 IMPC 数据上比较跨 MV 模型的交叉验证的对数似然。
显示论文中主要分析中使用的模型的行以粗体突出显示。CV 日志的最大可能性带有下划线。显示的结果是跨折叠c = 1,...,10的（34）的每个样本对数似然，以及±2个SEM区间。
https://doi.org/10.1371/journal.pbio.3001723.t008
其他数据集上的交叉验证似然。
我们在附加数据集上比较各种MV方法，以检查各种方法的相同定性性能是否仍然存在。要使用的自然数据集是Urbut及其同事的多组织eQTL研究[22]，其中首次开发了麦芽浆。该数据集包含16，069个样本，每个样本对应于一个（基因，单核苷酸多态性）对。在44个组织中测量这些表达数量性状位点（eQTLs）。所以，N托特= 16，069 和 P = 44，与 N 相反托特对于 IMPC 数据，= 4，584 和 P = 148。
我们使用方法-交叉验证和模型平均中列出的经验贝叶斯交叉验证框架来分析eQTL数据。我们遵循[22]的方法进行数据预处理和R的估计。我们在 eQTL 数据上的训练折叠大小为 5，000（与 IMPC 数据的 2，000 相比），但除此之外，参数设置是相同的。一个重要的区别是eQTL研究中没有缺失的数据。表9显示了eQTL数据集上的CV日志似然比较。
thumbnail 下载：
断续器幻灯片
巴新放大图像
断续器原始图像
表 9. 在Urbut及其同事的eQTL数据上比较MV模型的交叉验证的对数似然[22]。
CV 日志的最大可能性带有下划线。显示的结果是跨折叠c = 1,...,10的（34）的每个样本对数似然，以及±2个SEM区间。
https://doi.org/10.1371/journal.pbio.3001723.t009
讨论基准测试结果。
我们提出的主要结果基于MV模型（ComposeMV），其中S = 1和K = 20（结果 - 多变量模型中引入的符号），同时控制如方法 - 错误率控制中所述。我们在本节中引用此模型。现在，我们简要讨论了在各种错误率控制方式下对 12 个模型进行基准测试的结果（表 7 和表 8）。
以我们控制的表7A为中心，所有基准模型中测量数据的命中率最大，为10.5%（10.2，10.9）。与其他两种考虑的错误率控制方法相比，该命中率也是最佳的（表7B和7C）。值得注意的是，在将命中率与表7C进行比较时，表7C中监测的错误率通常高于表7A和7B，这归因于使用不同的错误率控制方法，名义上控制lfsr<5%。
在缺少数据时的最佳命中率方面，我们看到其他模型（例如，S = 2和K = 15的ComposeMV）的表现略好，为2.0%（1.8，2.1），而1.3%（1.2，1.4）为，但与 相比，估计错误率更高。
转到表 8 中的交叉验证似然比较，我们看到每个样本的交叉验证对数似然平均值为 = ?53.4 （?53.6， ?53.1）。ComposeMV对此略有改善，S = 2，K = 15，其= ?53.0（?53.2，?52.7）。和ComposeMV通常与现有方法相比是有利的，其性能比醪糟好一些，麦芽浆= ?54.2（?54.4，?53.9），并且比S = 2的XD好得多，后者= ?60.0（?60.7，?59.3）。
对于表9所示的附加eQTL数据集的基准测试，出现了类似的模式 - 当前论文的MV模型，标记为ComposeMV，与mash相比，在CV可能性方面表现更好，但性能明显优于XD。有趣的是，对于eQTL数据基准测试，S = 2的ComposeMV模型在表中表现最好，这表明混合了多个学习的协方差矩阵（除了ω的多尺度阶梯之外）m） 在某些情况下可能特别有用。
支持信息
单变量模型详细信息。
显示 1/13： pbio.3001723.s001.pdf
跳转到无花果共享导航
补充说明1.单变量模型详细信息。将转换分别应用于每个定量表型，并应用于同时跨越所有表型中心。对于任何具有一些观察值的定量表型婷婷≤0，则在变换之前向所有观测值添加一个常量，以满足：最小值（y） = （最大值（y)?最小值（y))/100.然后用指数对Box-Cox进行变换表型5λ被限制在λ∈ {?2,?1.5,...,1.5,2}并选择通过重新视乎λ在普通高斯线性模型下应用于来自有性别的基线动物的数据和日作为协变量。Box-Cox变换后，每个中心-表型对的数据为缩放到零中位数和单位中位数绝对偏差，然后在±20绑定极端数据点的影响。10转化的定量表型p用于 KO 线g在高斯响应下进行分析带日（α日），垃圾（α垃圾）、基因型（θ紫外线页码），性别（β性），应变(β应变），调查员（βinv） 和元数据组 （β元） 作为协变量，并带有惩罚的样条以解释基线动物测量的系统性时间趋势。受到惩罚的样条如 [1] 的第 16 章所述，纯三次多项式分量具有15系数β聚，以及具有系数的全三次样条的基函数αsplk哪些是通过具有方差分量的分层模型进行正则化σ2spl.日和产仔影响使用方差分量分层建模σ2日和σ2垃圾.表示残差方差由σ2残油.对于任何特定的突变系，分析仅限于该线的数据来自同一中心的所有基线动物的数据。模型是：20y我～N(μ我, σ2残油)μ我=θ紫外线页码我(动物我排队g) +α日d[我]+α垃圾l[我]+K+3∑k=1αsplkfk(td[我])+β性s[我]+β性别，基因s[我]，g[我]+β应变j[我]+βinvv[我]+β元m[我]+3∑p=1β聚ptpd[我]α日d|σ2日～N(0, σ2日),为d= 1,...,Dα垃圾l|σ2垃圾～N(0, σ2垃圾),为l= 1,...,升αsplk|σ2spl～N(0, σ2spl),为k= 1,...,千+ 3哪里g指数基因型，s性j应变v调查员，以及m元数据组;td是时间点对应于d第一天。该模型调整潜在的性别 - 基因型相互作用效应[2]，具有总和到零的对比度约束，即∑s∈{M，F}β性别，基因s，g≡0对于每个g- 表示主要基因型效应可解释为男性和女性基因型效应的平均值。这功能fk(·)表示三次样条的 B 样条基的基函数，其节点为规则25天的经验分布的间隔分位数，以及节数，K，向下舍入从唯一天数除以 10。非信息性先验被指定为β和σ2在共轭先验族内可用在所使用的软件包中 （MCMCglmm [3， 4]）。位置参数β已分配独立正态（平均值 = 0，方差 = 100）先验。方差参数σ2已分配30独立的逆伽马（形状 = 0.01，速率 = 0.01）先验。11非信息逆伽马（ε, ε） 先与小ε是差异的常见但实用的选择组件，我们以所使用的软件包中提供的内容为指导。众所周知，可以对特定选择的后向敏感度ε（例如ε从 0.01 到 0.001） [5]以后方法开发我们更喜欢[5]所建议的非信息性半柯西先验。1
引用[1] Ruppert， D.， Wand， M. P. & Carroll， R. J.半参数回归.剑桥系列《统计与概率数学》（剑桥大学出版社，2003年），第一版。[2] 卡普等。哺乳动物表型性状中性二态性的患病率。老挝语35通讯8, 15475 (2017).网址http://www.nature.com/doifinder/10.1038/ncomms15475.[3] Hadfield， J. D. MCMC方法的多响应广义线性混合模型：MCMCglmm R Package.统计软件学报33, 1–22 (2010).网址http://www。jstatsoft.org/v33/i02.40[4] R开发核心团队。R：统计计算的语言和环境.R统计计算基金会，奥地利维也纳（2018年）。网址http://www.r-project。组织.[5] Gelman， A. 分层模型中方差参数的先验分布。贝 叶 斯分析1, 515–533 (2006).452
1 / 13
下载
无花果份额
S1 注意。 单变量模型详细信息。
本说明详细介绍了结果-单变量模型中引入的 UV 模型的技术方面。
https://doi.org/10.1371/journal.pbio.3001723.s001
（英文）
S2 注意。 2. 电磁算法。
本说明包含方法–EM 算法中引入的 EM 算法的技术细节。
https://doi.org/10.1371/journal.pbio.3001723.s002
（英文）
S1 图 IMPC成体和胚胎表型管道。
科学目的，实验设计和每个程序的详细描述可在 www.mousephenotype.org/impress/pipelines 获得。还详细描述了每个程序中的每个表型。请注意，术语参数用于指代我们在本文中所说的表型。我们更喜欢使用表型来避免在统计推断中使用参数的任何术语歧义。
https://doi.org/10.1371/journal.pbio.3001723.s003
（断续器）
S2 图 缩放 z 统计量的热图，用于说明应用于 UV 结果的质量控制过滤器。
KO 线按中心内的时间水平排序。表型中心内的纵向趋势可以表明UV模型未捕获的实验人工制品。在这种情况下，用红色矩形勾勒出来，来自中心-过程对的数据将从下游 MV 分析中省略。用于生成此数字的数据和代码可在 [13，14] 中找到。KO，淘汰赛;MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.s004
（断续器）
S3 图例 针对有数据的基因-表型对的 MV 与 UV 模型的散点图。
轴延伸至 [?3， 3]，而计数适用于所有数据，包括超出图比例的数据。用于生成此数字的数据和代码可在 [13，14] 中找到。MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.s005
（断续器）
S4 图 按表型划分的功率和缺失数据百分比。
顶部面板显示每个表型的缺失百分比数据。下部面板显示表型特异性命中率（即，明显扰动的线的比例），用于UV方法，以及根据数据是否缺失或观察到的数据分层的MV方法。用于生成此数字的数据和代码可在 [13，14] 中找到。MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.s006
（断续器）
S5 图例 可复制性热图，比较不同表型中心的结果。
热图显示 UV 和 MV 模型下的参考线的缩放 z 统计量。显著的扰动 （） 用十字标记。白色方块表示 UV 模型下的缺失数据。在几个实验室中独立测量了七条KO线（标记为顶部），并使用UV和MV模型（从顶部标记的第三行）进行分析。每行对应于一个表型（右标记），按程序分组（左标记）。用于生成此数字的数据和代码可在 [13，14] 中找到。KO，淘汰赛;MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.s007
（断续器）
S6 图 缩放 z?统计量的杂合子/纯合子一致性散射图， .
每个点对应于特定基因的杂合子和纯合子KO谱系。叠加每个显著性组合的计数 （%）;当轴延伸至 [?3， 3] 时，计数适用于所有数据，包括超出图比例的数据。基于不一致程度的 Fsr 估计值（95% 置信区间）显示在面板顶部。用于生成此数字的数据和代码可在 [13，14] 中找到。Fsr，错误符号率;KO，淘汰赛。
https://doi.org/10.1371/journal.pbio.3001723.s008
（断续器）
S7 图 GO项（左）与 IMPC 表型（下图）的共富集，用于 UV 模型调用的命中。
GO项和IMPC表型之间具有统计学意义的共富集用粗体轮廓正方形表示（控制每种表型的家庭错误率<5%）。正方形的颜色表示GO项处明显扰动的KO基因的百分比，这些基因在正方向上改变了表型（参见顶部的比例尺）。IMPC表型由沿水平轴的GO项模式聚类，而BP GO项按表型模式垂直聚类。表型标签根据左下角图例中的程序进行着色。GO项的子集，由右侧的行（a-h）标记，在图7中更详细地检查。为了便于阅读，我们仅绘制至少具有 3 个显著共富集实例的 IMPC 表型和 GO 项。用于生成此数字的数据和代码可在 [13，14] 中找到。血压， 生物过程;GO，基因本体论;IMPC，国际小鼠表型分析联盟;KO，淘汰赛。
https://doi.org/10.1371/journal.pbio.3001723.s009
（断续器）
S8 图 相关结构的累积比例由特征向量Q解释（31）。
虚线表示超过 75% 的相关性由 20 个特征向量解释。用于生成此数字的数据和代码可在 [13，14] 中找到。
https://doi.org/10.1371/journal.pbio.3001723.s010
（断续器）
S9 图. 检查遮罩数据 （LOO-MV） 上的 MV 分析与 UV 模型一致性的散点图。
我们将LOO-MV结果（推断屏蔽数据的扰动）与应用于未屏蔽数据的UV模型的结果进行绘制;请参阅方法 – 预测屏蔽的数据。基于不一致程度的 Fsr 估计值（95% 置信区间）显示在面板顶部。用于生成此数字的数据和代码可在 [13，14] 中找到。Fsr，错误符号率;LOO-MV，一个程序出MV;MV，多变量;紫外线，单变量。
https://doi.org/10.1371/journal.pbio.3001723.s011
（断续器）
S10 图 10 CV 折叠的随机初始化和样本协方差初始化拟合之间的交叉验证 （CV） 对数似然比较。
随机初始化拟合在CV可能性方面表现得系统性较差，支持使用监督初始化来减轻优化的非凸性。用于生成此数字的数据和代码可在 [13，14] 中找到。
https://doi.org/10.1371/journal.pbio.3001723.s012
（断续器）
S11 图 因子载荷的灵敏度分析。
（a） 贝叶斯模型的可变最大旋转载荷在所有交叉验证折叠中平均协方差矩阵。（b） 折叠 c′ 协方差矩阵的可变最大旋转载荷，选择它是为了最大化折叠 c 之间和跨折叠 c 的对称 KL 散度。2 个荷载图在质量上相似，但存在一些小的差异。用于生成此数字的数据和代码可在 [13，14] 中找到。
https://doi.org/10.1371/journal.pbio.3001723.s013
（断续器）
引用
1.Oprea TI， Bologa CG， Brunak S， Campbell A， Gan GN， Gaulton A， et al.人类基因组中未开发的治疗机会。Nat Rev Drug Discov.2018;17(5):317–32.pmid：29472638
查看文章PubMed/NCBI谷歌学术搜索
2.布朗SDM，摩尔MW。国际小鼠表型分析联盟：小鼠表型的过去和未来观点。妈妈基因组。2012;23(9–10):632–40.pmid：22940749
查看文章PubMed/NCBI谷歌学术搜索
3.Brown SDM，Holmes CC，Mallon AM，Meehan TF，Smedley D，Wells S.用于表征哺乳动物基因功能的高通量小鼠表型组学。纳特·吉内特·2018;19(6):357–70.pmid：29626206
查看文章PubMed/NCBI谷歌学术搜索
4.Dickinson ME， Flenniken AM， Ji X， Teboul L， Wong MD， White JK， et al.高通量发现新型发育表型。自然界。2016;537(7621):508–14.pmid：27626380
查看文章PubMed/NCBI谷歌学术搜索
5.Karp NA， Mason J， Beaudet AL， Benjamini Y， Bower L， Braun RE， et al.哺乳动物表型性状中性二态性的患病率。纳特·科蒙。2017;8:15475.pmid：28650954
查看文章PubMed/NCBI谷歌学术搜索
6.Meehan TF， Conte N， West DB， Jacobsen JO， Mason J， Warren J， et al.从国际小鼠表型分析联盟的3，328个基因敲除中发现的疾病模型。纳特·热内特。2017;49(8):1231–8.pmid：28650483
查看文章PubMed/NCBI谷歌学术搜索
7.Bowl MR， Simon MM， Ingham NJ， Greenaway S， Santos L， Cater H， et al.大规模的听力损失筛查揭示了听觉功能障碍的广泛未开发的遗传景观。纳特·科蒙。2017;8(1):886.pmid：29026089
查看文章PubMed/NCBI谷歌学术搜索
8.Rozman J， Rathkolb B， Oestereicher MA， Schütt C， Ravindranath AC， Leuchtenberger S， et al.通过高通量小鼠表型鉴定代谢中的遗传元素。纳特·科蒙。2018;9(1):288.pmid：29348434
查看文章PubMed/NCBI谷歌学术搜索
9.英普瑞斯管道|标准化屏幕的国际小鼠表型分析资源;.可从以下位置获得： https://www.mousephenotype.org/impress/pipelines.
10.Karp NA， Speak AO， White JK， Adams DJ， Hrabé de Angelis M， Hérault Y， et al.时间变化对小鼠敲除表型研究设计与分析的影响.PLoS ONE.2014;9（10）：e111239.pmid：25343444
查看文章PubMed/NCBI谷歌学术搜索
11.de Angelis MH， Nicholson G， Selloum M， White JK， Morgan H， Ramirez-Solis R， et al.通过小鼠诊所联盟的广泛表型筛选分析哺乳动物基因功能。纳特·热内特。2015;47(9):969–78.pmid：26214591
查看文章PubMed/NCBI谷歌学术搜索
12.库尔巴托娃N，梅森JC，摩根H，米汉TF，卡普NA。PhenStat：用于标准化分析高通量表型数据的工具包。PLoS ONE.2015;10（7）：e0131274.pmid：26147094
查看文章PubMed/NCBI谷歌学术搜索
13.Nicholson G. 多变量表型分析能够对哺乳动物基因功能进行全基因组推断。2022. 可从以下日期获得： https://github.com/georgenicholson/multivariate_phenotype_data_and_code.
查看文章谷歌学术搜索
14.Nicholson G. georgenicholson/multivariate_phenotype_data_and_code：多变量表型分析能够对哺乳动物基因功能（数据和代码）进行全基因组推断。2022. 可从以下日期获得： https://zenodo.org/record/6787112.
查看文章谷歌学术搜索
15.Shriner D.在全基因组关联研究中通过多性状分析转向系统遗传学。前热内特。2012;3:1.pmid：22303408
查看文章PubMed/NCBI谷歌学术搜索
16.遗传关联研究中多元表型的分析方法.J 可能统计 2012;2012：1–13.pmid：24748889
查看文章PubMed/NCBI谷歌学术搜索
17.Galesloot TE， van Steen K， Kiemeney LALM， Janss LL， Vermeulen SH.多变量全基因组关联方法的比较.PLoS ONE.2014;9（4）：e95923.pmid：24763738
查看文章PubMed/NCBI谷歌学术搜索
18.Petretto E， Bottolo L， Langley SR， Heinig M， McDermott-Roe C， Sarwar R， et al.使用贝叶斯多组织方法对基因表达的遗传控制的新见解。PLoS计算机生物学. 2010;6（4）：e1000737.pmid：20386736
查看文章PubMed/NCBI谷歌学术搜索
19.Scott-Boyer MP， Imholte GC， Tayeb A， Labbe A， Deschepper CF， Gottardo R.用于多元 eQTL 映射的集成分层贝叶斯模型。2012;11（4）.pmid：22850063
查看文章PubMed/NCBI谷歌学术搜索
20.Flutre T， Wen X， Pritchard J， Stephens M.用于多个组织中联合eQTL分析的统计框架。PLoS Genet.2013;9（5）：e1003486.pmid：23671422
查看文章PubMed/NCBI谷歌学术搜索
21.李G，沙巴林AA，鲁辛一世，赖特FA，诺贝尔AB。用于多组织eQTL分析的经验贝叶斯方法。生物统计学。2018;19(3):391–406.pmid：29029013
查看文章PubMed/NCBI谷歌学术搜索
22.Urbut SM，Wang G，Carbonetto P，Stephens M.灵活的统计方法，用于估计和测试具有多种条件的基因组研究中的效果。纳特·热内特。2019;51(1):187–95.pmid：30478440
查看文章PubMed/NCBI谷歌学术搜索
23.Bovy J，Hogg DW，Roweis ST.极端反卷积：从嘈杂，异质和不完全观测中推断出完整的分布函数。Ann Appl Stat. 2011;5（2 B）：1657–77.
查看文章谷歌学术搜索
24.NC3R.提供自： https://www.nc3rs.org.uk/.
25.Box GEP， Cox DR.转换分析。J R Stat Soc B Methodol.1964;26(2):211–252.
查看文章谷歌学术搜索
26.鲁珀特D，魔杖MP，卡罗尔RJ。半参数回归。剑桥统计与概率数学系列丛书。剑桥大学出版社;2003.
27.Stephens M. False Discovery Rates： A New Deal.生物统计学。2017;18(2):275–94.pmid：27756721
查看文章PubMed/NCBI谷歌学术搜索
28.主教C.模式识别和机器学习|克里斯托弗·毕晓普·|施普林格。第1版纽约：施普林格出版社;2006. 可从以下日期获得： https://www.springer.com/gp/book/9780387310732.
29.Benjamini Y， Hochberg Y. 控制错误发现率：一种实用而强大的多重测试方法。J R Stat Soc B Methodol.1995;57(1):289–300.
查看文章谷歌学术搜索
30.Efron B，Tibshirani R.微阵列的经验贝叶斯方法和错误发现率。Genet Epidemiol.2002;23(1):70–86.pmid：12112249
查看文章PubMed/NCBI谷歌学术搜索
31.Westfall PH，Krishen A，Young SS.使用先验信息为多个端点分配显著性水平。统计医学 1998;17（18）：2107–2119.pmid：9789917
查看文章PubMed/NCBI谷歌学术搜索
32.Meinshausen N， Maathuis MH， Bühlmann P. Westfall-Young 排列程序在依赖性下进行多重检验的渐近最优性.统计年鉴。2011;39(6):3369–91.
查看文章谷歌学术搜索
33.鲁宾数据库。推理和缺失数据。Biometrika.1976;63(3):581–92.
查看文章谷歌学术搜索
34.小RJA，鲁宾DB。缺少数据的统计分析。美国新泽西州霍博肯：John Wiley & Sons， Inc.;2002. 可从以下日期获得： http://doi.wiley.com/10.1002/9781119013563.
35.Ashburner M， Ball CA， Blake JA， Botstein D， Butler H， Cherry JM， et al. Gene Ontology： tool for the unity of biology.纳特·热内特。2000;25(1):25–9.pmid：10802651
查看文章PubMed/NCBI谷歌学术搜索
36.基因本体学联盟。基因本体资源：丰富一个GOld矿。核酸研究 2021;49（D1）:D 325–D334.pmid：33290552
查看文章PubMed/NCBI谷歌学术搜索
37.听觉脑干反应协议。可从以下日期获得：https://www.mousephenotype.org/impress/ProcedureInfo?action=listprocID=149。
38.声学惊吓和预脉冲抑制 （PPI） 协议;.可从以下日期获得：https://www.mousephenotype.org/impress/ProcedureInfo?action=listprocID=148。
39.Crabbe JC， Wahlsten D， Dudek BC.小鼠行为遗传学：与实验室环境的相互作用。科学（纽约州纽约市）。1999;284(5420):1670–2.pmid：10356397
查看文章PubMed/NCBI谷歌学术搜索
40.Kafkafi N， Agassi J， Chesler EJ， Crabbe JC， Crusio WE， Eilam D， et al.临床前研究中啮齿动物表型的再现性和可复制性。Neurosci Biobehav Rev. 2018;87：218–32.pmid：29357292
查看文章PubMed/NCBI谷歌学术搜索
41.穆森博.异质种群中的潜在变量建模。Psychometrika.1989;54(4):557–85.
查看文章谷歌学术搜索
42.朗福德NT，穆森博。聚类观测值的因子分析。Psychometrika.1992;57(4):581–97.
查看文章谷歌学术搜索
43.Ansari A， Jedidi K. 贝叶斯因子分析，用于多级二元观测。Psychometrika.2000;65(4):475–96.
查看文章谷歌学术搜索
44.Goldstein H， Browne W. 使用马尔可夫链蒙特卡洛 （MCMC） 估计的多级因子分析建模。潜在变量和潜在结构模型。2002.
查看文章谷歌学术搜索
45.Goldstein H， Browne W. 连续和离散数据的多级因子分析模型。当代心理测量学：罗德里克·P·麦克唐纳的节日。2005. 第453–75页.
查看文章谷歌学术搜索
46.德什潘德 SK， 罗奇科娃 V， 乔治 EI.使用多变量穗板 LASSO 同时选择变量和协方差。计算图统计 2019;28（4）：921–31.
查看文章谷歌学术搜索
47.Bernardo JM，Smith A. Bayesian Theory。概率和统计中的Wiley系列。Chichester： John Wiley & Sons Ltd.;2000. 可从： https://cds.cern.ch/record/1319894.
48.Westfall PH， Young SS. 基于重采样的多重检验：P值调整的示例和方法。威利;1993.
49.Gelman A， Tuerlinckx F. Type S 错误率经典和贝叶斯单和多重复合过程。计算机统计 2000;15（3）：373–90.
查看文章谷歌学术搜索
50.本杰明尼，耶库蒂利D.依赖关系下多重测试中错误发现率的控制。Ann Stat. 2001;29（4）：1165–88.
查看文章谷歌学术搜索
51.R 核心团队。R：统计计算的语言和环境。奥地利维也纳：R统计计算基金会;2021. 可从以下日期获得： https://www.R-project.org/.