DNA Methylation Profile Associated with Fertility Trait in Goat Using Whole-Genome Bisulfite Sequencing

DNA methylation is explained as attachment of methyl group to DNA especially to the carbon atom of cytosine base C at 5th position (Blow e al., 2016). Also, the addition of the methyl group may be occurred to the adenine nucleotide (Greer et al., 2015). The C methylation in CpG context is the extensively examined one compared to other contexts; CHG or CHH (An et al., 2018). DNA methylation is considered an important epigenetic factor which can transmit through successive generations (Ming et al., 2021) and it has different important physiological and biological activities such as gene expression and cellular differentiation (Vento-Tormo et al., 2016; Adamkova et al., 2017) and can affect the phenotypes characters in domestic animals (Fan et al., 2020).

DNA methylation can be studied using different sequencing approaches including whole genome bisulfite sequencing (Zhang et al., 2017). The Next-generation sequencing of DNA treated with bisulfite became the most popular approach to investigate the DNA methylation profile for the whole genome (Jang et al., 2017) where it gives complete information about the methylation profile of cytosine nucleotides among the whole genome with high specificity.

As the levels of economic and living standard increasing, more meat food is required. Lamb meat is preferable and palatable for population in Mediterranean’s countries including Egypt (Teixeira et al., 2020). The demand for lamb is much more than the real lamb production in Egypt and this demand for lamb is continuously increasing due to the over populations. Goats are considered the important sources for meat and milk, particularly in dry and semidry regions (Abd-Allah et al., 2019).

Fecundity is an important economic trait for production of goat meat, therefore, the increasing of goat fecundity will dramatically increase the economic benefits (Ahlawat et al., 2015). The litter size is an important indicator for fecundity, and it contributes 74-93% for fecundity. The ovulation rate is the direct and most important factor which controls the litter size and fecundity where it is the major contributor for goat fecundity. In order to breed the higher fecundity goat, it is highly important to investigate the goat ovulation rate trait (Lai et al., 2016), so the present work aimed to investigate the DNA methylation pattern of Zaraibi goat ovaries and its association with fertility trait using the sequencing of bisulfite-treated DNA.

Two groups of Zaraibi goats were used in this study; high fertility group (HFG; goats had 3/litter) and low fertility group (LFG; goats had 1/litter) in the four previous generations. Three animals from each group were sacrificed under the normal conditions (so, it is not needing any ethic permission) and the ovaries were collected for the analysis of DNA methylation.

WGBS was carried out according to procedures described by Frattini et al. (2017) and An et al. (2018). Genomic DNA was extracted from the collected ovaries of two fertility groups. The genomic DNA was fragmented to 200-300bp and treated with bisulfite using DNA Methylation Kit. After treatment, PCR amplification of DNA was done to constitute DNA library. The sequencing of samples was done using Illumina HiSeq 2500 platform. After sequencing and filtration, the clean reads were recorded and compared to goat reference genome (Dai et al., 2010).

The methylation level analysis was done according to procedures described by Frattini et al. (2017) and An et al. (2018). The clean reads were compared with goat reference genome (Kasper et al., 2012) and the methylation level was recorded in the complete genome.

The DMRs were determined using methylKit software. Multiple samples data were analyzed together. Logistic regression model was used to analyze the DMR between groups. The DMR differences were expressed as the mean for each group. The significant different methylation levels between the two tested groups were determined (>10%; q value <0.05).

