Molecular Mechanisms Associated with the Development of the Metritis Complex in Dairy Cattle

The metritis complex (MC), a group of post-partum uterine diseases, is associated with increased treatment costs and reduced milk yield and fertility. The goal of this study was to identify genetic variants, genes, or genomic regions that modulate MC disease. A genome-wide association study was performed using a single-locus mixed linear model of 1967 genotypes (624,460 SNPs) and metritis complex records. Then, in-silico functional analyses were performed to detect biological mechanisms and pathways associated with the development of MC. The ATP8A2, COX16, AMN, and TRAF3 genes, located on chromosomes 12, 10, and 21, were associated with MC at p ≤ 0.0001. These genes are involved in the regulation of cholesterol metabolism in the stromal tissue of the uterus, which can be directly associated with the mode of transmission for pathogens causing the metritis complex. The modulation of cholesterol abundance alters the efficiency of virulence factors and may affect the susceptibility of the host to infection. The SIPA1L1, DEPDC5, and RNF122 genes were also significantly associated with MC at p ≤ 0.0001 and are involved in the PI3k-Akt pathway, responsible for activating the autophagic processes. Thus, the dysregulation of these genes allows for unhindered bacterial invasion, replication, and survival within the endometrium.


Introduction
During the peripartum period, dairy cows are particularly vulnerable to uterine infection caused by pathogenic bacteria including Escherichia coli, Fusobacterium necrophorum, and Trueperella pyogenes, resulting in a significant increase in the prevalence of uterine diseases [1,2].The metritis complex refers to a persistent group of post-partum uterine diseases that affect 20-40% of dairy herds, with a 25% average incidence rate per herd [3], and include all forms of endometritis, metritis, and pyometra diagnosed after parturition.These diseases result in an annual cost of approximately USD 900 million to the US dairy industry due to their high prevalence, treatment cost, and negative impact on milk production, uterine functionality, and reproductive efficiency [4][5][6][7].These infections often follow placental retention, abortion, and dystocia, and are accompanied by clinical signs of pyrexia, loss of appetite, depression, and reduced feed intake [1,7].In an effort to improve diagnostic accuracy, clinical endometritis (an infection of the endometrial tissue producing purulent or mucopurulent discharge), metritis (fetid reddish-brown discharge), and pyometra (pus accumulation in the uterus concomitant with a corpus luteum) have been consolidated into a single category, defined as the metritis complex (MC), due to the variations in clinical diagnosis and often ambiguous commercial dairy operations records [8].
Recent years have seen an increase in studies focused on genomic selection for calving ease, fertility, and health traits, leading to an increase in the genetic trends of functional traits, and subsequently bringing to light the negative impacts of MC in dairy cattle [9].Health traits, such as the development of MC, have been found to be highly complex due to their complex genetic architecture [10,11].In addition, these traits exhibit lower genomic heritability due to high environmental variance, making it difficult to detect genetic variants and genes prompting variations among animals in these traits.Nonetheless, with the advent of high throughput sequencing and genotyping, genome-wide association studies (GWASs) have been a powerful tool to identify genes for common diseases and quantitative traits [12].Recent studies using single-nucleotide polymorphism (SNP) data estimated the heritability of metritis from 0.004 to 0.07 [13][14][15], indicating that genetic variance could be captured [14,16].This does not mean that MC is not controlled by genes, but the genomic variances caused by genes are relatively small compared to the environmental variance.If genetic variation is detected for uterine disease traits, this allows for genetic progress at a slow rate and the detection of genes or molecular mechanisms with strong effects on these traits [15].Thus, a few GWASs have successfully identified genetic variants (e.g., single-nucleotide polymorphisms (SNPs)) involved in the genetic basis of health traits in cattle using small to medium-density SNP panels [13,15,17].Potential SNPs and genes on chromosomes BTA 1, 2, and 21 were associated with endometritis in German Holstein cows [17].Recently, another GWAS on German Holstein Friesian cattle revealed 24 potential genes for all uterine disease traits [13].Another GWAS discovered 51 genes located on nine different chromosomes that are associated with metritis in Canadian Holstein cows [15].Nonetheless, causative variants for MC remain relatively unknown, and missing heritability could be claimed [18].In addition, there is a need to detect genes and pathways from commercial dairy cattle that could be validated on the same herds using integrated omics approaches such as metabolomics and proteomics.It is also important to validate previously reported genes and genomic regions for uterine disease traits using higher-density SNP panels-used in the current study-that can be used in genomic selection.Genomic selection offers a promising opportunity to improve health traits including MC.Therefore, the objectives of this study were to (1) estimate genomic heritability for MC using HD SNP panels, (2) identify genetic variants and genes associated with the development of MC via a GWAS using high-density chip panels (HD SNPs), and (3) detect the biological mechanisms or pathways associated with MC development.

