iMetaOmics | 魏来/贾慧珏/何明光-多组学助力揭示塑造转录组的基因型-微生物组相互作用...

点击蓝字 关注我们

多组学整合揭示塑造结膜转录组的基因型-微生物组相互作用

研究论文

● 原文链接DOI: https://doi.org/10.1002/imo2.37

● 2024年10月14日，广州医科大学附属第二医院魏来、复旦大学贾慧珏、香港理工大学何明光等团队在iMetaOmics在线发表了题为“Multiomics integration unravels genotype-microbiome interactions shaping the conjunctival transcriptome”的文章。

● 本研究整合了基因组、转录组和微生物组数据，分析了基因组和微生物组及其相互作用与结膜转录组之间的关系，并开发了R包MicroGenix用于在多组学数据中识别潜在的基因型-微生物组相互作用，助力揭示分子表型塑造的机制。

●  第一作者：梁桥星、林国旺

●  通讯作者：魏来（weil9@mail.sysu.edu.cn）、贾慧珏（huijue_jia@fudan.edu.cn）、何明光（mingguang_he@yahoo.com）

●  合作作者：丁小虎、邹斌、刘小敏、李晶、张宇欣、文小凤、梁凌毅、贝锦新

●  主要单位：复旦大学生命科学学院、粤港澳大湾区精准医学研究院（广州）、中山大学中山眼科中心、南方医科大学珠江医院、广州医科大学附属第二医院等

 亮 点

●  宿主基因型和微生物组对宿主分子表型（例如基因表达）的影响可能存在相互依赖关系，但这些交互效应在很大程度上仍然未知；

●  我们通过分析来自120对双生子的多组学数据，探究了基因型和微生物组与结膜转录组的关联；

●  将基因型-微生物组交互项纳入基因表达模型，可提高大量基因的表达水平预测准确性，其中包括许多编码细胞粘附分子的基因;

●  我们开发了R语言工具包MicroGenix，用于在多组学数据中筛选可能参与塑造分子表型的基因型-微生物组相互作用。

摘  要

基因表达是由基因型和环境相互作用塑造的分子表型。微生物组是宿主的关键环境暴露。然而，塑造宿主转录组的基因型-微生物组交互作用（genotype-microbiome interaction, GMI）在很大程度上尚未被探索。在这里，我们通过整合分析来自120对双生子的结膜多组学数据来研究GMI。我们鉴定了与5946个基因相关的272,972个表达数量性状位点，以及基因表达与微生物丰度之间的241,073个基因型控制条件下的相关性。我们开发了一种建模方法并将其命名为MicroGenix，用于在宿主基因组、转录组和微生物组数据中筛选GMI，并通过基于基因的关联测试识别与疾病相关的GMI。我们将MicroGenix应用于结膜数据集，发现将GMI纳入基因表达预测模型显著提高了预测准确性。准确性提高的基因中包括许多编码细胞粘附分子的基因。我们进一步使用MicroGenix从结膜宏基因组预测转录组，并识别了与眼表疾病相关的GMI。我们的工作为研究结膜处宿主与微生物组的相互作用提供了资源，并提供了一种使用多组学数据研究GMI的计算方法。

视频解读

Bilibili：https://www.bilibili.com/video/BV1iSy3Y5ETM/

Youtube：https://youtu.be/JPJ_2yo6ysA

中文翻译、PPT、中/英文视频解读等扩展资料下载

请访问期刊官网：http://www.imeta.science/imetaomics/

全文解读

引  言

基因型-环境交互作用（genotype-environment interaction, GEI）是指基因型和环境对表型的影响因彼此而异的现象。GEI被认为是表型变异（包括对复杂疾病的易感性）的关键机制。微生物组可以被视为宿主的环境暴露，并可能与宿主基因型相互作用以塑造宿主表型。为简洁起见，我们将这种交互效应称为基因型-微生物组交互作用（genotype-microbiome interaction, GMI）。基因表达是一种分子表型，可能由宿主基因型和微生物组共同调节。例如，一项基于原代人结肠上皮细胞和肠道微生物组提取物共培养的研究提供了宿主基因型和微生物组在调节宿主转录应答中相互作用的证据。此外，微生物组应答表达数量性状位点（microbiome response expression quantitative trait locus）的概念已被提出，它是指宿主基因表达对特定微生物组暴露的响应的eQTL。然而，宿主基因表达调控中具体的GMI在很大程度上仍然未知。

研究GMI的实验方法仅限于针对少数基因和微生物，主要基于体外或动物模型研究，而多组学数据整合分析在生成GMI假设方面具有较大潜力，并且在人类研究中是可行的，特别是对于基因表达等分子表型。在GMI研究中，考虑分子表型至关重要，因为它们在复杂性状的出现中起着基础作用，并可以帮助阐明潜在的分子机制。现有的计算方法评估基因组-转录组或微生物组-转录组的关联，而没有明确建模两者加性效应或GMI。在多元回归中，交互项捕捉两个变量之间的乘法效应，使得可以建模一个变量对结果的影响可能取决于另一个变量的值。乘法交互项已被用于在模型中表示基因型和环境在表型塑造中的交互作用。然而，这种方法尚未被用于在多组学数据中研究GMI对分子表型的影响。

