INTRODUCTION
Over the past three decades, intensive genetic selection for high milk yield in dairy cattle has led to a pronounced energy deficit postpartum and a decline in fertility [1,2]. This occurs because negative energy balance (NEB) takes place when nutrient requirements for maintenance and lactation exceed dietary intake [3]. Cows experiencing NEB in late gestational and early lactation are particularly at risk of metabolic disorders due to utilization of body fat as a source of energy, leading to an elevation of circulating ketone bodies concentration in the bloodstream. Subclinical ketosis (SCK) is a prevalent metabolic disorder in high-producing dairy cows, characterized by increased concentrations of ketone bodies in the b1ood without clinical signs of ketosis [4]. The gold standard diagnostic test for SCK is the measurement of a beta-hydroxybutyric acid (BHBA) in serum/plasma [5,6] and cows with a BHBA concentration of 1.2 to 2.9 mmol/L are considered to have SCK [6]. The average incidence of SCK within the first weeks of lactation ranges from 26% to 56% [7,8]. Additionally, Cows with higher milk production were at increased risk for hyperketonemia (HYK), while cows with a body condition score (BCS) of 4 or higher before calving or those that lost more body condition during the transition period were more likely to develop HYK, emphasizing the importance of avoiding over-conditioning of cows during the dry period and excessive BCS loss during the transition period [9].
Interestingly, SCK has a negative association on reproductive performance, including an increased calving to first estrus, first insemination, and pregnancy intervals [2,10]. Specifically, SCK cows are less likely to become pregnant after first insemination and produced less milk [11-13].
Loor and colleagues reported that 2,415 genes were altered in liver from periparturient dairy cows undergoing nutrition-induced ketosis, and the genes were associated with oxidative phosphorylation, protein ubiquitination, cytokine signaling, fatty acid uptake/transport, fatty acid oxidation, cholesterol metabolism, growth hormone signaling, proton transport, and fatty acid desaturation [14]. Additionally, the proteome analysis in feed-deprived dairy cows revealed that altered proteins are involved in fatty acid oxidation, glycolysis, electron transfer, protein degradation, antigen processing, cytoskeletal rearrangement, and cholesterol transport [15]. However, our understanding of polymorphisms within genes linked to SCK, which could aid in the removal of cows with SCK-susceptible genotypes from the breeding population, remains limited. This is due to the low heritability range (0.01–0.16) of the trait and the challenges associated with detecting the phenotypic trait. Nonetheless, it’s noteworthy that a reference population, utilizing blood beta-hydroxybutyrate (BHB) concentration measurements, has the potential to estimate HYK breeding values for a genomic selection program [16,17].
In this study, we conducted a single nucleotide polymorphism (SNP)-based gene-set enrichment analysis (GSEA) to identify candidate genes associated with SCK. We first compared allele frequencies of each SNP between the SCK and non-SCK control groups using genome-wide association studies (GWAS). Subsequently, we used GSEA to better understand the association of genes with SCK susceptibility by identifying biological process or functional pathways. GSEA is a powerful analytical method for interpreting multiple contributing factors affecting complex traits and diseases that can be regulated by different combinations of mutations among different genes within biological groupings [18-20]. The identification of genetic markers and putative biomarker genes using GWAS-GSEA would contribute to successful breeding programs of dairy cattle and provide early and accurate diagnosis in terms of SCK susceptibility.
MATERIALS AND METHODS
All experimental protocols were approved by the Animal Care and Use Committee of the College of Agriculture and Life Sciences at the University of Wisconsin-Madison (ACUC no. A005802).
As previously described in other studies [9,17], Holstein cows were fed a corn silage and wheat straw based total mixed ration (TMR) during the dry period, and a TMR based on corn silage and alfalfa silage after calving. Blood samples were collected from a total of 112 Holstein cows two days per week between 5 and 18 days postpartum, during the morning feeding, to determine the incidence of SCK. Blood BHBA concentration was measured using an electronic hand-held biosensor, the Precision Xtra meter (Abbott Lab, Abbott Park, IL, USA), which is a useful cow-side ketone test with high sensitivity (91% to 94%) for detecting SCK in postpartum dairy cow blood samples, compared with laboratory assays [21]. Cows were diagnosed with SCK (hyperketonaemia) if blood BHBA concentration was ≥ 1.2 mmol/L. A binary SCK phenotype was assigned each cow. Of the 112 cows, 30 (26.7%) were diagnosed with SCK, and 82 (71.4%) were healthy.
Genomic DNA was extracted from hair follicles of cows and genotyped using the Illumina BovineSNP50K BeadChip according to the manufacturer’s instructions (Illumina, San Diego, CA, USA), which includes 54,609 Bos Taurus autosome SNPs. The Genome Studio Data Analysis software (Illumina) was used for visualizing SNP data and conducting preliminary analysis. Initially, we compared individual SNPs from SCK cows and healthy cows, and examined alternative alleles or variants at a given SNPs. SNP filtering was performed using PLINK 1.9 software with the following exclusion criteria: minor allele frequency (MAF) < 0.01; call rate < 0.10; and Hardy–Weinberg equilibrium (HWE) < 0.0001. After a quality control process, 43,552 SNPs were retained for further analyses. To identify genomic regions linked with SCK, we utilized a binary trait, assigning ‘1’ when blood BHBA concentration ≥ 1.2, and ‘0’ otherwise. The genetic relationship matrix (GRM) that we constructed was then employed in a GWAS using a following mixed linear model (MLM):
where yij is the binary SCK phenotype for a given cow with overall mean µ; β is the marker genotype effect, SNPi is the additive SNP substitution effect, aj is a random polygenic effect, and eij is a vector of random residual effects. Genome-wide p-values threshold of 10−5 was used for analyses of variants.
To assign SNPs to genes, we searched for SNPs that were detected within the genomic sequences of annotated genes (using dbSNP; https://www.ncbi.nlm.nih.gov/SNP/) or within 20 kb of the 5’ or 3’ ends of the first or last untranslated region (UTR), respectively. To explore functional enrichment in the gene sets and generate Kyoto Encyclopedia of Genes and Genomes (KEGG) biological pathways, we used Database for Annotation, Visualization, and Integrated Discovery (DAVID) Bioinformatics Resources 6.8 (http://david.abcc.ncifcrf.gov/). To confirm the predicted gene set function and interaction, we analyzed the listed genes in the KEGG pathways using the web-based tool Gene Multiple Association Network Integration Algorithm (GeneMANIA; https://genemania.org/). We targeted SNPs with a call rate of ≥ 95% and a MAF of ≥ 5% for prediction of gene sets in this study.
MAF was determined using the FREQ procedure of SAS version 9.2 (SAS Institute, Cary, NC, USA). The Benjamini-Hochberg (BH) false discovery rate (FDR) correction was applied to raw p-values (< 0.05). The distribution of genotypes was tested for deviation from HWE using a chi-square test. A p-value less than 0.05 was considered statistically significant.
RESULTS AND DISCUSSION
We conducted GWAS by analysing 54,609 SNPs and identified 8,360 loci that showed differences when comparing SNPs between the SCK and healthy control groups. We then eliminated SNPs with lower MAF frequencies, which were not considered to be genome-wide significant due to the sample size, and identified 194 SNPs associated with SCK in 163 functional candidate genes. This was achieved by evaluating all SNPs within each gene and the surrounding 1 Mb on each side of the gene for their suitability to serve as the gene’s proxy (Table 1 and Supplementary Table S1) [22].
To identify key genes and pathways associated with SCK, we performed KEGG pathway enrichment analysis using DAVID (https://david.ncifcrf.gov/). We found several pathways that may be involved in the incidence of SCK, including calcium signaling, starch and sucrose metabolism, cyclic adenosine monophosphate (cAMP) signaling, intestinal immune network for IgA production, and metabolic pathways (Table 2), although none of the evaluated KEGG pathways or Gene Ontology (GO) gene sets were significantly enriched for genes associated with diseases such as ketosis. Among these pathways, calcium signaling, starch and sucrose metabolism, and cAMP signaling pathways were statistically significant (p < 0.05). In the calcium signaling pathway, we identified several candidate genes including ATPase plasma membrane Ca2+ transporting 1 (ATP2B1), stromal interaction molecule 2 (STIM2), phospholipase C delta 3 (PLCD3), Ryanodine receptor 2 (RYR2), and prostaglandin F receptor (PTGFR). In the starch and sucrose metabolism pathway, three candidate genes were detected: UDP glucuronosyltransferase family 1 member A6 (UGT1A6), UGT2A1, and probable maltase-glucoamylase 2 (LOC1002966901). The cAMP signaling pathway included ATP2B1, T-cell lymphoma invasion and metastasis 1 (TIAM1), RYR2, phosphodiesterase 3A (PDE3A), and mitogen-activated protein kinase 10 (MAPK10). We then analyzed the genes listed in the KEGG pathways (Table 2) using GeneMANIA, which can extend the listed genes with functionally similar genes and predict their biological function/network. This approach showed that the listed genes are co-expressed (78.47%) and co-localized (13.73%). Interestingly, UGT1A6 and UGT2A1 may affect metabolism, such as glucuronidation.
Identification of genes affecting susceptibility to common disease is very difficult because each causal gene only makes a limited contribution to the diseases. Although ketosis is an important metabolic disease of dairy cows, and SCK occurs more often than clinical ketosis, little is known about the genes and signal pathways affecting SCK in the cattle. It has been reported that polymorphisms within the APOBR gene, which codes for the apoliporportein receptor, are associated with the milk level of a prognostic ketosis biomarker in dairy cows [23]. Another study has reported that a SNP in the 3́- UTR of protein kinase AMP-activated non-catalytic subunit gamma 1 (PRKAG1), a regulatory subunit of the AMP-activated protein kinase (AMPK) that plays a critical role in regulating cellular energy metabolism, can influence BHBA levels and milk production. This effects is due to the distortion of the target site of the highly expressed microRNA mir-423-5p in the bovine mammary gland, liver, and kidney [24]. These findings suggest that allele mutation in specific genomic regions may be related to the susceptibility of SCK. However, most studies have focused on evaluating the value of phenotypic traits such as the levels of BHBA and nonesterified fatty acids (NEFA) for predicting ketosis [16] rather than identifying specific variants associated with functional candidate genes that influence disease susceptibility.
A recent study has reported that a combined approach using GWAS, genome-wide interaction studies (GWIS), and metabolic pathway enrichment analyses identified SNPs and proximal candidate genes associated with susceptibility to HYK. This suggest that the combined analysis is an effective way to detect genetic factors contributing to HYK, and that the proposed genes are related to energy and lipoprotein metabolism, particularly insulin secretion or resistance [25].
In this study, we identified SNPs associated with SCK using GWAS and a pathway-based approach, leading to the discovery of functionally important genes. For example, SNPs closely located to ANAPC4 and SEC 23A are identified by using the combined analysis (GWAS-GSEA). This result is consistent with previous literatures which states that both genes are associated with nutrition-induced ketosis [14] and that expression level of these genes is altered in ketogenic diet rats, suggesting their involvement in mitochondrial biogenesis [26].
In addition, several genes listed in Table 1 are involved in metabolic dysfunction such as diabetes; regulator of calcineurin 1 (RCAN1) is highly expressed in type 2 diabetes in human and mice, and elevation of RCAN1 reduces β-cell mitochondrial function and ATP availability, resulting in a reduction of glucose-stimulated insulin secretion [27], ST8 alpha-N-acetyl-neuraminide alpha-2,8-sialyltransferase 1 (ST8SIA1), a membrane-bound glycosphingolipid, is reported to affect impaired glucose-stimulated insulin secretion [28], and an intronic SNP within leucine rich repeat and Ig domain containing 2 (LINGO2) genes is associated with body mass and adiposity in elder humans [29]. Moreover, analysis of integrated function/interaction using GeneMANIA supports our bioinformatics approach. For example, UGT1A6 and UGT2A1 are highly associated with metabolic dysfunction observed in the ketogenic diet in humans and high-fat diet-induced fatty liver in rat [30,31].
In this study, our analysis based on GWAS-GSEA suggested the involvement of calcium signaling, starch and sucrose metabolism, cAMP signaling, intestinal immune network for IgA production, and metabolic pathways in SCK. Currently, understanding of the interplay between the metabolic and immune systems is crucial to decipher the integration of metabolism and immunity in periparturient dairy cows. As previously mentioned, the periparturient dairy cow experiences a significant increase in nutrient requirements for lactation; roughly a threefold increase in glucose demand, a twofold increase in amino acids demand, and approximately a fivefold increase in fatty acids [32]. Furthermore, cows require a fourfold increase in calcium on the day of parturition, and early lactating cows experience negative calcium balance [33,34]. Intracellular calcium, in particular, must be maintained to avoid nerve and muscle dysfunction, as it is involved in many signaling pathways [35].
Additionally, NEB during early lactation may worsen periparturient immunosuppression, making dairy cows more vulnerable to infectious diseases [36]. HYK, in particular, appears to negatively affect immune function, impairing udder defense mechanisms related to leukocyte function in NEB cows [37]. Although the mechanisms of impairment have not been fully elucidated, the high concentration of ketone bodies has been associated with suppressed bovine lymphocyte blastogenesis, lowered chemotactic capacity of leukocytes, and decreased production of interferon (IFN)-gamma and tumor necrosis factor (TNF) alpha [38-40].
Moreover, we identified several candidate genes that are related to reproductive performance, including early embryo development, implantation, and germ cell development. For instance, expression of Ectodysplasin-A receptor-associated adapter protein (EDARADD) and LINGO2 affect early cleaving embryos development and oocyte maturation, respectively [41-43]. SMG6 Nonsense Mediated mRNA Decay Factor (SMG6) depleted blastocyst embryos did not form the ICM, leading to early embryonic lethality in mice [44]. Additionally, ATP2B1 and GA Binding Protein Transcription Factor Subunit Alpha (GABPA) may affect spermatogenesis via sperm calcium channel and folliculogenesis via granulosa cell interactions, respectively [45,46]. Importantly, a SNP located within an intron of Mucin glycoprotein 19 (MUC19) has been reported to be involved in early embryonic loss during the implantation window of heifer Holstein cows [42].
In summary, our study identified chromosomal regions associated with the incidence of SCK and their corresponding genes, highlighting the relationship between nutritional metabolism and immunity in periparturient dairy cows through GSEA. Our approach using SNP and GSEA provided integrated information for the identification of genes associated with SCK. Notably, we also found that the proximal genes of SNPs are related to gametogenesis, early preimplantation development, and the implantation window, suggesting potential links between SCK and reproductive performance. This genetic information can be utilized for optimizing breeding programs and managing high-performance cows with infertility, thus minimizing economic losses caused by SCK.
SUPPLEMENTARY MATERIALS
Supplementary materials are only available online from: https://doi.org/10.5187/jast.2023.e97.