Animals and Phenotypic Data
This study involved Holsteins and Jersey breeds at three dairy operations: Herd 1 (located in the coastal region of central California) and Herds 2 and 3 (located in central California).We extracted the historical records of primi-and multiparous dairy cattle born between 2015 and 2021 from a complied producer-recorded dataset of over 17,000 animals using DairyComp 305 TM Herd management software Release 27 (Valley Agricultural Software, Tulare, CA)).For the purpose of this study, diagnosed cases of metritis, endometritis, or pyometra within the transition period were collectively grouped as the phenotype referred to as the metritis complex (MC) and designated as cases [7,19,20].Thus, MC includes the diseases defined as clinical metritis, clinical endometritis, and pyometra, as described by Sandals et al. (1979) [8].All diseases were diagnosed by a trained herdsman or veterinarian.
In the case of one or multiple documented MC diagnoses, these cows were assigned as cases with MC.All cows without a diagnosis of MC or any other disease were assigned as controls (i.e., healthy).

DNA Extraction and Genotypic Data
A total of 2120 cows were sampled using different approaches.Blood collected from the coccygeal vein (n = 679) and ear notch samples (n = 1441) were sent for DNA extraction using the standard commercial kit and sequenced for low-pass skim sequencing by NEOGEN ® Genomics.The genomic coverage was from 0.5 to 3.0× and the concordance averaged 99.3% in genotype concordance to the illumina SNP chips [21].The assessment of imputation accuracies from low-pass sequencing to commercial SNP panels was performed using the Gencove pipeline [21] adopted by NEOGEN ® Genomics.For this study, the bovine high-density (HD; 782,138 SNPs) panel genotypes were extracted and utilized for GWAS analyses.Quality control was applied to genotypes using the PLINK v1.9 [22] software.Genetic variants with a minor allele frequency (MAF) of less than 0.01 and those deviating from the Hardy-Weinberg equilibrium were removed.Additionally, animals with duplicate genotypes and breeds other than Holstein and Jersey were excluded from the analysis.The final dataset included 624,460 SNPs and a total of 1967 cows consisting of 442 cases with MC and 1525 healthy controls.

Population Structure
To assess population stratification, a multidimensional scaling (MDS) analysis based on genome-wide identity-by-state pairwise distances was performed using PLINK v1.9 [22].Linkage disequilibrium pruning was also conducted to reduce the likelihood of obtaining principal components based on a few genomic regions, and the analysis was carried out with 66,245 genetic variants.The MDS plot was visualized using "scatterplot3d" and "ggplot2" in R [23][24][25].