在这项研究中，我们对广州双生子研究中的120对双生子的基因组、结膜转录组和结膜微生物组数据进行了整合分析。我们首先确定了宿主基因型与结膜转录组之间的关联，以及结膜微生物组与结膜转录组之间的关联。我们进一步开发了MicroGenix方法，并将其应用于结膜数据集，以筛选潜在的GMI。MicroGenix方法将基因表达建模为宿主基因型、微生物组及其交互作用的函数。这种方法还允许识别与疾病等复杂性状相关的GMI。

结  果

宿主遗传变异与结膜转录组的关联

迄今为止，大多数粘膜组织（包括结膜）的表达数量性状位点（eQTL）尚未确定。我们通过进行全基因组eQTL定位来表征遗传因素对结膜转录组的影响。分析仅纳入每对双生子中的一个个体，并平均左右眼的基因表达水平，因此分析数据集由不相关个体组成。为便于与GTEx数据集比较，我们采用与GTEx项目相同的分析方法和显著性阈值。我们鉴定了属于5946个基因的187,241个顺式eQTL（FDR < 0.05，表S1）。对于每个eGene（即具有eQTL的基因），检测到的eVariant（相应的显著变异）的数量的中位数为五个，73.4%的eGene有一个以上的eQTL（图1A）。随后，我们将结膜eQTL与GTEx项目中其他组织的eQTL进行比较。共695个eGene为结膜独有，例如IL5和FCRL4（图1B；表S1）。甲状腺与结膜共有的eGene数量最多，其次是日晒皮肤、胫神经、非日晒皮肤和食管粘膜（图1C）。对于结膜与其他组织共有的eVariant-eGene对，不同组织与结膜中作用方向相同的eVariant的百分比在90.77%到98.66%之间（图S1；表S1）。我们鉴定了编码细胞因子（例如IL5、TNF和VEGFC）和趋化因子（例如CCL22、CCL28、CXCL14、CXCL17和XCL1）及其受体（例如IL6R、IL12RB2、IL18R1、IL20RB、TNFRSF1A、CCR1和CCR2）的基因的顺式eQTL。这些结果还为研究结膜免疫性疾病的遗传基础提供了资源。

图1. 结膜eQTL定位分析

（A）具有一个或多个顺式eQTL的结膜基因数量。（B）维恩图表示结膜eGene（具有eQTL的基因）与GTEx项目所包括的所有组织的eGene的交集。（C）与GTEx项目中每个组织共有的eGene和eVariant-eGene对（eGene和相应的显著变异）的数量。

微生物组与结膜转录组的关联

接下来，我们在同卵双生子中分析了结膜基因表达与微生物丰度之间的关联，从而控制宿主遗传变异的影响。为了评估环境对基因表达表型的影响，我们确定了双生子对之间表达水平的倍数变化（FC）。我们观察到转录组中平均log₂转换的FC呈双峰分布（图2A）。具体来说，20.0%的基因显示FC > 2，而36.6%的基因显示FC < 1。受环境影响最大的基因（定义为平均FC > 2且中位数FC > 2）主要涉及免疫反应，如免疫球蛋白、细胞粘附分子、细胞因子和趋化因子（图2B）。

为了研究微生物组作为环境因素如何解释上述受环境影响较大的基因的双生子内表达水平变异，我们对双生子对中两个个体的基因表达的FC与β-多样性以及物种丰度的FC进行了相关性分析。我们鉴定了两个基因的FC与β-多样性呈正相关，分别为LCN1（Spearman ρ = 0.39，FDR = 0.0047，图2C）和LACRT（Spearman ρ  = 0.32，FDR = 0.059，图2D）。此外，我们鉴定了基因表达与物种丰度的FC之间的60个正相关关系（Spearman ρ > 0.2，FDR < 0.1）（图2E；表S2）。具体来说，编码免疫球蛋白（IGHA1、IGHG1、IGHG3、IGHM、IGKC）的基因，在双生子内部表达变化最大，其与马拉色菌属（Malassezia）的物种（Spearman ρ  > 0.2，FDR < 0.1，图2E）呈正相关。

此外，我们进行了双生子对内结膜基因表达与微生物组物种丰度之间的关联分析（图2F）。我们拟合了线性混合模型，将宿主基因对微生物物种进行回归，将家庭作为随机截距校正。总共鉴定了241,073个显著关联（FDR < 0.1，表S3）。我们针对每个物种对显示正相关的宿主基因（表S3）进行了通路富集分析。结果表明，与免疫反应相关的通路上调，包括FcγR介导的吞噬作用、趋化因子信号传导、MAPK信号传导、TNF信号传导和B细胞受体信号传导（图S2）。

图2. 结膜基因表达与微生物丰度在双生子内的相关性

（A）所有同卵双生子基因表达水平的双生子内倍数变化（FC）的分布。条形图的颜色表示log₂转换的FC范围。图中标注了log₂转换的FC > 2或 < 1的基因数量。（B）同卵双生子中表达水平平均log₂转换FC最高的前50个基因。数据表示为平均值±标准误差。误差条表示标准误差。（C, D）同卵双生子内微生物Bray-Curtis距离与LCN1（C）和LACRT（D）表达的log₂转换FC之间的Spearman秩相关性。（E）表达受环境影响的基因（中位数和平均FC > 2）的双生子内FC与微生物丰度的双生子内FC之间的Spearman秩相关性。标有星号的相关性满足Spearman ρ > 0.2和FDR < 0.1。（F）曼哈顿图展示了结膜基因表达与微生物丰度在双生子内的全转录组关联。线性混合模型的p值使用Satterthwaite方法计算。

