A large number of genome-wide association studies (GWAS) in livestock, especially in dairy cows, provide favorable conditions to integrate multiple independent studies. Methods such as meta-analysis provide the identification of effective QTLs with higher precision and power. A meta-analysis for milk production traits between different countries was conducted using the GWAS summary statistics (i.e., P-value, sample size, allele effects, and etc.) in Holstein cows. In the present study, METAL software was used for the weighted Z-score model. Gene network analysis was used as a complementary method to improve our knowledge of the genome structure of milk production traits and was implemented through the STRING plug-in in Cytoscape software. The Cytoscape ClueGO plug-in was also used for GO enrichment in order to identify biological process, cellular component, and molecular function associated with genomic regions. The aim of this study was to improve the power of QTLs detection and identify the biological mechanisms associated with milk production traits. Data were obtained from 26 published studies from 2010 to 2019. A total of 2,072 SNPs were identified for milk production traits, of which 1,583 SNPs were significant ($P<0.05$). Meta-analysis identified 9 QTLs for milk yield, 36 QTLs for fat percentage, and 10 QTLs for protein percentage. Some QTLs were confirmed on BTA14, e.g., BTA14:1801116 close to the DGAT1 gene (milk yield, $P=2.6\times {10}^{-131}$; fat percentage, $P=4.8\times {10}^{-347}$; protein percentage, $P=7.6\times {10}^{-24}$) and BTA14:1651311 close to the PPP1R16A gene (milk yield, $P=2.3\times {10}^{-162}$; fat percentage, $P=3.5\times {10}^{-153}$). We identified pleiotropic effects of lead SNPs for milk production traits, e.g., one SNP (rs109421300) at BTA14 had pleiotropic effects on milk yield, fat percentage, and protein percentage traits. The most important SNPs for studied traits across countries implicated to network scoring and visualization were including: rs109421300 (DGAT1 gene) for milk yield, fat percentage, and protein percentage; rs109146371 (PPP1R16A gene) for milk yield and fat percentage; rs109968515 (CYHR1 gene) for milk yield and fat percentage; rs134432442 (CPSF1 gene) for fat percentage; rs111018678 (TRAPPC9 gene) for protein percentage. Significant pathways involved in milk production traits through GO term enrichment analysis for biological process, cellular component, and molecular function included: regulation of cation channel activity ($P=1.6\times {10}^{-2}$), ion channel complex ($P=1.4\times {10}^{-2}$), and phosphoric diester hydrolase activity ($P=1.1\times {10}^{-3}$) for milk yield; negative regulation of organ growth ($P=8.2\times {10}^{-3}$), transmembrane transporter complex ($P=1.6\times {10}^{-3}$), and potassium ion transmembrane transporter activity ($P=8.8\times {10}^{-3}$) for fat percentage; mRNA polyadenylation ($P=1.2\times {10}^{-2}$), mRNA cleavage factor complex ($P=9.8\times {10}^{-4}$), and phosphoric diester hydrolase activity ($P=5.3\times {10}^{-3}$) for protein percentage, respectively. Thus, the combination of GWAS summary statistics through a powerful methodology such as meta-analysis will assist us to accurately identify QTLs, potential candidate genes, and biological mechanisms. This kind of studies will help us to have better understanding of QTL regions and genome structure for milk production traits and improve genomic evaluations in Holstein cows. To the best of our knowledge, this is the first meta-analysis of GWAS and GO enrichment across countries for milk production traits in Holstein cows.

家畜尤其是奶牛的大量全基因组关联研究（GWAS）为整合多项独立研究提供了有利条件。诸如荟萃分析之类的方法提供了具有更高精确度和功效的有效 QTL 的鉴定。使用荷斯坦奶牛的 GWAS 汇总统计数据（即 P 值、样本量、等位基因效应等）对不同国家之间的产奶性状进行了荟萃分析。在本研究中，METAL 软件用于加权 Z 分数模型。基因网络分析被用作补充方法，以提高我们对产奶性状基因组结构的了解，并通过 Cytoscape 软件中的 STRING 插件实施。Cytoscape ClueGO 插件也用于 GO 富集以识别生物过程，细胞成分，以及与基因组区域相关的分子功能。本研究的目的是提高 QTL 检测的能力，并确定与产奶性状相关的生物学机制。数据来自 2010 年至 2019 年发表的 26 项研究。共鉴定出 2,072 个产奶性状 SNP，其中 1,583 个 SNP 显着（$磷<0.05$）。Meta 分析确定了 9 个产奶量 QTL、36 个脂肪百分比 QTL 和 10 个蛋白质百分比 QTL。在 BTA14 上确认了一些 QTL，例如 BTA14:1801116 靠近DGAT1基因（牛奶产量、$磷=2.6\times {10}^{-131}$; 脂肪百分比，$磷=4.8\times {10}^{-347}$; 蛋白质百分比，$磷=7.6\times {10}^{-24}$) 和 BTA14:1651311 接近PPP1R16A基因（产奶量，$磷=2.3\times {10}^{-162}$; 脂肪百分比，$磷=3.5\times {10}^{-153}$）。我们确定了先导 SNP 对产奶性状的多效性影响，例如，BTA14 的一个 SNP (rs109421300) 对产奶量、脂肪百分比和蛋白质百分比性状具有多效性影响。涉及网络评分和可视化的国家间研究性状的最重要 SNP 包括： rs109421300（DGAT1基因），用于表示产奶量、脂肪百分比和蛋白质百分比；rs109146371（PPP1R16A基因）用于产奶量和脂肪百分比；rs109968515（CYHR1基因）用于产奶量和脂肪百分比；rs134432442（CPSF1基因）用于脂肪百分比；rs111018678 ( TRAPPC9基因）表示蛋白质百分比。通过对生物过程、细胞成分和分子功能的 GO 术语富集分析，涉及产奶性状的重要途径包括：阳离子通道活性的调节。$磷=1.6\times {10}^{-2}$), 离子通道复合物 ($磷=1.4\times {10}^{-2}$)，和磷酸二酯水解酶活性 ($磷=1.1\times {10}^{-3}$) 产奶量；器官生长的负调控（$磷=8.2\times {10}^{-3}$), 跨膜转运蛋白复合体 ($磷=1.6\times {10}^{-3}$) 和钾离子跨膜转运蛋白活性 ($磷=8.8\times {10}^{-3}$) 脂肪百分比；mRNA 聚腺苷酸化 ($磷=1.2\times {10}^{-2}$), mRNA 切割因子复合物 ($磷=9.8\times {10}^{-4}$)，和磷酸二酯水解酶活性 ($磷=5.3\times {10}^{-3}$) 分别表示蛋白质百分比。因此，通过荟萃分析等强大方法结合 GWAS 汇总统计数据将有助于我们准确识别 QTL、潜在候选基因和生物学机制。此类研究将有助于我们更好地了解产奶性状的 QTL 区域和基因组结构，并改进荷斯坦奶牛的基因组评估。据我们所知，这是针对荷斯坦奶牛产奶性状的跨国家 GWAS 和 GO 富集的首次荟萃分析。