Genome-Wide Association Study (GWAS) and Genomic Heritability
The GWAS analysis was conducted on animals (n = 1967) with genotypes (624,460 SNPs) and phenotypes that passed the quality control.Variance components and genomic heritability were estimated using the genomic best linear unbiased prediction (gBLUP) statistical method implemented in the Genome-wide Complex Trait Analysis (GCTA) software version 1.94.1 [26] through a restricted maximum likelihood algorithm (REML).To perform the GWAS, a single-SNP mixed linear animal model was constructed using GCTA software.The allele substitution effect and the association level of significance were estimated.The model included the fixed effects of breed, herd, and SNP and utilized the genomic relationships of all animals.The model is represented as follows: where Y ijk is the diagnosis of MC (i.e., phenotypes) in the kth animal of the ith breed, µ is the overall mean for the trait, Herd j is the fixed effect of the jth herd, Breed i is the fixed effect of the jth breed, β 1 is the allele substitution effect of the candidate SNP being tested, SNP represents the SNP genotype variable (coded as 0, 1, or 2), a k is the random additive genetic (polygenic) effect of the kth animal, and e jklm is the residual random effect associated with the kth animal record.Assumptions for this model were a k : a~N (0, G σ 2 a), where G is a genomic relationship matrix and σ 2 a is the additive genetic variance, and e ijk : e~(0, I σ 2 e), where I is the identity matrix and σ 2 e is the error variance.The expectations were that E(a k ) = 0 and E(e ijk ) = 0, and the variances were Var(a k ) = σ 2 a and Var(e ijk ) = σ 2 e. Gσ 2 a is the covariance matrix of the vector of genomic additive genetic effects and the genomic relationship matrix (G).
To declare a significant association between an SNP and MC, we considered a suggestive significance candidate threshold (p-value = 1 × 10 −4 ) [13,27,28].We also identified the top 100 significant (p ≤ 0.00025) SNPs for further functional analyses as the genomic inflation factor indicated a slight underestimation of p-values.We used a less stringent p-value threshold to declare the significance associated with MC for two main reasons: (1) the known polygenic nature of health traits and (2) to increase the power of quantitative trait loci (QTL) detection (e.g., genetic variants/genes/genomic region) with small but true effects on MC development.In addition, relaxing the threshold for declaring significant association is consistent with previous publications on GWASs for fertility and health traits in dairy cattle using SNP panels [29,30].To visualize the GWAS results, a Manhattan plot and Quantile-Quantile (Q-Q) plot were generated using the R package "qqman" [31].The Q-Q plot and the genomic inflation factor were used to evaluate the bias in the estimated p-values from the association analyses.The genomic heritability of MC was estimated using the additive genomic variance divided by the total phenotypic variance calculated via REML methods implemented in GCTA [26].The proportion of genetic variance explained by the top 100 SNPs for MC was estimated using GVCBLUP software version 3.9 [32].

In-Silico Functional Analyses
The SNPs significantly (p ≤ 0.00025) associated with MC were mapped to the bovine genome ARS-UCD1.2assembly using the BioMart tool in the Ensembl database [33] to identify the nearby gene(s) located within 10 kilobases upstream or downstream.Enrichment analyses of Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways and gene ontology (GO) terms were performed using the Database for Annotation, Visualization, and Integrated Discovery (DAVID) bioinformatics tool [34].To investigate if any of the identified candidate genes in the current study had previously been reported for QTLs related to uterine infection, the Cattle QTL database Release 50 [35] was accessed on 18 April 2023.

Results
The descriptive statistics of the observed phenotypes and incidence rates separated by herd are presented in Table 1.A total of 1967 out of 2124 multiparous cows were used in this study after removing duplicate animals, low-quality genotypes, and other breeds.The records for Jersey (13.5%) were fewer than those of the Holstein (86.5%) breed.Nonetheless, Jersey (22.8%) and Holstein displayed similar diseased phenotype diagnoses, with 22.8% and 21.8% for Holstein and Jersey, respectively.Of the three herds, the incidence rate of cattle diagnosed with a uterine disease within the metritis complex (endometritis, metritis, and pyometra) was 33.7%, 16.1%, and 24.5% with respect to Herds 1, 2, and 3. 1 Case = the diagnosis of at least one disease from the metritis complex; Control = no prior diagnosis from the metritis complex.

Population Structure
The MDS analysis, which is similar conceptually to principal component analysis, revealed no stratification due to herd and breed, as illustrated in Figure 1a,b.The results from plotting PC1 against PC3 showed no difference in genetic structure stratification, illustrated in Figure 1c,d.Nonetheless, to account for the slight stratification in family structure, shown in Figure 1, we included the cluster groups obtained by MDS analysis as fixed effects in the MLM model, but this caused a severe underestimation of observed p-values in comparison with the expected p-values.In addition, the genomic inflation factor from the model with MDS decreased, indicating an underestimation of SNP effects and overparameterization of the model.Thus, the cluster groups were eliminated from the model in further analyses.The value of the genomic inflation factors from the GWAS evaluated models ranged between 0.97 and 0.999, indicating the absence of bias in the estimated SNP effects or p-values due to family structure.

FOR PEER REVIEW
5 of 14 the model.Thus, the cluster groups were eliminated from the model in further analyses.The value of the genomic inflation factors from the GWAS evaluated models ranged between 0.97 and 0.999, indicating the absence of bias in the estimated SNP effects or p-values due to family structure.

Genome-Wide Association Study (GWAS) and Genomic Heritability
The results from the GWAS, presented as a Manhattan plot displaying the -log10 pvalue of each polymorphism with respect to their autosomal positions, for the metritis complex are illustrated in Figure 2. In this study, the genomic heritability was 0.04 (±0.02) and the calculated genomic inflation factor (λ) was 0.98.The Q-Q plot confirmed the findings from the genomic inflation factor, with a slight underestimation of p-values (Supplementary Figure S1).