将GMI纳入模型可提高基因表达预测的准确性

为了鉴定可能参与塑造结膜转录组的GMI，我们将宿主基因表达建模为基因型和微生物组的主要效应及其交互效应的组合。对于主要效应，我们包括了位于基因起始和结束位置± 1 Mb的单核苷酸多态性（SNP）的剂量，以及所有微生物物种的丰度。应用了两阶段拟合策略来确定要包括作为预测因子的SNP-微生物交互项。具体来说，首先拟合一个没有交互作用的模型，以选择重要的SNP和微生物预测因子，定义为在正则化模型中有非零系数或在标准线性回归模型中有名义显著性。然后通过将每个选定的SNP的剂量和微生物的丰度相乘来构建交互项。使用所有主要效应预测因子（不进行重要性筛选）加上交互项来拟合包含交互的模型。预测因子的高维性和多重共线性对标准线性回归构成挑战。对于包含和不包含交互的两种模型，通过决定系数R²进行评估（配对t检验，p < 2.2×10^-16；图S3），弹性网络在预测准确性方面优于多元线性回归和线性混合模型，线性混合模型与多元线性回归相比改进不显著（配对t检验，p > 0.05）。因此，我们在后续分析中使用了弹性网络。

在分析中包括的所有蛋白质编码基因（n = 14,197）中，18.9%的基因（n = 2681）在没有交互的模型中至少有一个SNP和一个微生物物种具有非零系数。在这些基因中，92.3%（n = 2475）至少有一个交互项具有非零系数。在包含交互的模型中预测良好的基因包括编码HLA类II分子的基因，如HLA-DQB1、HLA-DQA2和HLA-DOA（R² > 0.2，图S4）。值得注意的是，包含交互的模型的预测准确性高于不包含交互的模型（配对t检验，p < 2.2×10^-16；图3A），并且2158个基因在包含交互的模型中显示出准确性提高（表S4）。准确性提高的基因在细胞-基质粘附这一生物过程富集（图3B）。在细胞粘附分子中，ANTXR1（ANTXR细胞粘附分子1）的准确性提升最多，R²从0.012提升到0.194（图3C、D；表S4）。准确性提高的细胞粘附分子编码基因（表S4）的SNP-微生物交互形成了一个网络，包括217个顺式eQTL和35个与基因表达相关的微生物物种（后者在双生子对内分析中所识别）。例如，ANTXR1的表达水平与rs6760455和Corynebacterium ureicelerivorans之间的交互作用相关（图3E）。我们检查了2158个基因的最重要的交互项（定义为系数的绝对值最大的项）涉及的物种（图3F；表S4）。Porphyromonas gingivalis参与最多基因的最重要交互项，其中LRFN2是R²提升最多的基因。

图3. 基于基因型、微生物组及其相互作用的基因表达预测

（A）比较包含和不包含SNP-微生物交互项的弹性网络模型的预测准确性。散点图显示了两种模型类型的决定系数（R²），每个点代表一个基因。呈现了每次交叉验证中的平均R²。p值是使用配对t检验计算的。密度图说明了基因准确性提高的分布，定义为包含交互的模型中的R²减去不包含交互的模型中的R²。（B）显示在包含交互的模型中预测准确性提高的基因的富集通路。点的大小表示基因数量，点的颜色表示通路中预测准确性提高的基因比例。（C, D）使用包含和不包含交互的模型预测的ANTXR1表达水平与观测值的对比。点的颜色表示交叉验证中每个折叠的测试集。图中标注多个折叠的平均R²。线条表示使用线性模型的最佳拟合，阴影区域表示95%置信区间。（E）网络图展示编码细胞粘附分子的基因（在包含交互的模型中准确性提高）与模型中交互项所涉及的微生物。（F）参与构成最多基因的最重要交互项（弹性网络中系数绝对值最大的交互项）的前30种微生物物种。

识别与复杂性状相关的GMI

为了便于使用我们的多组学数据筛选GMI的方法，我们开发了一个名为MicroGenix的R包，它由三个模块组成：使用基因型、微生物组及其交互作用对参考转录组进行训练；使用基因型和微生物组数据预测基因表达水平；以及将预测的表达与复杂性状相关联（图4A）。预测的表达水平可以被视为由基因型、微生物组及其交互作用解释的转录变异的估计。将解释的基因表达组分与复杂性状相关联有助于优先考虑GMI，这些GMI在关联基因的建模中发挥作用。基于基因的关联测试使用逻辑和线性回归分析定性和定量性状，使用单个基因的预测表达水平作为变量。

