科研|SCI ADV：对线粒体异质性的核全基因组关联分析发现mtDNA复制、产能等相关基因附近的核变...
编译：微科盟Nelly，编辑：微科盟悄咪咪、江舜尧。
微科盟原创微文，欢迎转发转载。
目前，关于核基因组在维持线粒体基因组（mtDNA）稳定性中扮演什么样的角色尚不清楚。mtDNA序列变异可以异质性状态存在，表明具有不同序列的细胞器基因组可以共存。损害线粒体产能的异质性变异会导致疾病发生，异质性状态其本身也是有害的。然而，线粒体异质性可能为新型线粒体单倍群的出现提供过渡状态。研究使用来自982,072名欧洲血统个体的全基因组基因分型数据，来评估线粒体异质性中的变异情况，并鉴定影响该变异的核基因组区域。年龄、性别和线粒体单倍群与异质性程度相关联。GWSA鉴定出20个超出全基因组显著性与异质性相关的基因座。其包括一个与线粒体转录因子A（TFAM）重叠的区域，该区域在mtDNA组装、复制和转录中具有多种作用。这些结果表示，线粒体异质性具有可遗传的核组分。
论文ID
原名：Nuclear genome–wide associations with mitochondrial heteroplasmy
译名：线粒体异质性的核全基因组关联分析
期刊：Science Advances
IF：14.136
发表时间：2021年3月
通讯作者：David Hinds；Neal Sondheimer
通讯作者单位：23andMe基因检测公司；多伦多大学儿科学与分子遗传学院
DOI号：10.1126/sciadv.abe7520
实验设计
结果
1 数据集和分析策略
研究使用的全基因组基因分型数据来自982,072名欧洲血统个体的唾液样本，这些个体是个人基因组学和生物技术公司—23andMe公司研究项目的参与者。MtDNA基因分型密集，检测出3287个单核苷酸多态性（SNP）。应用质量控制（QC）检测，去除由于缺乏杂交或等位基因不准确区分而导致基因分型不佳的检测，最终评估出326个mtDNA SNP（附图1）。用较小等位基因强度与该位置总强度的比值作为MtHz值，因此可能的最大异质性值为0.5。研究使用母体-后代组合来检测母体异质性值>5%的成对点（附图S1；n=28,963对）。正如对于所遗传的异质性所预期的那样，母体-后代值存在相关性。
表1|按平均MtHz四分位数划分的受试者性别和年龄群体表征。
图1|在23andMe v4阵列上评估线粒体SNP。该阵列包括整个线粒体基因组的3287个位置。这些位置经过修剪以保持双等位基因变异，在整个测试群体中检出率> 99%，次要等位基因频率（MAF）>0.001。若个体参与者在单个SNP的检出强度（LRR、log2 R比）比该位置所有个体的平均强度低>3个标准差（SD），将该检出移除。如果一个位置有>1%的样本达不到这些标准，则在所有个体中都将该位置移除。（A）在两个同质簇之间的异质样本和少量强度低于平均值（最低线）3个SD（个体检出已排除）的样本的性能良好分析。（A）（绿色）和（B）（红色）等位基因检测均显示SD线。那些鉴定出差同质簇的过度异质性的检测已被移除。（B）性能不佳的检测，该检测中数据集中>10%个体的MtHz >20% [B等位基因频率（BAF）值，0.2到0.8]。所有个体关于此位置的数据都从分析中移除。
MtHz广泛存在，所有被检测位置都存在>25% MtHz的个体（附图2A）。由于本研究的方法没有对线粒体中所有位置进行预估，因此将每个个体的MtHz进行定量，作为检测的326个位置的平均值。所评估的所有个体的平均MtHz为0.00744（四分位距=0.0046至0.012；附图2B）。值得注意的是，平均异质性值不是由线粒体位置的子集驱动的，因为每个位置的个体间平均异质性是紧密分布的（附图2C）。
图2|MtHz数据定量。（A）MtHz按群体中每个线粒体位置评估的百分位数绘制。每个位置处均观察到高度异质性个体。（B）数据集中个体的平均线粒体异质性，在评估位置上取平均值。（C）QC子集中检测出的326个线粒体位置的平均MtHz值分布。（D）平均MtHz的自然对数转换，异常值去除率为99.5%。
2 年龄、性别和线粒体单倍群对MtHz的影响
为检测遗传和非遗传因素对平均MtHz的影响，研究者计算了自然对数转换的平均MtHz值，在99.5%处进行缩尾处理，以消除异常值的影响（附图2D）。建立初始模型以计算对平均MtHz的影响（表2）。纯合常染色体SNP等位基因强度方差（“常染色体方差”）与MtHz呈强正相关。这是预料之中的，因为在纯合位置具有高常染色体方差的样本可能具有表明样本噪点更大的技术特征。这种噪点使MtHz被高估，在全基因组关联研究（GWAS）中使用常染色体方差可校正样本质量的微小差异。
表2|协变量对平均MtHz的影响。
每个群体至少有10,000个个体，针对存在于所有群体中的19个最常见单倍群（即总群体的约1%），评估线粒体单倍群的影响（表2）。不同单倍群之间的平均MtHz不同。H1是数据集中最常见的单倍群（n=203,003；20.7%）。尽管许多单倍群对MtHz的影响是显著的，但MtHz差异的幅度很小，具有T1单倍群的个体（n=21,921；2.2%）的平均MtHz仅比H1个体大0.0010（在空模型中，P=1.1×10−248）。群体中的平均MtHz随年龄增加而显著下降。此外，女性的平均MtHz显著低于男性。将年龄和性别作为协变量包含在GWAS模型中。
3消除线粒体假基因
对平均MtHz进行GWAS。初始时，在调整基因组膨胀因子λ=1.077后，37个基因座达到了全基因组显著性，P<5×10−8（附图3和附图S2）。然而，核mtDNA片段（NUMTs）的存在是重要的潜在混杂因素。NUMTs包括已逆转录到核基因组中的mtDNA的部分的、片段的、复合的或完整的拷贝。如果NUMT序列包含与个体真实mtDNA序列不同的等位基因，那么表面线粒体探针与来自NUMTs的核DNA的直接杂交将造成MtHz被高估。因此，NUMT有可能在处于连锁不平衡（LD）的SNP处产生假阳性GWAS信号。
图3|平均MtHz初始关联的曼哈顿图。
例如rs1951197（AKAP6附近）所示，它与对数平均MtHz相关联。该位置与NUMT474重叠，后者在欧洲人中呈现多态性（esv3633987）（附图S3）。NUMT474与修订的剑桥参考序列（rCRS）mtDNA（NC_012920.1）的第5583至6606位有93%一致。从MtHz的定量中去除mt.5583和mt.6606之间的任何位置后，通过重复分析来测试rs1951197关联的有效性（表S1）。正如预期的那样，rs1951197和平均MtHz之间的关联不再显著。使用类似方法，本研究还消除了rs7728823（与NUMT228重叠）和rs571982832（比对至最近鉴定出的复合物NUMT）的关联。
因为并非所有NUMT都是已知的或精确定位的，研究者试图确定其他由未知NUMT驱动的明显关联的基因座。将线粒体位置分为三组后，分别估算平均MtHz。如果NUMT驱动了伪关联，则推断在核位置与这些线粒体位置子集的关联上将显示出显著差异。当异质性分析基于子区域时，当关联降低至全基因组显著性以下时，将SNP从后续分析中移除。这鉴定出四个另外的基因座，它们的关联具有很强的区域依赖性（表S1），也将这些基因座从随后的分析中移除。
4 全基因组关联鉴定与MtHz关联的基因座
在排除NUMT依赖性位置后，30个基因座最初有至少一个具有全基因组显著性的位置。在其中的10个基因座上，只有一个单独推定的SNP超过全基因组显著性，为了避免估算出支持不足的关联，从后续分析中将这10个基因座移除，于是得到最后与平均MtHz相关的20个基因座（表3和附图S4）。从阵列估算，MtHz观测尺度遗传力为0.65%（SE=0.1%，平均χ2=1.1914，截距=1.032）。20个基因座占观察遗传力的32%（表S2）。通过在模型中包含单倍群的情况下重新检测峰值关联，确认了这些关联不依赖于单倍群（表S3）。
表3|与平均MtHz关联的基因座。
在这些候选基因中，三个基因座临近四个具有明确线粒体功能和已发表的线粒体定位的基因（TFAM、TWNK、MRPL43和NDUFS4）。在rs1049432处发现了最强的关联（P=1.7×10−223），rs1049432邻近线粒体转录因子A（TFAM）。TFAM是个有意思的候选基因，它与mtDNA有许多相互作用。TFAM最初被定性为线粒体转录因子，但随着时间的推移，TFAM在mtDNA的维持和包装中的多种作用被发现。先前的组织特异性表达数量性状基因座分析表明，rs1049432处的T等位基因（与升高的MtHz关联）与较低的TFAM表达显著关联（附图4和附图S5）。rs1049432是非编码的，与rs1937（p.Ser12Thr）（r2=0.44）关联，后者是TFAM中最常见的编码SNP（MAF=0.08）。
先前的研究表明，rs11006126在两项研究中与唾液和血液样本中的mtDNA拷贝数关联，rs11006126与在rs1049432处鉴定出的TFAM峰值接近且处于强LD，并包含在可信集合（rs11006126）中。研究者对mtDNA拷贝数的改变是否可能混淆在MtHz上的发现进行检测。使用来自基因分型阵列的强度数据对相对mtDNA拷贝数表型进行定量。重新检测MtHz空模型以及平均MtHz与来自GWAS的候选基因座之间的关联，发现所鉴定的关联与mtDNA拷贝数无关（表S4和表S5）。
除了TFAM，其他邻近相关SNP的基因有可能基于其已知的细胞活动影响MtHz。这包括靠近C10orf42的rs58678340，它编码线粒体解旋酶Twinkle（TWNK）和线粒体核糖体蛋白MRPL43。值得注意的是，该SNP与TWNK编码变异rs17113613处于强LD中（r2=0.84；p.Val368Ile）。TWNK是mtDNA复制机制的一部分，遗传缺陷会导致mtDNA耗竭或缺失综合征。TWNK可防止出现mtDNA变异，这暗示了TWNK活性变化与MtHz相关联的潜在机制。
MtHz也与rs10063311相关，rs10063311邻近NDUFS4，NDUFS4是电子传输链复合体I的亚基。尽管电子传递链的运行与mtDNA的完整性之间不存在明确的联系，但之前的研究已经表明复合物V亚基缺失会影响mtDNA的数量。
除此之外，尽管CLEC16A和PRKAB1编码的蛋白质不存在线粒体定位，但是它们与mtDNA复制的保真度和稳定性具有潜在的功能联系。从1型糖尿病的GWAS中鉴定出的CLEC16A，通过与NDRP1和PARKIN的相互作用调节线粒体自噬。PRKAB1编码一磷酸腺苷活化蛋白激酶的亚基，该亚基与促进线粒体内生物发生和能量产生的一系列途径有关。这表明对线粒体QC很重要的细胞质和核过程在MtHz的调节中发挥作用。
几个在免疫中起作用的基因（包括HLA-DQB1、IL1RN、IFNL4和FUT2）位于与MtHz关联的SNP的近端。免疫系统功能的变化可能对MtHz有直接影响，但这些变异也有可能影响DNA样本中细胞类型之间的比例。需要对其他样本类型进行进一步检测以确认这种关联。
5 基于通路和基因的分析
使用基因集分析（MAGMA [v1.07]）对GWAS结果进行评估，鉴定出尿酸盐转运的基因本体生物学途径富含关联（P=6.7×10−7），在对15,484条通路进行Bonferroni校正后，这种关联仍然显著。该基因集中有五个基因（SLC22A13、SLC2A9、SLC17A1、SLC17A3和SLC22A12），但没有一个基因位于符合GWAS显著性的基因座中。SLC17A1和SLC17A3的基于基因的结果均符合显着性标准（分别为P=4.1×10−11和8.4×10−9）并且在基因组中是相邻的，而该基因集的其他三个成员没有在对多个基因的校正中存在，表明信号可能由单个基因座驱动。此外，尿酸转运途径与控制MtHz的机制之间的相似性并不直观。
6 PheWAS分析
使用来自23andMe的全表型组关联研究（PheWAS），检测19个SNP与1123个性状的关联。在对评估的性状和SNP应用Bonferroni校正后，287个SNP-表型对符合显著性标准（表S6）。研究者重点关注与MtHz最密切相关的SNP的关联。对于rs1049432（TFAM），与较高异质性相关的T等位基因与多囊卵巢综合征风险降低相关（61,181例病例和839,824例对照；优势比=0.96，P=6×10−7），这种关联先前未在10,074例多囊卵巢综合征的meta-GWAS中被鉴定出。
讨论
本研究使用通过阵列得到的大型Mt变异的基因分型群组，评估了影响MtHz水平的参数。本研究工作存在以下几个限制。首先，先前通过与等位基因特异性定量聚合酶链式反应进行比较，已确认阵列可以用于MtHz，但本研究中使用的这个阵列进行并未经过确认，并且MtHz的估值中可能存在噪点。本研究的另一个限制是，在估算MtHz时并未评估所有线粒体位置，而是侧重于选择高检出率和可观的BAF的SNP子集。最后，所使用的组织类型（唾液）可能无法代表影响mtDNA的性状的所有组织。
本研究鉴定出年龄、性别、线粒体单倍群和MtHz值之间的关联。MtHz随着年龄的增长而降低。这一发现出乎意料，因为之前的研究鉴定出MtHz随着年龄的增长而增加。来自弗雷明汉心脏研究的356名个体的血液DNA测序发现，随着年龄的增长，基因组中多个位置的MtHz升高。同样，一项使用白细胞DNA对2077名撒丁岛人的研究也发现，随着年龄的增长，MtHz增加，拷贝数减少。对于这种差异，一种可能的解释是，所研究的组织类型可能会影响MtHz和衰老的动态，本研究的观察可能仅特定于唾液。另一种可能性是，本研究对大量位置的MtHz的定量可能与之前的研究不同，之前的研究在分析中使用了异质性峰值，或比本研究检测的mtDNA位点少。
MtHz的变异并不经常按性别进行评估。一项对235名具有致病性mt.3243A>G变异异质性患者尿液样本的研究发现，男性尿液样本中的MtHz高于女性。然而，使用1035名没有线粒体疾病的个体的白细胞DNA进行的测序研究没有发现显著的性别差异。与年龄类似，这种差异可能反应了组织的特异性。此外，本研究具有更高的能力来检测MtHz中较小的年龄和性别依赖性差异。
从线粒体遗传的角度来看，女性MtHz通常低于男性的这种可能性很有趣，因为男性MtHz的影响应该仅限于个体，而不具有将MtHz传播给后代的风险。然而，体细胞异质性的研究不能轻易扩展到对雌性生殖系的理解。
本研究鉴定出20个与MtHz水平关联的基因座。其中两个基因座临近DNA结合TFAM和线粒体解旋酶TWNK，这两个基因编码直接参与mtDNA复制的蛋白质。表达数据表明TFAM处的变异rs1049432与许多组织中TFAM表达的差异有关。一个明显的问题是，TFAM的单体剂量不足是否与mtDNA MtHz增加有关。最近在小鼠模型中的研究显示，TFAM表达减少会导致致病性MtHz水平下降。值得注意的是，这项研究表明，无论异质性如何变化，致病性变异造成的后果在较低拷贝数下比在高拷贝数下更大，这表明野生型等位基因的拷贝数可能是疾病表型的控制特征。在TFAM致病性变异患者的单一报告中，罕见错义变异的纯合性导致TFAM缺失、mtDNA耗竭和一种线粒体表型，但同样，这种TFAM缺失对线粒体序列保真度的影响未知。
TFAM是mtDNA的多功能结合蛋白。顾名思义，它最初被鉴定出在线粒体转录中发挥作用，但TFAM在mtDNA复制和基因组的整体压缩中也具有多种作用。在两项研究中，一个SNP与唾液和血液样本中的mtDNA拷贝数均相关，这个SNP邻近本研究鉴定出的TFAM峰值（r2=0.96）并具有强LD，并且包含在本研究的可信集中（rs11006126）。对拷贝数的观测影响的方向性很有趣，因为在这些研究中，等位基因关联的拷贝数越大，相关联地，MtHz增加越多。之前观察到的TFAM对mtDNA拷贝数的影响可能与本研究中TFAM对异质性影响的发现有关。蔡等人观察到不稳定二核苷酸重复mt.514到mt.523处的MtHz也与mtDNA拷贝数相关，并且可能不足以在其他基因座观察到这种关联。
在本研究的分析中，研究者试图排除NUMT的影响。如果核基因组内的线粒体序列块的序列与个体的mtDNA不同，那么他们可以使异质性估值变大。这反应出本研究存在一个潜在的问题，因为NUMT与邻近核基因座处于LD中。本研究鉴定出驱动与核基因座的错误关联的已知NUMT，但也证实额外的基因座可能与未知NUMT邻近。鉴定这些错误关联的策略可能对包含全长线粒体基因组的非常大的“mega-NUMT”无效；然而，它们在群体中显然不常见，因此它们的影响存在限制。
MtHz的状态是细胞器基因组的一个重要特性，它影响新的多态性组合的出现，并在由致病变异引起的疾病的外显率和严重程度中起着至关重要的作用。异质变异既可以发生在体细胞水平，也可以通过种系传播改变。在本研究中，研究者评估了一个大样本，以鉴定影响异质性的核编码变异和附近基因。结果表明，与mtDNA复制所需的基因（TFAM和TWNK）以及其他与线粒体产能和QC相关的基因（CLEC16A和PRKAB1）邻近的核变异与MtHz相关，这可能是基于这些基因在维持线粒体内复制保真度方面的作用。