The functional enrichment of Differentially Methylated Genes was done using Metascape online (https://metascape.org/gp/index.html#/main/step1) (Zhou et al., 2019).

To validate the sequencing results of these DMGs, ten genes of them were selected for the assessment of their mRNAs expression level using qRT-PCR. The reaction volume of 20 μl containing 10 μl SYBR Green master mix, 0.6 μl of 10 μM of both primers (Zhang et al., 2017), 1 μl cDNA and the volume was completed with RNase-Free water. The Data were analyzed using the equation described by Livak and Schittgen (2001). GAPDH was used as a House-Keeping gene to normalize the results of gene expressions and the statistically significant was p<0.05

Whole genome DNA bisulfite sequencing was done for six samples: three from each groups, high (HFG) and low (LFG) fertility groups. Averages for total numbers of sequenced basses and reads number were 177,265,899,051 and 1,173,946,351 for HFG and for LFG 177,437,595,615 and 1,175,083,415, respectively. After trimming process, the total numbers of sequenced bases in HFG and LFG were 153,975,467,691 and 153,874,094,771, respectively and total read›s number were 1,168,585,920 and 1,169,067,815 for HFG and LFG, respectively. The means of Q20% and Q30% for clean and full-length reads were 97.73 and 93.72 for HFG and 97.56 and 93.35 for LFG, respectively (Table 1).

Table 1: DNA methylation data

Table 2: Estimated bisulfite conversion rate by sample

Table 3: Methylation ratio of each cytosine in CG, CHG and CHH contexts

Bisulfite conversion rate was calculated at CG, CHG, CHH sites within the lambda genome (Table 2). The bisulfite conversion rates were estimated and ranged from 99.79% to 99.81% among samples using lambda phage DNA as a spike-in control.

The results declared that there is a small difference in the DNA methylation profiles between the two groups, where the frequency of methylated C to total C in context CG is 74.43% for HFG and 73.75% for LFG. The frequency of methylated C to total C in context CHG is 0.39% for HFG and 0.38% for LFG and the frequency of methylated C to total C in context CHH is 0.38% for HFG and 0.37% for LFG (Table 3).

The distribution ratios of methylated C in the three contexts: CG, CHG and CHH were 89.89%, 2.39%, 7.72% and 90.19, 2.34%, 7.47% in high and low fertility groups, respectively (Figure 1). This finding showed that the level of methylated C in CG context is in directionally opposite with fertility trait where this level is lower in HFG compared with that in LFG. In contrast, the levels of methylated C in contexts CHG and CHH (where H = A, C or T) are high in HFG compared with those in LFG groups

The differentially methylated regions (DMR) with different methylation levels in the two fertility groups (HFG and LFG) were identified. Then, differentially methylated region-related genes (DMG) were determined; where DMR exists within the gene and its upstream 2 kb range, it is considered as a DMR-related gene. Totally, 4640 DMRs were observed, and their analysis showed that most of the

Table 4: List of Differentially Methylated Genes (DMGs)

Table 5: List of fertility-associated DMGs

DMRs (4615 DMRs) were found in CG context whereas little DMRs were found in CGH (15 DMRs) and CHH (10 DMRs) contexts. The quantity of hypermethylation in CG islands was determined in both HFG and LFG (Figure 2). For HFG, about 65.56% of hypermethylation CGI was present in distal intergenic regions followed by exons and interons with 16.20% and 12.40% of hypermethylation CGI, respectively. The fewer hypermethylation CGI were distributed in promoters (4.12%), 3›UTR (1.20%) and 5›UTR (0.52%). With no significant difference with HFG, the hypermethylation CGI were distributed in LFG as follows: 65.45% in distal intergenic regions, 15.85% in exons, 12.90 in introns, 4.01 in promoter, 1.30 in 3›UTR and 0.49% in 3›UTR.

The enrichment analysis declared that some of DMGs are involved in the reproductive system development and many functional pathways which are related to follicle development and other important functional pathways responsible on regulation of growth, hormone production and development. The data declared that DMGs in goat ovaries are associated with the fertility phenomena in both fertility groups. Five-hundred and thirty-six DMGs (Table 4) from these DMRs were further analyzed in both fertility groups into HM (hyper-DMGs in HFG group) and LM (hyper-DMGs in LFG group where there were 1505 and 1135, respectively. The protein-coding regions (CDS) in these DMGs were determined (Figure 3) and enrichment analysis of protein coding regions declared their involvement in different important pathways (Figure 4).

One-hundred and seventy fertility-related genes with different methylation levels were selected for further functional enrichment analysis and the results declared that these genes with their methylation patterns are responsible for the difference in the fertility phenomena between high and low fertility groups (Table 5).

Ten DMGs were randomly selected to assess their gene expression for validation the association between their methylation patterns and the fertility trait. Data were represented as the fold change in target gene expression normalized to GAPDH gene in HFG relative to LFG. The mean cycle threshold (Ct) values of triplicate samples for each gene were calculated and used in data analysis. The results declared that the expressions of five tested DMGs: SCYL1, mTOR, ABCG2, STK3 and PSMD7 were higher in HFG than those in LFG, with significantly level (p<0.01) for SCYL1, mTOR, ABCG2 genes whereas the expressions of STK3 and PSMD7 genes were insignificant increasing (Figure 5). The expressions of other five genes were lower in HFG than those in LFG: GPNMB, ELK4, BACH1 at significant level p<0.01, CDIPT with significant level p<0.05 and ACVR1 with insignificant value. Generally, the qRT-PCR results were concordance with the sequencing data and confirming the reliability of sequencing data.

The domestic goat is an important small ruminant which distribute all over the world, especially in desert and semi-desert regions where goat is easier to breed with less expensive than other large livestock (Nowier et al., 2019). The improvement of goat productivity has a respectable role for covering the gap between over-growing populations and adequacy of milk and meat supplies (Miller and Lu, 2019). Production improvements of farm animals can be achieved by using new genetic technology for better selection of heritable traits. Selection has done depending on economically important productivity traits such as disease resistance, fertility performance as well as meat and milk production (Gutierrez-Reinoso et al., 2021).

The litter size represents one phenomenon associated with fertility performance in goats and it is related to the healthy reproductive organs and high ovulation rate (Assan et al., 2021). DNA methylation is an epigenetic factor which has an important role in gene expression regulation. In farm animals, DNA methylation regulates sexual, ovarian maturation and ovulation rate (Breton-Larrivée et al., 2019). Genome-wide bisulfite sequencing was used to explore the methylation profiles of genes which are related to complex trait such as litter size or ovulation rate trait (Lai et al., 2016).

The genome-wide methylation levels were assessed in different farm animals which discussed the association between DNA methylation profiles, gene expression and phenotype characters related to production and reproduction traits (Hao et al., 2016; An et al., 2018; Fan et al., 2020). Recently, the DNA methylation profiles in sheep ovaries and their associations with different parameters of fertility trait were described (Barboni, et al., 2011; Russo et al., 2013; Zhang et al., 2017). Regarding to goat, An et al. (2018) studied the change in the whole-genome methylation profiles between oestrous and dioestrous stages. The present study was focused on the assessment of methylation pattern in the whole genome goat ovaries and its relationship with the difference in litter size and ovulation rate traits. Our results showed a slightly non-significant increase in DNA methylation level in goat ovaries in high fertility group over than that in low fertility group despite the detection of many DMRs and DMGs among the two fertility groups. The similar finding was reported by Kang et al. (2022) who studied the genome-wide DNA methylation pattern in goat ovaries of Chinese Ji›ning grey goats with different litter size groups.

In the present study, most of methylated C was present in CG context and this result agrees with other reports in different species including human (Lister et al., 2014), pig (Hao et al., 2016), sheep (Zhang et al., 2017) and Chinese Ji›ning grey goats (Kang et al., 2022). Our finding declared that the major hypermethylation CGI in Zaraibi goat breed was present in distal intergenic regions followed by exons and interons whereas the lowest frequencies of CGI were in 3›UTR and 5›UTR. These proportions of methylated CGI were nearly like those which were reported in Chinese sheep and goat (Zhang et al., 2017; Kang et al., 2022).

One hundred and seventy DFGs were selected for the enrichment analysis and the result showed their role in different biological functions related to fertility. The directionally opposite association between methylation levels of ten selected DMGs and their gene expression was confirmed where SCYL1, mTOR, ABCG2, STK3 and PSMD7 were hypomethylated with highly expression in HFG compared with LFG. On the other hand, GPNMB, ELK4, BACH1, CDIPT and ACVR1 were hypermethylated with decreased expression in HFG than those in LFG. This directionally opposite association between methylation levels of DMGs and their gene expression was also reported by Kang et al. (2022) in goat and by Zhang et al. (2017) in sheep.

It is concluded that DNA methylation patterns of Zaraibi goat ovaries may be responsible for the fertility trait in goats through their important roles in folliculogenesis, oocyte ovulation rate and finally the fertility phenomena. Also, this study confirmed the opposite association between the methylation levels and expressions of differentially methylated genes which are related to fertility phenomena in goats.

The authors are grateful for the financial funding by Science, Technology & Innovation Funding Authority (STDF) under grant (STDF 30352).

There is absolutely no conflict of interest between the authors of this manuscript and any other scientists or producers.

This work was financial funded by Science, Technology & Innovation Funding Authority (STDF) under grant (STDF 30352).

To our knowledge, this study is the pioneer one for providing a comprehensive analysis of whole-Genome DNA methylation patterns in goat ovaries which helps for understanding the relation between ovarian DNA methylation and Egyptian goat fertility.

All authors certify that they have participated sufficiently in contributing of the manuscript. Othman E. Othman: Work design, Methodology, Data analysis and Manuscript writing. Lingjiang Min: Work design, Primer providing and Data analysis. Amira M. Nowier: Providing animals and all data associated with their fertility records.