我们通过识别可能在过敏性结膜炎中起作用的GMI来举例说明MicroGenix的应用。为此，我们预测了一个包含健康个体和过敏性结膜炎或干眼病患者样本的宏基因组数据集中的结膜基因表达水平。我们使用根据纳入分析的总基因数量的Bonferroni校正来确定显著关联阈值（p < 2.4×10^-5）。我们鉴定了148个与过敏性结膜炎相关的基因，其中39个基因显示出疾病特异性，即不与干眼病相关联。此外，我们使用不包含交互的模型、仅用SNP或微生物预测因子拟合的模型进行了相同的分析。39个基因中有22个是仅通过包含交互的模型鉴定的，并且所有这些基因在包含交互的模型中的R²都高于其他三种模型（表S5）。值得注意的是，22个基因中有11个与过敏性结膜炎呈正相关（图4B），其中5个位于先前已与过敏相关联的遗传区域，包括PTTG1IP、SCIN、TMEM18、CFLAR和LAMC2，其中LAMC2是一种细胞粘附分子。其他参与免疫反应的基因包括CD276和CXCR2。CD276在过敏性哮喘儿童中升高，并调节Th2/Th17细胞分化。基因表达预测模型将这些基因与涉及24种微生物物种和31个SNP的SNP-微生物交互作用联系起来（图4C；表S6）。这些发现为进一步研究眼表疾病中的GMI提供了基础。

图4. MicroGenix用于筛选可能与疾病发病机制有关的SNP-微生物交互作用

（A）MicroGenix方法包括使用参考数据集训练基因表达预测模型，这些数据集包括宿主基因型、微生物组和转录组数据。这一过程生成了一个包含预测模型和SNP-微生物交互的数据库。这些模型将被应用在包括宿主基因型和微生物组数据的新数据集上，以预测基因表达。随后进行关联测试，将预测的基因表达水平与特定目标性状相关联。这有助于优先考虑基因，并在数据库中识别关键的SNP-微生物交互作用。G×M表示基因型（G）和微生物组（M）之间的交互作用，捕捉两者之间的乘法效应。宿主基因表达（Y）被建模为截距（β₀）、基因型（β₁）和微生物组（β₂）的主要效应以及交互效应（β₃）的线性组合。（B）与过敏性结膜炎正向且特异性相关的11个基因，在健康个体和过敏性结膜炎患者中的预测基因表达水平。（C）过敏性结膜炎相关基因的弹性网络模型中，拥有非零系数的SNP-微生物交互项的系数绝对值。对于具有多个交互项的基因-微生物对，系数绝对值最大的SNP在图中展示。

讨  论

生物体的遗传构成与其环境之间的相互作用在表型塑造过程中起着重要作用。一个引人注目的假设是宿主基因型和微生物组的相互作用塑造了转录组，后者可以被视为一种分子表型。本研究建立了结膜基因表达与宿主基因型或微生物组之间关联的资源。利用这一资源，我们开发了一种方法，使用多组学数据探索可能参与塑造转录组的GMI。

在以往的研究中，我们已经表征了眼表疾病（包括过敏性结膜炎、干眼病和慢性眼移植物抗宿主病）中结膜微生物组的失衡。我们分析了结膜微生物组和转录组之间的关联，发现马拉色菌属（一种在过敏性结膜炎患者结膜中富集的真菌物种）与编码各种免疫球蛋白的基因表达呈正相关。M. furfur与CD79A、MS4A1和FCRL5的正相关为B细胞的功能角色提供了进一步的证据。值得注意的是，马拉色菌已被发现通过引发IgE反应在过敏中发挥促炎作用。此外，我们的发现揭示了多种微生物物种的丰度与编码两种主要泪液蛋白的LCN1和LACRT的表达呈正相关，这些蛋白质维持稳态。这些蛋白质在干眼病患者中被认为下调。鉴于干眼病患者的结膜微生物组以共生菌的减少为特征，确定共生菌是否在维持这些蛋白质的表达中发挥作用，从而预防干眼病的发作，将是有利的。

几种经常出现在多个基因最重要的交互项中的微生物物种已与眼表疾病相关联。例如，Porphyromonas gingivalis已被证明在眼移植物抗宿主病中富集，它出现在LRFN2和CCL28的最重要交互项中。LRFN2是一种调节T细胞反应的细胞粘附分子，对移植物抗宿主病的病理生理学至关重要。CCL28是一种趋化因子，促进T细胞在粘膜的归巢，并表现出抗菌活性。此外，M. furfur出现在SERPINE1的交互项中，SERPINE1是一种与过敏性疾病遗传相关并广泛参与过敏发病机制的细胞粘附分子。其他值得注意的物种包括与干眼病相关的Rothia mucilaginosa和Neisseria subflava。鉴定出的潜在GMI可能代表眼表疾病遗传和微生物发病机制中的关键交汇点。

与单个SNP或微生物与基因之间的成对关联不同，使用多个SNP和微生物对基因表达进行建模捕捉了表型决定的多变量性质。与该领域应用正则化的常见做法一致，我们发现弹性网络方法在大多数基因上的表现优于标准多变量线性回归。弹性网络对处理高维性和多重共线性有效。它通过交叉验证进行变量选择，并倾向于将相关性高的预测因子同时包含或排除。这些特点避免了人为设定构建SNP-微生物交互项的显著性阈值和任意丢弃相关SNP和微生物，减小了错过重要发现的可能性，因而优于一般的线性回归。然而，我们承认对于一小部分基因，一般线性回归模型的准确性更高。为了保持我们分析的一致性，我们对所有基因使用了弹性网络方法。结合使用不同建模方法拟合的预测模型可能有助于揭示单一方法错过的发现。