Genome-Wide Association Study (GWAS) and Genomic Heritability
The results from the GWAS, presented as a Manhattan plot displaying the −log10 p-value of each polymorphism with respect to their autosomal positions, for the metritis complex are illustrated in Figure 2. In this study, the genomic heritability was 0.04 (±0.02) and the calculated genomic inflation factor (λ) was 0.98.The Q-Q plot confirmed the findings from the genomic inflation factor, with a slight underestimation of p-values (Supplementary Figure S1).

In-Silico Functional Analyses
The enrichment analyses for gene ontology and the KEGG pathway for the 20 identified genes within or flanking the regions associated with the metritis complex were not significant at a 5% false discovery rate.Nonetheless, GO terms Regulation of Dendrite Morphogenesis, Regulation of Axonogenesis, Ephrin Receptor Signaling Pathway, Actin Cytoskeleton Reorganization, and GTPase Activator Activity were enriched at p < 0.05, (Supplementary Tables S2-S5).The summary of all GO terms and KEGG pathways linked to the gene list associated with the metritis complex are given in Supplementary Tables S2-S5.

Discussion
Diseases such as metritis, endometritis, and pyometra have been the topic of numerous studies within the last decade due to their negative effects on animal welfare, animal efficiency, and the economic welfare of the dairy industry [36][37][38].Thus, one of the main goals of the current study was to understand the genetic architecture of uterine disease diagnoses in the post-partum period.In the present study, the estimated genomic heritability was 0.04 for the metritis complex.The proportion of genetic variance explained by the top 100 SNPs for MC was 0.1792 ± 0.04.While some studies have reported pedigree-based heritability of metritis, there are few published estimates for genomic or SNP-based heritability for MC.A previous study using imputed and actual 44,747 genotypes and a linear model estimated genomic heritability as 0.008, 0.01, and 0.004 for metritis, endometritis (mucopurulent), and pyometra, respectively, and gave an estimate of 0.02 when assessing the genomic heritability of the metritis complex [13].Their genomic heritability was smaller than what was reported based on pedigree-based estimates due to the inability of the imputed 50K SNP chips to fully capture the genetic variation in the whole genome [13].In the current study, we utilized high-density SNP chips (624,460 SNPs), which resulted in similar heritability to that when using pedigree-based analysis with first-parity records but slightly smaller heritability than that with genomic-based analysis (h 2 = 0.06-0.07)reported from producer-recorded metritis event data in US dairy cattle [14].Thus, these results support that genomic improvement via selection for MC is feasible in the United States using producer-recorded data and combining uterine diseases as part of the MC [14].
In this genome-wide association study, the top 100 significant SNPs were within or flanking 20 genes.These genomic regions were previously reported to be associated with other QTLs, including placenta retention (n = 1), calving ease (n = 4), stillbirth (n = 4), and dystocia (n = 5), and each has also been reported for their roles in the development of uterine infection or because of uterine infection in dairy cattle.With this, as many as 54.8% of cows that have placental retention are soon diagnosed with at least one of the diseases that are part of the metritis complex, and dystocia was seen preceding 40%, 48%, and 10%, of metritis, endometritis, and pyometra cases, respectively [39,40].Here, calving ease, or the scale of calving difficulty experienced during parturition, is interchangeable with the true definition of dystocia and closely related to stillbirth [41].The incidence of stillbirth has been previously associated with the pathogenic invasion of the maternal and fetal tissues by Arcanobacterium pyogenes and E. coli, two prominent bacterial species closely tied to the development of metritis and endometritis [1,42].Further, three QTLs for non-return rates, the proportion of cattle that is not re-bred within a certain time frame proceeding insemination, corresponded to the three separate significant regions on autosomes 12, 15, and 17, suggesting their relevance in post-partum uterine disease traits [43].Metritis has been found to impede ovarian functionality, prolong ovulatory periods, and disrupt breeding timeframes, thereby contributing to an increased non-return rate [44].Interestingly, cows diagnosed with placental retention and/or a disease within the metritis complex have a reduced likelihood of becoming pregnant post-artificial insemination and a heightened chance of pregnancy loss [45].With this, a transparent linkage is seen between multiple uterine infections, the precursory complications in the antepartum to post-partum period, and their resulting impact on future fertility traits, forming a multiplex of reproductive setbacks and economic deficiencies.This transpicuous linkage shown between the recognized QTLs and uterine disease in dairy cattle validates the power of our GWAS and its ability to identify potential genes and genetic variants associated with the metritis complex.The current findings confirmed how complications throughout pregnancy increase vulnerabilities within the vaginal canal, cervix, placenta, and endometrium and establish an environment destined for infection [46,47].
The most significant SNP, rs133231370, corresponded to a nearby hub gene, ATP8A2, known for its involvement in lipid metabolism.Previously, ATP8A2 has been discussed for its role in milk fat synthesis in Holstein cattle [48].This gene produces an ATPase protein responsible for transferring, or flipping, phosphatidylserine (PS) and phosphatidylethanolamine from the ectoplasmic to the cytoplasmic layers of the cell membrane lipid bilayer.The nature of this function results in the asymmetric partitioning of lipids across the membrane, which is vital for vesicle trafficking, cell signaling, apoptosis, cell survival, and cholesterol/bile homeostasis [49].The ATP8A2 gene is also expressed within the uterus, where it plays a role in lipid metabolism and cholesterol abundance.Interestingly, higher cholesterol abundance within the uterine endometrial tissue increases pathogenic invasion through pyolysins, cholesterol-dependent cytolysins from T. pyogenes, a prominent Gram-positive bacteria in purulent infections, such as the metritis complex [50].
Importantly in the current study, five other genes expressed within the uterus corresponding to significant SNPs (SLC10A1, COX16, AMN, TRAF3, and POR) were found to regulate cholesterol levels.The TRAF3 and COX16 genes work to regulate intracellular cholesterol indirectly through their involvement in neighboring pathways.The TRAF3 gene, an immunity-related gene in cattle, was previously found to be a negative regulator of the nuclear factor-κB (NF-κB) pathway, a pathway known to increase cholesterol accumulation within the cell to promote atherosclerosis and macrophage foam cell formation in recent culture studies [51,52].Similarly, COX16 has been shown to inhibit the signaling activity of p53, which cooperates with the Hippo pathway to regulate the downstream activity of the sterol regulatory element-binding proteins (SREBPs; [53,54]).These proteins activate a multitude of genes responsible for the uptake and metabolism of phospholipids, triglycerides, fatty acids, and, cholesterol [55].Conversely, the activity of SLC10A1, AMN, and POR works directly to alter cholesterol levels within the endometrium.The POR gene encodes the cytochrome p450 protein, previously seen in differing levels between healthy and metritis cows, and is known to produce superoxides during the oxidation process, contributing to the reactive oxygen species (ROS) pool [56,57].This ROS activity increases cholesterol and glucose influx and the synthesis of cholesterol from glucose [58].Likewise, SLC10A1, once translated, produces a protein that co-transports sodium and bile acids, the catabolic product of cholesterol, in the dairy cow post-partum endometrium, showing the importance of SLC10A1 in intracellular cholesterol homeostasis [59,60].The AMN gene, found to be expressed within the cattle uterine epithelial tissue at all stages of estrous, encodes the amnionless protein, which functions to anchor another protein, cubilin, to the membrane [61].Cubilin is an endocytic lipoprotein receptor known to mediate high-density lipoprotein (HDL) cholesterol endocytosis [62].Our study observed a functional relationship between the genes corresponding to SNPs significantly associated with the metritis complex in dairy cows and the importance of cholesterol regulation in uterine tissues.This study confirms that there is crosstalk between cholesterol abundance, pathogenic invasion, and susceptibility to uterine infection; however, details should be further elucidated [50,63].
Moreover, in the current study, 9 of the 20 genes found flanking the top significant SNPs were discovered to be linked to the phosphatidylinositol-3-kinase/protein kinase B (PI 3 K/Akt) pathway.The PI 3 K/Akt pathway has been recently discussed in terms of its essential role in signal transduction, cell proliferation, apoptosis, and autophagy and its correspondence to pathogenesis within bovine endometrial epithelial cells [64].The activated pathway phosphorylates Akt downstream, and active Akt has been found to deter the fusion of the autophagosome to the lysosome, allowing for the survival of pathogenic materials and promoting infection [65,66].The SIPA1L1, DEPDC5, and RNF122 genes encode proteins that repress the RAS, mTOR, and Rig-1 signaling pathways, respectively [67][68][69].These secondary pathways work to inhibit, activate, and compete with the activity of PI 3 K/Akt and, consequently, regulate autophagosome elongation, maturation, and termination and the autophagy process [70][71][72].The PI 3 K/Akt pathway was considerably discussed for its role in bovine physiology including oocyte competence, immune responsiveness, milk fat synthesis, bacterial resistance, vascular homeostasis, and angiogenesis [64,[73][74][75][76].Our study reports the link between the development of uterine infections in dairy cattle and the PI 3 K/Akt pathway, considering its major role in pathogenesis and the immune response.Little is known about the biological function of the remainder of the identified genes, CNTN5, SRSF5, SLC24A2, SLC8A3, and KIAA1671, regarding B. taurus.
This study revealed unique genetic variants and corresponding candidate genes associated with the metritis complex.This study acknowledges that there is no advantage in accuracy when using an HD SNP panel compared to that of a medium-density panel for the GBLUP method in Holstein and Jersey cows.With this, the power of associated polymorphism detection may be affected by including data from Jersey breeds compared to those of the Holstein breed, given that detecting QTLs for a multi-breed population requires a larger sample population size [77].Nonetheless, utilizing a multi-breed population for reference permits an increase in the accuracy of genomic estimated breeding values (GEBVs) for smaller breeds and the potential use of the detected genetic variant for a crossbreed genomic evaluation.
Previously identified fertility-related genes, such as POR and SCRN1, were also revealed in this study as candidate genes corresponding to uterine disease in the transition period [57,78].A number of the putative genes for the metritis complex have already been associated with immunity and reproductive efficiency in past publications [79,80].The underlying common biological processes amongst the identified candidate genes were further confirmed using the correspondence of past literature in both human and bovine disease studies.As the metritis complex is a complex polygenic trait and develops as a collective infection involving the detriment of multiple species, studies looking to further analyze how these candidate genes may influence the ability of dairy cows to resist infection are required.In addition, an unrelated population is required to validate the detected significantly associated SNPs and genes from this study.Furthermore, these SNPs must be tested for genomic prediction prior to application in breeding programs.In this case, the Jersey breed should be analyzed in a larger sample size to corroborate previous results.The novel SNPs and their corresponding genes identified in this study should be considered candidates associated with the development of uterine disease in dairy cattle.Even so, these results and genetic variants should be further studied in regard to their influence on bacterial resistance and susceptibility to diseases defined under the metritis complex in the transition period.However, future studies using larger sample populations and more dense genotypes (e.g., sequence genotypes with millions of genetic variants) may increase the resolution of detecting the causal mutation and candidate genes for the metritis complex.Further validation studies are required to confirm the detected genetic variants and genes before use in genomic prediction for the metritis complex.

Conclusions
This study provided a novel perspective on the genetic architecture of uterine infection, identifying genetic variations and genes associated with the metritis complex in Jersey and Holstein cattle.We detected 40 SNPs and 20 candidate genes flanking the significant genomic regions.Several genes expressed in the reproductive tract were found to be linked to two processes in the endometrial tissue: the PI 3 K/Akt pathway and cholesterol homeostasis.These findings provide new insights into the proportion of genetic variability elucidated by SNPs for uterine disease and shed light on the physiological vulnerabilities that can lead to the pathogenesis of metritis, endometritis, and pyometra after calving.The results of this study should be taken into consideration when selecting index Holstein and Jersey dairy cattle.

Supplementary Materials:
The following supporting information can be downloaded at https://www.mdpi.com/article/10.3390/genes15040439/s1: Figure S1: Quantile-Quantile (QQ) plot of the expected vs. the observed null distribution of the p-values for the association between the 624,460 SNPs and 1967 dairy cattle across three herds diagnosed with at least one disease in the metritis complex.Table S1: Identified single nucleotide polymorphisms (SNPs) with flanking candidate genes and associated quantitative trait loci (QTL).Table S2: Enriched gene ontology (GO) biological process (BP) terms for the metritis complex in Jersey and Holstein dairy cattle.Table S3.Enriched gene ontology (GO) cellular component (CC) terms for the metritis complex in Jersey and Holstein dairy cattle.Table S4.Enriched gene ontology (GO) molecular function (MF) terms for the metritis complex in Jersey and Holstein dairy cattle.Table S5.Enriched Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways for the metritis complex in Jersey and Holstein dairy cattle.

Figure 1 .
Figure 1.The multidimensional scaling (MDS) analysis.Representative of the 624,460 single-nucleotide polymorphisms (SNPs) and 1967 dairy cows separated according to (a) breed and (b) herd.The 3D format shows the population structure of the first three principal components (PCs), representing the distribution of genetic variation across breeds and herds.The 2D format shows a visualization for (c) breed and (d) herd between PC1 and PC3.

Figure 1 .
Figure 1.The multidimensional scaling (MDS) analysis.Representative of the 624,460 singlenucleotide polymorphisms (SNPs) and 1967 dairy cows separated according to (a) breed and (b) herd.The 3D format shows the population structure of the first three principal components (PCs), representing the distribution of genetic variation across breeds and herds.The 2D format shows a visualization for (c) breed and (d) herd between PC1 and PC3.

Figure 2 .
Figure 2. Genome-wide association study (GWAS) results from 624,460 SNPs and 1967 dairy cattle across three herds diagnosed with at least one disease in the metritis complex.Manhattan plot for −log10 p-values of SNP effects for the metritis complex in Holstein and Jersey dairy cows.The red horizontal line corresponds to a p-value threshold of −log10 ≥ 4.00.The blue line denotes the greatest p-value within the top 100 significant single-nucleotide polymorphisms (SNPs) at −log10 = 3.59.The region around 8,150,383-81,931,753 bp on BTA 10 held three genes: Solute Carrier Family 10 Member A1 (SLC10A1), Serine and Arginine-rich Splicing Factor 5 (SRSF5), and Solute Carrier Family 8 Member A3 (SLC8A3).Along with this, the regions around 10,026,821 bp and 29,252,315 bp on BTA 18 and 27, respectively, were linked to two genes: Cadherin-13 (CDH13) and Ring Finger Protein (RNF122).These regions on BTA 10, 18, and 27 have previously been reported as QTLs related to dystocia in Holstein cattle.The regions around 2,461,770 bp, 83,245,103 bp, 82,066,559 bp, and 92,480,213 bp on BTA 8 (n = 1), BTA 10 (n = 2), and BTA 11 (n = 1), respectively, were linked to four genes: Solute Carrier Family 24 Member A2 (SLC24A2), Cytochrome C Oxidase Assembly Factor (COX16), Signal Induced Proliferation Associated1 Like 1 (SIPA1L1), and DAB2 Interacting Protein (DAB2IP), all reported regions containing previous QTLs associated with stillbirth.We identified regions around 66,416,290 bp, 65,237,605 bp, and 34,170,987 bp on BTA 4, 17, and 25, respectively, linked to four genes: WAS/WASL Interacting Protein Family Member 3 (WIPF3), Secerin-1 (SCRN1), Uncharacterized Protein KIAA1671 (KIAA1671), and Cytochrome p450 Oxidoreductase (POR), each of which is QTLs previously reported relating to calving ease in dairy cows.The regions around 33,569,134 bp, 8,809,053 bp, and 70,481,417 bp on BTA 12, BTA 15, and BTA 17, respectively, were linked to three genes: ATPase Phospholipid Transporting 8A2 (ATP8A2), Contactin-5 (CNTN5) Signal, and DEP Domain-containing 5 Protein (DEPDC5), reported previous as QTLs associated with non-return rates in dairy cows.Ten significant SNPs spanning 87,619,311-87,650,109 on BTA 10 corresponded to one gene, Estrogen-related Receptor β (ESRRB), and this location housed a primary QTL for the occurrence of placental retention in the post-calving period.The remaining SNPs, rs110783124 and rs135585624, were not previously reported to any fertility QTL specific to Holstein or Jersey cattle.

Figure 2 .
Figure 2. Genome-wide association study (GWAS) results from 624,460 SNPs and 1967 dairy cattle across three herds diagnosed with at least one disease in the metritis complex.Manhattan plot for −log10 p-values of SNP effects for the metritis complex in Holstein and Jersey dairy cows.The red horizontal line corresponds to a p-value threshold of −log10 ≥ 4.00.The blue line denotes the greatest p-value within the top 100 significant single-nucleotide polymorphisms (SNPs) at −log10 = 3.59.

Table 1 .
Descriptive statistics for the metritis complex in Holstein and Jersey multiparous cows.