我们发现在模型中纳入GMI可以提高基因表达的预测准确性，这与GEI在许多性状的表型变异中解释很大比例的观察结果一致。值得注意的是，准确性提高的基因中，包括许多编码细胞粘附分子的基因。这些分子在免疫系统的各个方面发挥关键作用，包括白细胞迁移、炎症和淋巴细胞发育。基于结肠上皮细胞和肠道微生物组提取物共培养的研究，鉴定了遗传变异影响微生物组暴露反应的基因。这些基因涉及细胞粘附和迁移，包括在结直肠癌中表达增加的LASP1。尽管该共培养研究没有探索参与交互作用的特定微生物类型，我们的数据表明Fusobacterium nucleatum，一种已知促进结直肠癌的细菌，参与了LASP1表达模型的交互项。

MicroGenix方法可以利用多组学数据并行筛选数千种分子表型中潜在的GMI。除了基因表达，这种方法可以应用于各种分子过程，包括剪接和甲基化。此外，如果有足够的宿主读取量，可以在宏基因组测序数据上进行基因分型。这使得同时提取基因型和微生物组数据成为可能，这些数据是MicroGenix预测模块的输入，因而有可能从宏基因组数据推断宿主转录组。一旦模型通过一次性训练，它们就可以被用来鉴定与疾病和其他复杂性状相关的GMI。然而，在解释预测值和关联结果时需要谨慎，因为基因型和微生物组对分子表型的影响在不同的生理条件下可能会变化。

预测过程可以被视为估计由基因型、微生物组和GMI决定的分子表型组分。因此，预测值提供了超出真实表达水平的额外信息层面。此外，当应用于包括分子表型测量的疾病队列的多组学数据时，这种方法可以用来鉴定疾病特异性的GMI。这种分析将分子表型的真实值与潜在的GMI联系起来，为分子表型与疾病之间的关联提供了更深入的生物学见解。

本研究的局限性包括无法直接将结膜基因表达与基因型、微生物组及其交互作用之间的关联解释为因果关系。需要进一步采用实验方法来阐明这些关联和交互作用的潜在机制。

结  论

综上所述，本研究鉴定了结膜基因表达的eQTL以及基因型控制条件下结膜基因表达与微生物丰度之间的相关性，并证实GMI可解释一部分基因表达表型的变异。我们开发了MicroGenix，这是一个可以利用多组学数据识别GMI 的R语言工具包，并且可以根据基因型和微生物组数据优先考虑与复杂性状相关的SNP-微生物交互作用。

方  法

多组学数据

结膜多组学数据来自广州双生子研究（GTES）。从2006年起，居住在广州的双生子被招募参加GTES。所有参与者在入组时提供了书面知情同意。研究遵循赫尔辛基宣言的人体研究原则，并获得了中山大学中山眼科中心伦理委员会的批准。本研究纳入110对同卵双生子和10对异卵双生子的数据，包括117名男性和123名女性，平均年龄为18.9岁（18.9 ± 2.9，平均值±标准差）。本研究分析了480个结膜样本（来自240名个体的双眼）和240个血液样本。所有双生子都是健康受试者，没有系统性疾病或眼病的证据。至少有一名个体在过去6个月内有抗生素治疗或隐形眼镜佩戴史的双生子被排除在外。230名健康受试者、89名过敏性结膜炎患者和166名干眼病患者的宏基因组数据在我们的既往研究中收集。

样本收集和处理

按照我们之前建立的方法收集结膜印迹细胞学样本。简而言之，在采样前用1-2滴Alcaine眼药水（Alcon, Fort Worth, TX, USA）对眼睛进行局部麻醉。将一个无菌的半圆形膜滤器（REF: HAWP01300, 0.45 µm; Merck Millipore, Burlington, MA, USA）放置在下球结膜上10秒钟。然后将膜立即放入含有300 µl组织和细胞裂解液（Epicentre, Ambleside, UK）的无菌管中，并储存在-80°C。未用于采样并暴露于与用于采样的膜相同环境的膜被收集作为阴性空白对照。使用MasterPure Complete DNA和RNA纯化试剂盒（Epicentre, Ambleside, UK）按照制造商的说明从结膜样本中提取RNA。阴性空白对照与结膜样本一起使用相同的程序进行处理，包括储存、RNA提取、文库制备和测序。血液样本被收集并储存在-80°C。使用TIANamp Blood DNA试剂盒（DP348, TianGen Biotech, Beijing, China）按照制造商的说明提取DNA，并进行SNP分型。

测序

使用VAHTS Total RNA-seq文库制备试剂盒（Vazyme, Nanjing, China）按照制造商的说明为RNA-seq准备文库。这个试剂盒可以耗尽核糖体RNA和DNA，留下mRNA和其他非编码RNA。使用KAPA SYBR FAST qPCR试剂盒（Kapa Biosystems, Wilmington, MA, USA）通过qPCR定量文库。在NovaSeq 6000仪器（Illumina, San Diego, CA, USA）上进行2×150-bp的双端测序。

SNP阵列数据的插补

使用Illumina Infinium Asian Screening Array进行基因分型。我们使用PLINK在SNP和样本水平上对基因型数据进行质量控制。在SNP水平上，我们排除了次要等位基因频率（MAF）< 1%、基因型缺失率> 5%和偏离哈代-温伯格平衡（HWE）的SNP（p < 10^-5）。在样本水平上，我们移除了基因型缺失率> 5%和杂合度值超出平均值±3标准差的样本。我们进行了同系身份分析，并在存在一级或二级关系的情况下保留了基因型缺失率较低的样本。我们进行了主成分分析（PCA），并在前五个主成分（PC）中偏离平均值±6标准差的样本被排除。使用SHAPEIT4和IMPUTE5以及1000 Genomes Project数据集（第三阶段）的参考面板进行质量控制的基因型数据的阶段和插补。排除了MAF < 1%或INFO值< 0.8的SNP。插补后概率< 0.9的基因型被视为缺失，基因型缺失率> 5%的SNP被移除。

宏基因组测序数据的插补

使用内部管道将宏基因组数据作为浅全基因组测序数据进行基因分型。简而言之，质量控制后的读取与人类参考基因组（GRCh37）使用BWA进行比对。消除PCR重复。使用GLIMPSE进行基因型插补。用于插补的单倍型参考面板来自1000 Genomes Project数据集。排除了MAF < 1%和INFO分数< 0.8的SNP。基因型概率< 0.9的位点被视为缺失。

人类基因表达分析

RNA-seq数据的比对和定量是遵循GTEx项目的RNA-seq管道进行的。简而言之，RNA-seq数据与人类参考基因组GRCh37使用STAR进行比对。使用Picard套件标记重复读取。使用RNA-SeQC进行基因水平表达定量。定量基于GENCODE Release 41注释。

微生物组成分析

在质量过滤后，使用KneadData去除来自人类基因组和核糖体RNA的读取。使用Kraken 2将读取映射到自定义数据库。总共从国家生物技术信息中心（NCBI）下载了29,943个完整的微生物基因组，其中19,362个是细菌，368个是古菌，9346个是病毒，867个是真菌。使用置信度分数0.20过滤分类结果。只有在至少一个样本中读取数超过10个且在> 10%的样本中存在的物种才被保留。使用decontam R包的isContaminant函数的基于频率和基于流行度的方法进行污染物识别。在空白样本中检测到的物种被进一步移除。

eQTL定位

我们遵循GTEx项目的流程进行顺式eQTL定位。读取计数使用修整均值的M值（TMM）进行标准化。将在≥ 20%的样本中表达量≥ 0.1 TPM且在≥20%的样本中≥ 6个读取的基因纳入分析。每个基因在样本中进行逆正态转换。我们使用在peer R包中实现的PEER方法估计基因表达水平的潜在协变量。在顺式eQTL定位中，每个家庭只包括双生子的其中一个个体。对于每对双生子，选择基因分型率较高的双生子。个体的左眼和右眼的基因表达水平被平均。使用修改版的FastQTL进行顺式eQTL定位，参数与GTEx项目相同：定位窗口设置在转录起始位点（TSS）的±1 Mb，应用FDR阈值≤ 0.05。FDR由FastQTL内部使用置换方法估计。我们包括了前15个PEER因子、前五个基因型PC、年龄和性别作为协变量。GTEx项目的eQTL数据从GTEx门户网站下载（https://www.gtexportal.org）。

转录组和微生物组的双生子对内相关性

我们计算了同卵双生子内基因表达和微生物丰度的倍数变化之间的Spearman秩相关。通过计算双生子中一个个体的左眼和右眼的基因表达或微生物丰度的平均值，将一个双生子的平均值除以共生双生子的平均值，并进行log₂转换来获得倍数变化。对于全转录组关联分析，我们使用lme4 R包拟合线性混合模型来分析所有宿主基因表达与微生物丰度的关联。在每个模型中，我们将宿主基因对微生物物种进行回归，将家庭作为随机截距校正。这个测试可以被视为配对设计（双生子对中的两个个体）的配对t检验的泛化，线性混合模型的p值使用lmerTest R包计算。与之前的研究一致，使用Benjamini-Hochberg FDR < 0.1来确定显著的基因-微生物对。

通路富集分析

我们使用clusterProfiler R包基于京都基因与基因组百科全书（KEGG）数据库和基因本体（GO）数据库进行通路富集。

基于基因型和微生物组的基因表达建模

基因表达预测模型以表达为响应变量，以SNP剂量、微生物物种丰度和交互项为预测因子进行拟合。模型定义为：

其中Y是基因的表达，w_k是SNP k对基因的效应大小，S_k是SNP k的等位基因剂量；w_j是微生物j对基因的效应大小，M_j是微生物j的log₂转换丰度；w_i是SNP-微生物交互i对基因的效应大小，X_i是SNP-微生物对i的交互项，ε是假设与遗传和微生物组分独立的其他因素的贡献。

基因表达使用TMM进行标准化，并使用线性模型调整年龄、性别、前五个基因型PC和前15个PPER因子，残差用于拟合模型。纳入从基因起始上游1 Mb到基因结束下游1 Mb的SNP和所有微生物物种。SNP-微生物交互项，即SNP剂量和微生物丰度的乘积，用于表示GMI。在拟合包含交互的模型之前，先拟合不包含交互的模型，以识别用于构建交互项的重要SNP及微生物：

弹性网络使用glmnet R包进行拟合。对于在不包含交互的模型中具有非零系数的每个SNP和微生物，构建了交互项。使用stats和lme4 R包分别拟合了多元线性回归和线性混合模型，其中线性混合模型将家庭作为随机截距校正。对于未正则化的模型，为在不包含交互的模型中具有名义显著性（p < 0.05）的每个SNP和微生物构建了交互项。

基因预测准确性的交叉验证

嵌套交叉验证的执行方式如下：我们首先从100对双生子（n = 400）中随机选择样本作为训练数据，剩余的样本作为保留的测试数据。通过在家庭水平选择样本，我们确保了同一对双生子的样本在相同的集合中（无论是作为训练还是测试数据），以避免数据泄露。该过程针对每个折（fold）重复，测试集中的决定系数（R²）的平均值用于衡量准确性。

MicroGenix方法

MicroGenix方法已经实现为一个R包，包含三个模块：训练、预测和关联。训练模块接受SNP剂量、微生物丰度和分子表型（如基因表达）作为输入，并执行两阶段拟合以构建包含交互的模型。预测模块使用这些构建的模型从SNP剂量和微生物丰度预测分子表型。关联模块使用逻辑或线性回归测试预测的分子表型与复杂性状（如疾病）之间的关联，这取决于性状是定性的还是定量的。它还从与性状相关的分子表型模型中检索SNP-微生物交互项。

识别与过敏性结膜炎相关的GMI

利用在双生子数据集中拟合的模型，预测来自健康个体和过敏性结膜炎或干眼病患者的宏基因组数据集的基因表达。分析中包括了年龄和性别匹配的健康个体和患者的样本。使用逻辑回归分析表达量预测值与疾病的关联。分析纳入在包含交互的模型中R² > 0.01的2086个基因。采用了全基因显著性阈值p < 2.4×10^-5（0.05/2086），遵循全转录组关联研究的常见做法。

代码和数据可用性

测序数据已存储在国家基因组科学数据中心（National Genomics Data Center）（https://ngdc.cncb.ac.cn/），登录号为HRA005329。R包MicroGenix可在GitHub获取：https://github.com/MicrobeLab/MicroGenix。作图数据及脚本已存于GitHub：https://github.com/MicrobeLab/conjunctival-multi-omics。补充材料（图、表、脚本、图文摘要、幻灯片、视频、中文翻译版本和更新材料）可通过DOI或iMeta Science（http://www.imeta.science/）在线获取。

引文格式：

Qiaoxing Liang, Guo-Wang Lin, Xiaohu Ding, Bin Zou, Xiaomin Liu, Jing Li, Yuxin Zhang, Xiaofeng Wen, Lingyi Liang, Jin-Xin Bei, Mingguang He, Huijue Jia, Lai Wei. 2024. “Multiomics integration unravels genotype-microbiome interactions shaping the conjunctival transcriptome.” iMetaOmics 2: e37. https://doi.org/10.1002/imo2.37

作者简介

梁桥星（第一作者）

● 复旦大学生命科学学院博士后，博士毕业于中山大学中山眼科中心。

● 研究方向为微生物组；以第一作者在Cell、Allergy等杂志发表论文8篇。

林国旺（第一作者）

● 南方医科大学珠江医院副教授。

● 研究方向为复杂疾病遗传学；以第一作者在Lancet Oncology发表文章一篇。

魏来（通讯作者）

● 广州医科大学附属第二医院教授，广州医科大学眼科研究所所长，广东省过敏反应与免疫学重点实验室主任，入选国家高层次人才特殊支持计划（领军人才）、国家海外引才计划（青年项目）。

● 研究方向为人类疾病中host-microbe的相互作用；在包括Cell, Immunity, Cancer Cell等杂志发表论文130篇，总引用超过15000次。

贾慧珏（通讯作者）

● 复旦大学生命科学学院暨粤港澳大湾区精准医学研究院（广州）双聘研究员。

● 研究方向为人体共生微生物组；第一作者或通讯作者（含共同）成果包括Cell、Nature Biotechnology、Nature Medicine、Nature Genetics、Nature Reviews Microbiology等杂志论文，科睿唯安高引学者（1.8万次）；受邀出版首部人体共生微生物组与疾病研究英文专著。

何明光（通讯作者）

● 香港理工大学杰出创科学人计划讲席教授，国家杰出青年基金获得者、国家“万人计划”领军人才、广东特支计划杰出人才（“南粤百杰”）。

● 研究方向主要集中在眼科疾病的预防和治疗，特别是在近视控制、高度近视并发症的防治以及人工智能在眼科的应用方面。研究成果多次以论著方式发表在Lancet、JAMA等国际权威期刊。

iMetaOmics

更多资讯

●  iMeta姊妹刊iMetaOmics(定位IF>10)欢迎投稿！(2024.2.27)

●  iMeta姊妹刊iMetaOmics编委招募 (定位IF>10) (2024.3.2)

●  iMeta姊妹刊iMetaOmics电子版和印刷版ISSN申请获批(2024.4.1)

●  iMeta姊妹刊iMetaOmics投稿系统正式上线(2024.4.17)

●  iMeta姊妹刊iMetaOmics主编正式官宣(2024.4.22)

● 出版社iMetaOmics主页正式上线！(2024.4.28)

● iMetaOmics | 浙江大学宗鑫组揭示两猪种宿主-肠道菌群互作差异

● iMetaOmics | 罗鹏/袁硕峰/苗凯/程全发表STAGER: 生成式人工智能可靠性的标准化测试和评估推荐

● iMetaOmics | 徐州医科大杨欢组揭秘沙门氏菌-宿主-微生物群在免疫与代谢中的相互作

● iMetaOmics | 中科院动物所金坚石组综述16S rRNA基因扩增子测序技术的“前世今生”

● iMetaOmics | 浙大张天真组完成二倍体棉种泛基因组构建

● iMetaOmics | 张勇/李福平-先进糖蛋白组学在男性生殖研究中的潜在应用

● iMetaOmics | 暨南大学潘永勤/杨华组-炎症蛋白联合检测利于诊断甲状腺乳头状癌和结节性甲状腺肿

● iMetaOmics | 张开春组利用多组学方法揭示甜樱桃加倍后果色变化的候选基因

● iMetaOmics | 杜娟/林婷婷-慢性泪囊炎患者眼部菌群类型和纵向菌群变化

● iMetaOmics | 陈汉清/陈俊综述有关肝细胞癌治疗的新兴纳米医学策略

● iMetaOmics | 基因组所刘永鑫/卢洪评述微生物在提高杂种优势中的作用

● iMetaOmics | 上科大刘雪松组开发基于通路的肿瘤细胞鉴别工具TCfinder

● iMetaOmics | 中山大学刘鹏/邹宇田-整合人工智能实现HER2阳性乳腺癌精准管理

● iMetaOmics | 安徽农大李晓玉组-丛枝菌根真菌对玉米内生菌群的影响

● iMetaOmics | 徐涛/黄蓉/苏国海-急性冠脉综合征纵向多组学队列建设

● iMetaOmics | 通过整合宏组学促进人类与环境健康发展

● iMetaOmics | 苏州大学林俊组-揭示活性微生物及益生元/益生菌与关节炎联系

● iMetaOmics | 中国药科大学徐文波开发叶绿体基因组数据分析软件

● iMetaOmics | 清华刘晓组和复旦王久存组揭示特定细菌在皮肤老化中的作用

●iMetaOmics | 中南大学夏晓波团队揭示青光眼和SLE发病机制新关联

●iMetaOmics | 庐山植物园刘芬组揭示了自噬在植物-根微生物互作机制中的调控作用

●iMetaOmics | 杨瑞馥/袁静综述微生物组与“同一健康”的联系

●iMetaOmics | 同济/上海交大-开发支持群体分组分析的宏基因组测序综合分析软件

●iMetaOmics | 陈绍鸣-关于靶向NF-κB的潜伏逆转剂及其在HIV潜伏期的表观遗传和突变影响的评论

● iMetaOmics | 甘肃农大刘自刚组-强抗寒甘蓝型冬油菜的基因组组装和基因组特征解析

● iMetaOmics | 南京农大朱伟云组-外周血清素在结肠稳态中的作用

更多推荐

（▼ 点击跳转）

高引文章 ▸▸▸▸

iMeta | 引用14000+，海普洛斯陈实富发布新版fastp，更快更好地处理FASTQ数据

高引文章 ▸▸▸▸

iMeta | 德国国家肿瘤中心顾祖光发表复杂热图(ComplexHeatmap)可视化方法

高引文章▸▸▸▸

iMeta | 高颜值绘图网站imageGP+视频教程合集                                        

1卷1期

1卷2期

1卷3期

1卷4期

2卷1期

2卷2期

2卷3期

2卷4期

3卷1期

2卷2期封底

2卷4期封底

3卷2期

3卷3期

3卷3期封底

3卷4期

3卷4期封底

1卷1期

期刊简介

“iMeta” 是由威立、肠菌分会和本领域数百千华人科学家合作出版的开放获取期刊，主编由中科院微生物所刘双江研究员和荷兰格罗宁根大学傅静远教授担任。目的是发表所有领域高影响力的研究、方法和综述，重点关注微生物组、生物信息、大数据和多组学等。目标是发表前10%(IF > 20)的高影响力论文。期刊特色包括视频投稿、可重复分析、图片打磨、青年编委、前3年免出版费、50万用户的社交媒体宣传等。2022年2月正式创刊发行！发行后相继被Google Scholar、ESCI、PubMed、DOAJ、Scopus等数据库收录！2024年6月获得首个影响因子23.7，位列全球SCI期刊前千分之五(107/21848)，微生物学科2/161，仅低于Nature Reviews，同学科研究类期刊全球第一，中国大陆11/514！

“iMetaOmics” 是“iMeta” 子刊，主编由中国科学院北京生命科学研究院赵方庆研究员和香港中文大学于君教授担任，是定位IF>10的高水平综合期刊，欢迎投稿！

iMeta主页：

http://www.imeta.science

姊妹刊iMetaOmics主页：

http://www.imeta.science/imetaomics/

出版社iMeta主页：

https://onlinelibrary.wiley.com/journal/2770596x

出版社iMetaOmics主页：

https://onlinelibrary.wiley.com/journal/29969514

iMeta投稿：

https://wiley.atyponrex.com/journal/IMT2

iMetaOmics投稿：

https://wiley.atyponrex.com/journal/IMO2

邮箱：

office@imeta.science