Abstract
The way in which living organisms mobilize a combination of long-term adaptive mechanisms and short-term phenotypic plasticity to face environmental variations is still largely unknown. In the context of climate change, understanding the genetic and epigenetic bases for adaptation and plasticity is a major stake for preserving genomic resources and the resilience capacity of livestock populations. We characterized both epigenetic and genetic variations by contrasting 22 sheep and 21 goats from both sides of a climate gradient, focusing on free-ranging populations from Morocco. We produced for each individual Whole-Genome Sequence at 12X coverage and MeDIP-Seq data, to identify regions under selection and those differentially methylated. For both species, the analysis of genetic differences (FST) along the genome between animals from localities with high vs. low temperature annual variations detected candidate genes under selection in relation to environmental perception (5 genes), immunity (4 genes), reproduction (8 genes) and production (11 genes). Moreover, we found for each species one differentially methylated gene, namely AGPTA4 in goat and SLIT3 in sheep, which were both related, among other functions, to milk production and muscle development. In both sheep and goats, the comparison between genomic regions impacted by genetic and epigenetic variations suggests that climatic variations impacted similar biological pathways but different genes.
Introduction
The evolution and distribution of species is driven by the variation of their environment. Optimization of the individual’s phenotype to fit the environment may occur on the long-term, as populations adapt to the local values of environmental drivers (i.e., native individuals having on average a higher fitness than migrants, Savolainen et al., 2013). The access to whole genome sequences now allows to study the genomic bases of local adaptation by identifying genes and genomic regions under selection and the environmental parameters responsible for their selection (e.g., honeybees in the Iberian climate Henriques et al., 2018; dogs with the elevation Gou et al., 2014; or goats in the Moroccan climate Benjelloun et al., 2015). In the meantime, the context of global changes lead to an increasing interest on the ability of populations to develop short-term responses to face e.g., climate variations, as well as the interaction between short- and long-term mechanisms. Indeed, populations might quickly react to an environmental change through migration towards more favorable environments, which is well documented in the global warming context (Chen et al., 2011; McDonald et al., 2012). In the short-term, phenotypic plasticity can also increase the ability of organisms to cope with environmental changes. The ability of a genotype to produce, under different conditions, different phenotypes that are best fitted to the environment involves non-genetic mechanisms by which favorable variations can be acquired and even sometimes transmitted. They consist in epigenetic variations, parental effects, ecological and cultural variations (Danchin et al., 2011). One of such mechanism is epigenetics, notably molecular modifications of chromatin without modification of the DNA sequence itself (including e.g., DNA methylation or histone modification) which affects gene expression (Gibney and Nolan, 2010). Until now, most of the studies pointing out the role of molecular epigenetic mechanisms were performed in stress-controlled conditions on plant or animal models. In this case, the effect of one varying parameter on the physiological or epigenetic responses of the organisms is assessed. For example, the glyphosate herbicide injury on Arabidopsis thaliana (Kim et al., 2017), or the maternal diet in mice (Cooney et al., 2002) were shown to affect the methylome.
Since the end of the 2000s, the role of these epigenetic processes in ecology and evolution is increasingly being studied (Bossdorf et al., 2008). In this context, epigenetic marks were analyzed in natural populations. For example, the search for epigenetic differentiation by contrasting natural populations of Lilium bosniacum living in different habitat conditions (Zoldoš et al., 2018), or in wild baboons (Papio cynocephalus) populations with different food resources (natural fodder in a savanna environment or human food scraps, Lea et al., 2016).
Like wild species, farm animals will have to face future climate change, and understanding the genetic and epigenetic mechanisms responding to environmental changes, especially for free-ranging populations, is a key issue for the conservation of Farm Animal Genomic Resources (Bruford et al., 2015). Our study focuses on small domestic ruminants, namely sheep (Ovis aries) and goats (Capra hircus), which are a very good model for studying local adaptation. Indeed they spread all over the world for about 10,000 years from the domestication center in the Middle East, through human migration and commercial trades (Zeder, 2008). Until recent centuries, they were traditionally managed in a sustainable way leading to populations well adapted to the large variety of environments worldwide (Taberlet et al., 2008). Nowadays, the effect of the environment on farm animals strongly depends on breeding conditions. Intensively farmed animals are rather preserved from environmental variation because they are generally kept indoors with access to quality supply, in contrast to extensive livestock which are raised outdoors with little additional water and feed. Since the rise of the breed concept about 200 years ago, the selection for morphological and productivity traits improved the productivity, but resulted in the loss of genetic variability in many breeds, compromising their ability to adapt to future environmental changes (Taberlet et al., 2008). Thus, in the context of climate change gathering “epigenomic information to be integrated with phenotypic and genomic data to scrutinize the biological basis for adaptation and plasticity/resilience in livestock populations” is a key question in farm animal genomics research (Bruford et al., 2015). Until now, epigenetic analyses on livestock have been performed almost exclusively in stress-controlled conditions in relation with agronomic traits such as reproduction (Lan et al., 2013), production (Peñagaricano et al., 2014; X. Wang et al., 2015) and disease sensibility (Doherty et al., 2014). To our knowledge, only the study by Sevane et al. (2018) searched for DNA methylation in relation with climate variation by comparing cattle adapted to tropical and temperate climates.
In this context we aimed at characterizing both epigenetic and genetic variations with regards to climate variation in a single framework. Our study focused on sheep and goats sampled in Morocco where farming is mainly carried out in traditional ways with a strong impact of environmental conditions upon animals. Within each species, our goal was to identify candidate genomic regions involved in both genetic (selection on Single Nucleotide Polymorphism, SNP) and epigenetic (differential DNA methylation) mechanisms related to environmental variations, using whole genome comparisons of groups of individuals from locations with contrasted annual temperature range.
Materials and Methods
Sampling
Morocco has very contrasted climates including hot/cold desert and semi-arid, warm-summer Mediterranean, dry-summer subarctic (Born et al., 2008). In this context, we sampled sheep (Ovis aries) and goats (Capra hircus) from locations with contrasted annual temperature variations (Figure 1, bioclimatic variable BIO7 from the WorldClim dataset, www.worldclim.org/bioclim, Hijmans et al., 2005). This variable is the difference between the maximum temperature of the warmest month and the minimum temperature of the coldest month. It is correlated with several temperature variables, sunshine and elevation (Supplementary Table S1). Twenty-one unrelated goats (10 and 11 from low and high annual temperature ranges, respectively) and 22 unrelated sheep (12 and 10 from low and high annual temperature ranges, respectively) were selected (Figure 1; Supplementary Table S2) from the samples collected between January 2008 and March 2012 in the frame of the NextGen European project (Grant Agreement no. 244356), in accordance with ethical regulations of the European Union Directive 86/609/EEC, and for which Whole Genome Sequences were already available (projects.ensembl.org/nextgen/). The average annual temperature range between the two groups is 20 vs. 39.5°C, which is almost twice for the animals living in Eastern Morocco.
FIGURE 1
Variation of the annual temperature range and sampling of small ruminants in Morocco. Distribution of sheep and goats sampled at each end of the temperature gradient in Morocco. The color gradient corresponds to the variation of the annual temperature range (BIO7).
Whole Genome Sequences
Dataset
We extracted the states of SNPs from the variant files generated by the NextGen project. Whole Genome Sequences (WGS) at 12–14 X coverage corresponding to the studied individuals were retrieved from the European Nucleotide Archive (www.ebi.ac.uk/ena, accession code PRJEB7436). The sequences, produced and filtered following the protocol in Alberto et al. (2018), were aligned on the reference genomes CHIR v1.0 (GenBank assembly GCA_000317765.147, Dong et al., 2013) for goats and OAR v3.1 (GenBank assembly GCA_000298735.1, Jiang et al., 2014), for sheep.
Detection of Selection Signatures
We kept SNPs with a missing data rate inferior to 30% and a minor allele frequency above 1%. Then we calculated the FST value (Weir and Cockerham, 1984) between the two groups (low vs. high temperature annual range for each species) for each SNP using the weir-fst-pop function in vcftools version 0.1.16 (Danecek et al., 2011). Candidate genomic regions were defined by 1) merging SNPs with the highest FST values (i.e., top 0.025% and top 0.001% for the pathway and gene approaches, respectively) closer than 5,000 bases from each other and then 2) adding 2,500 bases at each side of the region obtained. At the end 3) the genes overlapping these regions were extracted.
To annotate the genomic regions, two gene ontology analyses were performed with the GO Ontology database released on 2020-08-10 (
Carbon and Mungall, 2018). The
Bos taurusreference was used to annotate the candidate genes. We do not expect any major issue due to the use of
B. taurusas a reference as the genome synteny between these three Bovidae species is very high (see
Jiang et al., 2014).
1) Enrichment analysis: it was performed on the list of genes extracted with the 0.025% threshold. For that, we used the statistical over-representation test in PANTHER (version 15.0, Mi et al., 2019) and the Bos taurus reference. Each biological process with a significant Fisher’s exact test and a False Discovery Rate (FDR) lower than 0.05 was annotated. We kept the more basal significant biological pathway and searched the GO terms with the highest hierarchical level (i.e., child term of “biological process”) associated.
2) Candidate gene identification: we extracted the GO term corresponding to the top genes using a more stringent threshold (0.001%) to avoid the false positive and highlight the most impacted genes. In addition, a specific gene bibliography was done using Google Scholar with the gene’s name and “goat/sheep” or “livestock” or “mammals” as keywords, to determine what could be their effect on livestock phenotypes. These were grouped into large phenotype categories such as “Environmental perception,” “Immunity,” “Production” (genes involved in milk or meat production), “Reproduction” and “Other.”
DNA Methylated Regions
DNA Extraction, Library Preparation and Sequencing
DNA was extracted from ear biopsies collected by the Nextgen consortium (see above) using the isolation of genomic DNA from tissues protocol from QIAamp® DNA micro kit (QIAgen, Germantown, MD, United States). Methylated DNA immunoprecipitation sequencing (MeDIP-Seq) was performed at the GeT-PlaGe core facility, INRAE Toulouse (http://www.get.genotoul.fr). Methylated DNA libraries were prepared according to Bioo Scientific’s protocol using the Bioo Scientific NEXTflex™ Methyl-seq Library Prep Kit for Illumina Sequencing. Briefly, DNA was fragmented by sonication on a covaris M220, size selection was performed using AMPure XP beads and adapters were ligated to be sequenced. 3 µg of each library was diluted in 450 µl of TE buffer, denatured in boiling water for 10 min and immediately cooled in ice for 10 min. Then 50 µL of 10X concentrated IP buffer was added to the mixture as well as 1 µg of anti-5-methylcytosine monoclonal antibody (clone 33D3, Diagenode®). After overnight incubation at 4°C with agitation, the DNA-antibody complexes were purified using 40 µL of “Dynabeads-ProteinG” previously washed in PBS. The DNA-antibody complexes were washed twice in 700ul of IP buffer and then the DNA was purified using the iPure kit (Diagenode®). Control or input DNA samples were not coupled to the antibodies and were therefore not enriched in epigenetic markers. Then, 12 cycles of PCR were performed. Library quality was assessed using an Advanced Analytical Fragment Analyzer and libraries were quantified by QPCR using the Kapa Library Quantification Kit. MeDIP-Seq was performed on an Illumina® HiSeq 2500 with the Illumina Reagent Kits. At the end, the sequenced reads correspond to fragment of methylated DNA.
Dataset Preparation
Illumina paired-end reads (mean size of 125 bp) from all animals were aligned on the same reference genomes as WGS data by using the default parameters of BWA-MEM version 0.7.12 (Li and Durbin, 2009). The function multicov from bedtools version 2.29.2 (Quinlan and Hall, 2010) was used to report the MeDIP-Seq coverage of alignment files inside windows of 100 non-overlapping base pairs along the genome. The same function was used to extract the WGS coverage. In these windows we also extracted the number of CpG from the reference genomes by using the MEDIPS package (Lienhard et al., 2014) on R software version 3.5.3 (R Core Team, 2019).
We compared genome wide coverage profiles of MeDIP-Seq by using the MEDIPS.correlation function from the MEDIPS package with the uniq, extent, and shift parameters equal to 0 and the window size parameter equal to 100, to keep all MeDIP-Seq reads in each 100 bases window on the genome. We verified with a pairwise Pearson test that the number of read by windows for one animal had a correlation coefficient greater than 0.6 with all other individuals of the same species.
The numbers of WGS and MeDIP-Seq reads in each window were normalized between individuals with default parameters of calcNormFactors function (edgeR package, Robinson et al., 2010). Then, normalized MeDIP-Seq read counts were divided by normalized WGS read counts in order to get rid of the effect of sequencing coverage variations (e.g., due to repeated regions). Normalized values have been rounded to two decimals. Four outlying windows (2 for each species) with 10 times more MeDIP-Seq reads than the others were removed as this may be due to PCR or alignment errors.
Whole Epigenome Correlations
For each species, correlations between the number of CpG and that of MeDIP-Seq reads found in the 100 bp window were tested using a Pearson’s correlation test. We used the same test to test the correlation between species for each of these two variables. For that we used the correspondence between Capra and Ovis orthologous regions resulting from the cross-alignment between the two reference genomes CHIR_1.0 and OAR v3.1 performed by Alberto et al. (2018). Results were represented as density plots performed with the hexbin package (Carr et al., 2019).
Detection of Differentially Methylated Regions
We used the edgeR package to detect DMRs. Starting from non-standardized MeDIP-Seq read counts, we kept the 100 bp windows for which at least 10% of the samples had at least 1 count-per-million (1 read per million windows, this accommodates differences in library sizes between sample) and we normalized the read counts with the calcNormFactors function with default parameters. We estimated common dispersion with the estimateGLMRobustDisp function which mitigates the influence of outliers (X. Zhou et al., 2014) and we tested the occurrence of differential methylation between groups with the default parameters of the exactTest function. The p-values obtained were converted into q-values with the Bioconductor package qvalue (Storey et al., 2015) and the ones showing q-values inferior to 0.01 were selected for the rest of the analysis. Finally, we merged significant windows distant than less than 5,000 bases into regions.
In order to delineate the peaks of methylation occurring for some of the individuals within these regions, we used the changepoint package (Killick and Eckley, 2014) which detects breakpoints in time series. We looked for all possible breakpoints within a region using the cpt.meanvar function that identify change in mean and variance of the sum of coverages for all animals of each group. We kept the region between two breakpoints that contained at least one window detected as significant with edgeR, thus, the delimited regions were merged with the bedr package (Waggott et al., 2017) when they overlapped. Then, we used a Wilcoxon test to detect peaks with different methylation levels between the two groups. Among those peaks, we discarded those where the number of MeDIP-Seq reads was related to the DNA coverage in the WGS data, because here high MeDIP-Seq counts would result from the high DNA coverage (e.g., due to repeated sequences). The regions containing the remaining peaks were mapped on the reference genomes and included or overlapping genes were identified. We used the Functional classification viewed in gene list analysis in PANTHER to extract the Biological Pathways annotated in the Gene Ontology Database associated to these genes.
Results
Detection of Selection Signatures
The FST values at each genome position and their overall distribution are presented for goats and sheep in Supplementary Figures S1, S2, respectively. Almost 94,5% of the FST values are less or equal to 0.1 for both species. The top 0.025 and 0.001% thresholds correspond to FST value of 0.4 and 0.55 respectively.
Enrichment Analysis
The 5,000 SNPs for goat and 6,000 SNPs for sheep which harbored the highest 0.025% FST values (Supplementary Figures S3, S4) corresponded to 442 and 489 genes, respectively, putatively under selection.
The enriched pathways were associated with 8 and 10 GO terms in goat (Figure 2, light green color; Supplementary Table S3) and sheep (Figure 2, dark green color; Supplementary Table S4), respectively. Eight were common to both species.
FIGURE 2
GO terms putatively involved in the adaptation and acclimation of small ruminants to temperature variations. Number of genes or child GO terms per basal GO terms for candidate gene identification (under selection -oranges- and differentially methylated -yellow/brown-) and enrichment analysis (greens), respectively.
Among the 0.025% top genes (i.e., 931 genes selected for the enrichment analysis), 41 were orthologous between sheep and goats. When we did the same enrichment analysis (described above) with the 41 genes shared by sheep and goats we obtained 21 GO terms. Four of them had the “Nervous system development” in their parental terms; and the last three were associated with cellular processes, one of which was related to synapses organization.
Candidate Gene Identification
The top 0.001% FST values corresponded to 194 SNPs for goat and 146 SNPs for sheep distributed in 56 and 42 genomic regions, respectively. Their location in the species genome is detailed in Supplementary Tables S5, S6 for goat and sheep, respectively. In goats, 34 of the 56 regions were intergenic and the others were related to 22 genes (Table 1). For sheep, we obtained 24 intergenic regions and 18 genes (Table 1), which differed between the two environmental groups. No genes were common between the two species.
TABLE 1
| Species | Phenotype categories | Genes |
|---|---|---|
| Goat | Environmental perception | TMTC2, EDIL3, SASH1 |
| Immunity | — | |
| Production | CP, SLC9A9, PPFIA2, TMTC2, GALNTL6, CTSB, SASH1, PLCG2, MYADM | |
| Reproduction | ADGRB3, USH1C, PLCG2, MYADM, SERPINB7 | |
| Other | TRABD2B, MVB12B, KLF12, TNFSF9, CALHM3, SORCS1, RAB30 | |
| Unknown | LOC102184299, LOC106503718 | |
| Sheep | Environmental perception | NOX3, KSR2 |
| Immunity | SAMD12, SEMA5A, SEZ6L | |
| Production | BMPER, KSR2, PTPRE | |
| Reproduction | GCSAML, KSR2, BMPER | |
| Other | AGO3, IQCJ, KCNG3, ANO6, CASP8AP2, CNKSR3, ENGASE, KCNA4 | |
| Unknown | LOC101107868, LOC101119001 |
Phenotype categories in which the candidate genes (top 0.001%) are involved.
The 22 genes were linked to 15 GO terms in goat (Figure 2, light orange color; Supplementary Table S5) and the 18 genes in sheep were associated with 14 GO terms (Figure 2, dark orange color; Supplementary Table S6). Twelve of these GO terms were common to both species.
In livestock, several of these genes are known to be involved in phenotypic changes in relation with the perception of the environment (3 genes in goat, 2 in sheep), immunity (1 in goat, 3 in sheep), production (9 in goat, 2 in sheep), reproduction (5 in goat and 3 in sheep). The “Other” category brings together the genes whose effect on the phenotype was not documented in livestock species. We also detected genes that were uncharacterized (Table 1, Supplementary Tables S5, S6).
Distribution of CpG and Methylated Regions in Both Species
The majority of the 22 million windows (i.e., 96% for each species) carried one methylated read, 3 and 2% more than one, and 1 and 2% none, for goat and sheep, respectively. Concerning the number of CpG, 47% of the windows did not have any, 50% had between 1 and 4, and 3% had more than 4 in both species.
Genomic windows of 100 bp that were orthologous between Capra and Ovis showed a correlation between species for the number of methylated read (R = 0.81, p-value < 2.10-16, Figure 3A) and the number of CpGs (R = 0.67, p-value < 2.10-16, Figure 3B) A large number of windows had few CpG and methylated reads (25% carried less than 3 CpGs and 1 methylated read). Within species, we found a correlation between the number of CpGs and that of methylated reads, R = 0.49 (p-value < 2.10-16, Figure 3C) for goats and R = 0.44 (p-value < 2.10-16, Figure 3D) for sheep.
FIGURE 3
Characterization and comparison of sheep and goats methylomes. Relationship between the number of (A). CpG, or (B). MeDIP-Seq reads per 100 bp window between the two species across the genome. Relation between the number of MeDIP-Seq reads and the number of CpG per 100 bp widow for (C). goat, and (D). sheep. The colors represent the number of points and the black line (y = x) is represented.
For the detection of DMRs, 20,256 windows in goat and 16,281 windows in sheep were kept out of more than 24 millions. The edgeR analysis revealed 20 and 7 windows with significant differences in methylation levels between the two groups (i.e., low vs high annual temperature variations) in goat and sheep, respectively. Merging candidate windows closer than 5,000 bp, we obtained 8 and 4 regions in which we detected 8 and 5 methylation peaks for goat and sheep, respectively. From these, 4 and 2 peaks respectively showed significant differential MeDIP-Seq coverage between the 2 groups (low vs. high temperature annual range). For 4 peaks, the differential MeDIP-Seq coverage was concomitant to a differential WGS coverage, indicating that the differentiation would probably result from repetition polymorphism. Then only one region per species was kept as resulting from a differential methylation. They were associated to AGPTA4 in goat and SLIT3 in sheep (Figure 4).
FIGURE 4
Differentially methylated regions putatively involved in temperature acclimation in small ruminants. DMRs between genomes from locations with low (blue) and high (red) temperature annual range for (A). goat or (B). sheep. Each line represents the variations of the level of DNA methylation along a gene for one individual; the grey bar represents part of the gene; the two black lines delineate the peak on which the differential level of methylation was tested. Boxplots represent the level of methylation (MeDIP-Seq coverage) for each group compared with a Wilcoxon test (Significance codes: ** = p-value ≤ 0.01, * = p-value ≤ 0.05).
The GO terms associated with these genes were “Cellular process,” and “Metabolic process” for AGPAT4 (Figure 2, yellow color); and “Biological regulation,” “Cellular process,” “Developmental process,” “Growth,” “Locomotion,” “Metabolic process,” “Multicellular organismal process,” “Reproduction,” “Reproductive process,” “Response to stimulus,” “Rhythmic process,” and “Signaling” for SLIT3 (Figure 2, brown color).
Discussion
To our knowledge, this is the first study to investigate simultaneously the genetic and epigenetic variations putatively involved in the adjustment of phenotypes to field conditions in small ruminants. Until now, few publications jointly addressed both mechanisms. Some studies highlighted significant association between local genetic variation and the presence of DMRs like that of Eichten et al. (2013) who found that 51% of the DMRs between maize populations were associated with local SNPs. Foust et al. (2016) and Herrera and Bazaga (2010) found a correlation between AFLP and MS-AFLP marks in two salt marsh perennials (Spartina alterniflora and Borrichia frutescens) and violet populations (Viola cazorlensis), respectively. In our study, we investigated both genetic and epigenetic variations in sheep and goats originating from field populations. The sampling was carried out in order to avoid a breed effect as much as possible. However, given the specific distribution and adaptation of some breeds, several Ghazalia goats and D’man sheep were part of the groups submitted to high temperature annual range (c.a. 39.5°C) as opposed to that with low temperature annual range (c.a. 20°C). They represented less than half of the individuals (3/11 and 4/10 in goats and sheep, respectively) and had similar methylation and FST patterns to the other breeds from the same group. Thus, the mechanisms inferred would be shared by the variety of Moroccan breeds represented in our sampling. Our study points out several genes which polymorphism or methylation level is contrasted between the two groups and may contribute to the adaptation and acclimatization to these specific environments. The candidate genes putatively under selection are primarily involved in environmental perception, immunity, production and reproduction, and their variations could be driven by different environmental factors. First, a direct effect of high temperatures might operate through thermoregulation by reducing body temperature, through the decrease of metabolic heat and the increase of heat dissipation. This might be obtained by avoiding solar radiations (Al-Tamimi, 2007), increasing respiration or sweating and reducing food intake (Marai et al., 2007; Sejian et al., 2010). These adjustments might reduce the energy allocated to other biological processes including productivity (Tao et al., 2011; Dado-Senn et al., 2019) or reproduction (Meyerhoeffer et al., 1985; Monteiro et al., 2016; Krishnan et al., 2017). In relation with temperature annual range, we found candidate genes involved in production, including milk production for CP, MYAMD, PTPRE, SASH1 and SLC9A9 (Hussein and Staufenbiel, 2012; Dong et al., 2015; Laodim et al., 2017; Shin et al., 2017; Xiang et al., 2017), mammary gland development for PPFIA2 (Mota et al., 2020) and morphology for TMTC2 (Abo-Ismail et al., 2017; Z.-H. Fang and Pausch, 2019). In addition, CTSB, GALNTL6, and KSR2 genes were related to meat production (Russo et al., 2002; Doran et al., 2014; Puig-Oliveras et al., 2014; Lukić et al., 2020), and BMPER and KSR2 with body size (Puig-Oliveras et al., 2014; Zhao et al., 2015; Lukić et al., 2020). Genes involved in reproduction were also impacted. ADGRB3, GCSAML, and KSR2 affect fertility (Pimentel et al., 2011; X. Wang et al., 2015; Mohammadi et al., 2020), while PLCG2 and SERPINB7 affect gestation (Vallée et al., 2003; Ponsuksili et al., 2012). Moreover offspring survival can also be impacted through litter size via BMPER and USH1C (Kwon et al., 2016; Bolormaa et al., 2017; Mota et al., 2020) and offspring weight with MYADM (Gonzalez et al., 2020). Also, we found two candidate genes related to environmental perception. NOX3 is involved in the response to stresses in general (Lambeth, 2004) and has been related with heat stress in chicken (Aggrey et al., 2018), and KSR2 plays a role in the energy balance by adapting feeding behavior and thermogenesis in mice (Guo et al., 2017). Second, environmental parameters covarying with the annual temperature range (BIO7) may be prominent in driving the selection on the genes which polymorphism was contrasted between the two groups. TMTC2 and EDIL3 are involved in adaptation to high altitude in sheep and cattle (Ji Yang et al., 2016; Pierce et al., 2020), this later variable being correlated with BIO7 in Morocco. In addition, BIO7 variable is correlated with sunshine which impacts directly skin pigmentation, which would explain the detection of SASH1, which regulates melanocytes migration in relation with the production and transfer of melanin to protect skin from UV radiations (D. Zhou et al., 2013). Third, the temperature regime may have indirect impacts on phenotypes by driving other environmental factors. Especially, it determines plant resources (i.e., food quality and quantity) and is related to water availability, and thus strongly impacts the diet of free ranging livestock. The lack of resources has a negative impact on, e.g., reproduction (Martin et al., 2004; J. J. Robinson et al., 2006) and lactation (Razzaghi et al., 2016). In contrast, supplementation, which has a positive effect on milk yield and composition in sheep, could temperate this negative impact on lactation (Caroprese et al., 2011). As a candidate gene, we found PLCG2 which is regulated in the cattle liver according to feed intake (Salleh et al., 2017). Likewise, genes involved in immune functions may be selected by multi-factorial drivers related to nutrition or the pathogenic context in relation with climate (Bett et al., 2017). This could be the case for SAMD12 playing a role in pancreatitis in human (Giri et al., 2016), SEZ6L associated with foot-and-mouth disease (Lee et al., 2015) or SEMA5A with mastitis in cattle (Sugimoto et al., 2006).
Besides the candidate genes putatively involved in local adaptation through allelic selection, we identified two strong candidates in relation with an epigenetic response. The global epigenomic context was similar in sheep and goats, and the level of DNA methylation along the genome was correlated to the CpG content. This is consistent with what is known about mammalian epigenomes where the CpG methylation generally covers the entire nuclear genome, with the exception of CpG rich regions near the promoters of active genes (Suzuki and Bird, 2008; Feng et al., 2010). We found two genes differentially methylated according to the provenance of the individual, namely AGPAT4 (1-acylglycerol-3-phosphate O-acyltransferase 4) in goats, and SLIT3 (Slit Guidance Ligand 3) in sheep. AGPAT4 is part of the 1-acyl-sn-glycerol 3-phosphate acyltransferases family which is involved in triglyceride synthesis (Lu et al., 2005; Takeuchi and Reue, 2009), with a role in the milk composition of ruminants (Bionaz and Loor, 2008; He et al., 2011; Jie Yang et al., 2016). It has an organ-related level of expression in human (Lu et al., 2005; Prasad et al., 2011) and plays a role in muscle development and meat quality in cattle (X. Fang et al., 2017). In addition, the variation of AGPAT genes expression in the mammary gland has been related to milk fat depression in dairy cows (Peterson et al., 2003). On its side, SLIT3 takes place into the SLIT-ROBO interactions which are involved in the development of the nervous system (K. H. Wang et al., 1999; Andrews et al., 2007), cardiac development (Liu et al., 2003; Mommersteeg et al., 2015) and reproduction, especially menstrual cycle (Dickinson and Duncan, 2010) and fertility (Amundson et al., 2015). SLIT3 was also shown to be differentially methylated and upregulated in goat mammary glands during lactation compared to the dry period (Zhang et al., 2017). Such DMRs related to milk quality and lactation might result directly from heat stress. Indeed, it was shown that heifers submitted to a thermic shock during fetal development had smaller mammary gland alveoli, resulting in lower milk production (Skibiel et al., 2018; Dahl et al., 2019). The DMRs observed could also result from indirect effects of the climate through a variation of the diet resulting from a shift in plant resources. Peterson et al. (2003) showed that the concentrate forage diet influences the mRNA abundance of AGPAT genes in mammary glands in dairy cows, involving a modification of milk fat. The regulation of the SLIT-ROBO pathway is also affected by the diet (Amundson et al., 2015). Indeed, the diet can directly influence the levels of DNA methylation (Anderson et al., 2012) either due to direct exposure of the individuals or to parental effects. For example, a heat-stress-induced perturbation of mammary glands development in the fetus, related to differential methylation and differential gene expression (Skibiel et al., 2018), results in effects persisting over the lactations and in the next generation too (Ouellet et al., 2020). Maternal nutrition also affects DNA methylation and gene expression of lambs, especially in muscle and adipose tissues (Lan et al., 2013; Peñagaricano et al., 2014; Namous et al., 2018). Obviously, we found differential methylation in genes related to lactation in ear biopsies, not mammary glands. However, sweat glands localized in the dermis and mammary glands have a common tissue origin, and recent studies showed that the most similar tissues had comparable CpG methylation patterns (Lokk et al., 2014). Given the epigenetic memory, embryonic patterns may be retained in adult tissues (Hon et al., 2013). In this case, an early stress (e.g., in utero) could modify the methylation pattern of embryonic cells, which could be retained in structures with a common origin. Further characterization of the methylation patterns in different organs, including mammary glands, is needed to confirm the possible role of these mechanisms.
We found that candidate genes for local adaptation in response to different temperature regimes were related to a large variety of biological pathways which were mainly involved in functions such as immunity, perception of the environment and reproduction, or related to productivity traits. While these functions were affected in both sheep and goats, the best candidate genes (i.e., top 0.001%) identified were different, confirming that alternative genes might be selected under the same environmental constraints for close species (e.g., Benjelloun, 2015; Raeymaekers et al., 2017) or even populations of the same species (Manceau et al., 2010; Benjelloun et al., 2015). However, out of the 0.025% top genes (i.e., 931 genes selected for the enrichment analyses), 41 were orthologous between sheep and goats and mainly related to the nervous system development. Interestingly, the present study also investigated epigenetic variations in the individuals displaying candidate genes for local adaptation. Our deliberately stringent approach to avoid false positive, which identified only one candidate DMR per species, prevents a global view but already points out similarities between traits impacted by genetic and epigenetic mechanisms. In sheep, genes impacted by selection such as BMPER, GCSAML, KSR2, PTPRE affect embryo quality, litter size and lactation, which are also influenced by SLIT3 that was differentially methylated. Likewise, in goats, AGPAT4 that was differentially methylated has an effect on lipid metabolism, meat quality and the presence of lipids in milk, these traits being influenced by several candidate genes such as CTSB, MYADM, SERPINB7, USH1C, GALNTL6 and CALHM3.
Conclusion
This study showed differentially methylated genes and candidate genes under selection related to a differential regime of temperatures underwent by sheep and goats. Genetic and epigenetic mechanisms share several target functions including lactation and milk composition. These adjustments of phenotypes to the two contrasted environments may result from drivers related to temperature either directly or indirectly through e.g., the variation of plant resources conditioning the diet. We do not yet have enough evidence to say that these results reflect the joint effects of selection and epigenetic regulation on the same traits. However, such combined effects might play an important role in evolution (Richards, 2011), and we can wonder whether DMRs might cushion phenotypic variations due to prior selection of genes or reflect the persistence of adaptive plasticity as long as adaptation is not complete? Deciphering the relative impact of climate and e.g., differential diets on the genetic and epigenetically-induced variations of phenotypes will be a major stake for developing future breeding strategies in the context of global change.
Statements
Data availability statement
The genetic dataset (fastq of Whole Genome Sequencing) analyzed in this study can be found in the European Nucleotide Archive (www.ebi.ac.uk/ena, accession code PRJEB7436). The epigenetic dataset (fastq of the MeDIP-seq) generated for this study can be found in the European Nucleotide Archive (www.ebi.ac.uk/ena, accession code PRJEB47531).
Ethics statement
The animal study was reviewed and approved by European Union Directive 86/609/EEC.
Author contributions
FP and FB supervised the project. FP and BB contributed to the design project. BB conducted the sampling of Moroccan small ruminants. FA, FP, and LD performed DNA extraction and managed the DNA samples. LD, CG, MK, PD, and FB performed the analyzes. LD and FP wrote the paper.
Funding
We thank the Ministère de l’Enseignement Supérieur, de la Recherche et de l’Innovation for financing LD. The LECA is part of the Labex OSUG@2020 (ANR10 LABX56). This study was partly funded by the projects Climgen (Program FACE ERA-NET PLUS 2014, ANR- 14-JAFC-0002) and MethylGoat (AGIR-2014, Univ. Grenoble Alpes). Computations were performed using the CIMENT infrastructure, which is supported by the Auvergne-Rhône-Alpes region (CPER07_13 CIRA).
Acknowledgments
We thank Florence Plisson and Sylvie Veyrenc for their help during the DNA manipulations (extraction, MeDIP-Seq); and Hervé Acloque for help in setting up the MeDIP-Seq approach. This work was performed in collaboration with the GeT core facility, Toulouse, France (http://get.genotoul.fr), and was supported by France Génomique National infrastructure, funded as part of “Investissement d’avenir” program managed by Agence Nationale pour la Recherche (contract ANR-10-INBS-09). All the computations presented in this paper were performed using the GRICAD infrastructure (https://gricad.univ-grenoble-alpes.fr), which is partly supported by the Equip@Meso project (reference ANR-10-EQPX-29-01) of the programme Investissements d’Avenir supervised by the Agence Nationale pour la Recherche.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2021.745284/full#supplementary-material
References
1
Abo-IsmailM. K.BritoL. F.MillerS. P.SargolzaeiM.GrossiD. A.MooreS. S.et al (2017). Genome-wide Association Studies and Genomic Prediction of Breeding Values for Calving Performance and Body Conformation Traits in Holstein Cattle. Genet. Sel Evol.49 (1), 82. 10.1186/s12711-017-0356-8
2
AggreyS. E.HabashyW. S.MilfortM. C.AdomakoK.FullerA. L.RekayaR. (2018). “Molecular and Cellular Mechanisms that Underlie Genes and Antioxidant Enzyme Activities in Meat-type Birds during Heat Stress,” in World Congress on Genetics Applied to Livestock Production, Auckland, New Zealand.
3
Al-TamimiH. J. (2007). Thermoregulatory Response of Goat Kids Subjected to Heat Stress. Small Ruminant Res.71 (1), 280–285. 10.1016/j.smallrumres.2006.04.013
4
AlbertoF. J.BoyerF.Orozco-terWengelP.StreeterI.ServinB.de VillemereuilP.et al (2018). Convergent Genomic Signatures of Domestication in Sheep and Goats. Nat. Commun.9 (1), 813. 10.1038/s41467-018-03206-y
5
AmundsonO. L.FountainT. H.LarimoreE. L.RichardsonB. N.McNeelA. K.WrightE. C.et al (2015). Postweaning Nutritional Programming of Ovarian Development in Beef Heifers1,2. J. Anim. Sci.93 (11), 5232–5239. 10.2527/jas.2015-9067
6
AndersonO. S.SantK. E.DolinoyD. C. (2012). Nutrition and Epigenetics: An Interplay of Dietary Methyl Donors, One-Carbon Metabolism and DNA Methylation. J. Nutr. Biochem.23 (8), 853–859. 10.1016/j.jnutbio.2012.03.003
7
AndrewsW. D.BarberM.ParnavelasJ. G. (2007). Slit?Robo Interactions during Cortical Development. J. Anat.211 (2), 188–198. 10.1111/j.1469-7580.2007.00750.x
8
BenjellounB.AlbertoF. J.StreeterI.BoyerF. d. r.CoissacE.StuckiS.et al (2015). Characterizing Neutral Genomic Diversity and Selection Signatures in Indigenous Populations of Moroccan Goats (Capra hircus) Using WGS Data. Front. Genet.6, 107. 10.3389/fgene.2015.00107
9
BenjellounB. (2015). “Diversité des génomes et adaptation locale des petits ruminants d’un pays méditerranéen: Le Maroc,”. doctoral thesis (Grenoble (France): Université Grenoble Alpes). http://www.theses.fr/2015GREAV011.
10
BettB.KiungaP.GachohiJ.SindatoC.MbothaD.RobinsonT.et al (2017). Effects of Climate Change on the Occurrence and Distribution of Livestock Diseases. Prev. Vet. Med.137, 119–129. 10.1016/j.prevetmed.2016.11.019
11
BionazM.LoorJ. J. (2008). ACSL1, AGPAT6, FABP3, LPIN1, and SLC27A6 Are the Most Abundant Isoforms in Bovine Mammary Tissue and Their Expression Is Affected by Stage of Lactation. J. Nutr.138 (6), 1019–1024. 10.1093/jn/138.6.1019
12
BolormaaS.BrownD. J.SwanA. A.van der WerfJ. H. J.HayesB. J.DaetwylerH. D. (2017). Genomic Prediction of Reproduction Traits for Merino Sheep. Anim. Genet.48 (3), 338–348. 10.1111/age.12541
13
BornK.ChristophM.FinkA. H.KnippertzP.PaethH.SpethP. (2008). “Moroccan Climate in the Present and Future: Combined View from Observational Data and Regional Climate Scenarios,” in Climatic Changes and Water Resources in the Middle East and North Africa. Editors ZereiniF.HötzlH. (Springer Berlin Heidelberg), 29–45. 10.1007/978-3-540-85047-2_4
14
BossdorfO.RichardsC. L.PigliucciM. (2008). Epigenetics for Ecologists. Ecol. Lett.11, 106–115. 10.1111/j.1461-0248.2007.01130.x
15
BrufordM. W.GinjaC.HoffmannI.JoostS.Orozco-terWengelP.AlbertoF. J.et al (2015). Prospects and Challenges for the Conservation of Farm Animal Genomic Resources, 2015-2025. Front. Genet.6, 314. 10.3389/fgene.2015.00314
16
CarbonS.MungallC. (2018). Gene Ontology Data Archive [Data Set]. Zenodo. 10.5281/zenodo.3980761
17
CaropreseM.AlbenzioM.BrunoA.FedeleV.SantilloA.SeviA. (2011). Effect of Solar Radiation and Flaxseed Supplementation on Milk Production and Fatty Acid Profile of Lactating Ewes under High Ambient Temperature. J. Dairy Sci.94 (8), 3856–3867. 10.3168/jds.2010-4067
18
CarrD.Lewin-KohN.MaechlerM.SarkarD. (2019). Hexbin: Hexagonal Binning Routines. (R package version 1.28.0) [Computer software]. https://CRAN.R-project.org/package=hexbin.
19
ChenI.-C.HillJ. K.OhlemüllerR.RoyD. B.ThomasC. D. (2011). Rapid Range Shifts of Species Associated with High Levels of Climate Warming. Science333 (6045), 1024–1026. 10.1126/science.1206432
20
CooneyC. A.DaveA. A.WolffG. L. (2002). Maternal Methyl Supplements in Mice Affect Epigenetic Variation and DNA Methylation of Offspring. J. Nutr.132 (8), 2393S–2400S. 10.1093/jn/132.8.2393S
21
Dado-SennB.SkibielA. L.FabrisT. F.DahlG. E.LaportaJ. (2019). Dry Period Heat Stress Induces Microstructural Changes in the Lactating Mammary Gland. PLOS ONE14 (9), e0222120. 10.1371/journal.pone.0222120
22
DahlG. E.SkibielA. L.LaportaJ. (2019). In Utero Heat Stress Programs Reduced Performance and Health in Calves. Vet. Clin. North America: Food Anim. Pract.35 (2), 343–353. 10.1016/j.cvfa.2019.02.005
23
DanchinÉ.CharmantierA.ChampagneF. A.MesoudiA.PujolB.BlanchetS. (2011). Beyond DNA: Integrating Inclusive Inheritance into an Extended Theory of Evolution. Nat. Rev. Genet.12 (7), 475–486. 10.1038/nrg3028
24
DanecekP.AutonA.AbecasisG.AlbersC. A.BanksE.DePristoM. A.et al (2011). The Variant Call Format and VCFtools. Bioinformatics27 (15), 2156–2158. 10.1093/bioinformatics/btr330
25
DickinsonR. E.DuncanW. C. (2010). The SLIT-ROBO Pathway: a Regulator of Cell Function with Implications for the Reproductive System. Reproduction139 (4), 697–704. 10.1530/REP-10-0017
26
DohertyR.FarrellyC. O.MeadeK. G. (2014). Comparative Epigenetics: Relevance to the Regulation of Production and Health Traits in Cattle. Anim. Genet.45, 3–14. 10.1111/age.12140
27
DongY.XieM.JiangY.XiaoN.DuX.ZhangW.et al (2013). Sequencing and Automated Whole-Genome Optical Mapping of the Genome of a Domestic Goat (Capra hircus). Nat. Biotechnol.31 (2), 135–141. 10.1038/nbt.2478
28
DongY.ZhangX.XieM.ArefnezhadB.WangZ.WangW.et al (2015). Reference Genome of Wild Goat (capra Aegagrus) and Sequencing of Goat Breeds Provide Insight into Genic Basis of Goat Domestication. BMC Genomics16 (431). 10.1186/s12864-015-1606-1
29
DoranA. G.BerryD. P.CreeveyC. J. (2014). Whole Genome Association Study Identifies Regions of the Bovine Genome and Biological Pathways Involved in Carcass Trait Performance in Holstein-Friesian Cattle. BMC Genomics15 (1), 837. 10.1186/1471-2164-15-837
30
EichtenS. R.BriskineR.SongJ.LiQ.Swanson-WagnerR.HermansonP. J.et al (2013). Epigenetic and Genetic Influences on DNA Methylation Variation in Maize Populations. Plant Cell25 (8), 2783–2797. 10.1105/tpc.113.114793
31
FangX.ZhaoZ.YuH.LiG.JiangP.YangY.et al (2017). Comparative Genome-wide Methylation Analysis of Longissimus Dorsi Muscles between Japanese Black (Wagyu) and Chinese Red Steppes Cattle. PLoS ONE12 (8), e0182492. 10.1371/journal.pone.0182492
32
FangZ.-H.PauschH. (2019). Multi-trait Meta-Analyses Reveal 25 Quantitative Trait Loci for Economically Important Traits in Brown Swiss Cattle. BMC Genomics20 (1), 695. 10.1186/s12864-019-6066-6
33
FengS.CokusS. J.ZhangX.ChenP.-Y.BostickM.GollM. G.et al (2010). Conservation and Divergence of Methylation Patterning in Plants and Animals. Proc. Natl. Acad. Sci.107 (19), 8689–8694. 10.1073/pnas.1002720107
34
FoustC. M.PreiteV.SchreyA. W.AlvarezM.RobertsonM. H.VerhoevenK. J. F.et al (2016). Genetic and Epigenetic Differences Associated with Environmental Gradients in Replicate Populations of Two Salt Marsh Perennials. Mol. Ecol.25 (8), 1639–1652. 10.1111/mec.13522
35
GibneyE. R.NolanC. M. (2010). Epigenetics and Gene Expression. Heredity105 (1), 4–13. 10.1038/hdy.2010.54
36
GiriA. K.MidhaS.BanerjeeP.AgrawalA.MehdiS. J.DhingraR.et al (2016). Common Variants in CLDN2 and MORC4 Genes Confer Disease Susceptibility in Patients with Chronic Pancreatitis. PLOS ONE11 (1), e0147345. 10.1371/journal.pone.0147345
37
GonzalezM.VillaR.VillaC.GonzalezV.MontanoM.MedinaG.et al (2020). Inspection of Real and Imputed Genotypes Reveled 76 SNPs Associated to Rear Udder Height in Holstein Cattle. J. Adv. Vet. Anim. Res.7 (2), 234. 10.5455/javar.2020.g415
38
GouX.WangZ.LiN.QiuF.XuZ.YanD.et al (2014). Whole-genome Sequencing of Six Dog Breeds from Continuous Altitudes Reveals Adaptation to High-Altitude Hypoxia. Genome Res.24, 1308–1315. 10.1101/gr.171876.113
39
GuoL.Costanzo-GarveyD. L.SmithD. R.NeilsenB. K.MacDonaldR. G.LewisR. E. (2017). Kinase Suppressor of Ras 2 (KSR2) Expression in the Brain Regulates Energy Balance and Glucose Homeostasis. Mol. Metab.6 (2), 194–205. 10.1016/j.molmet.2016.12.004
40
HeC.WangC.ChangZ. H.GuoB. L.LiR.YueX. P.et al (2011). AGPAT6 Polymorphism and its Association with Milk Traits of Dairy Goats. Genet. Mol. Res.10 (4), 2747–2756. 10.4238/2011.November.4.8
41
HenriquesD.WallbergA.Chávez-GalarzaJ.JohnstonJ. S.WebsterM. T.PintoM. A. (2018). Whole Genome SNP-Associated Signatures of Local Adaptation in Honeybees of the Iberian Peninsula. Sci. Rep.8, 11145. 10.1038/s41598-018-29469-5
42
HerreraC. M.BazagaP. (2010). Epigenetic Differentiation and Relationship to Adaptive Genetic Divergence in Discrete Populations of the Violet Viola Cazorlensis. New Phytol.187 (3), 867–876. 10.1111/j.1469-8137.2010.03298.x
43
HijmansR. J.CameronS. E.ParraJ. L.JonesP. G.JarvisA. (2005). Very High Resolution Interpolated Climate Surfaces for Global Land Areas. Int. J. Climatol.25 (15), 1965–1978. 10.1002/joc.1276
44
HonG. C.RajagopalN.ShenY.McClearyD. F.YueF.DangM. D.et al (2013). Epigenetic Memory at Embryonic Enhancers Identified in DNA Methylation Maps from Adult Mouse Tissues. Nat. Genet.45, 1198–1206. 10.1038/ng.2746
45
HusseinH. A.StaufenbielR. (2012). Variations in Copper Concentration and Ceruloplasmin Activity of Dairy Cows in Relation to Lactation Stages with Regard to Ceruloplasmin to Copper Ratios. Biol. Trace Elem. Res.146 (1), 47–52. 10.1007/s12011-011-9226-3
46
JiangY.XieM.ChenW.TalbotR.MaddoxJ. F.FarautT.et al (2014). The Sheep Genome Illuminates Biology of the Rumen and Lipid Metabolism. Science344 (6188), 1168–1173. 10.1126/science.1252806
47
KillickR.EckleyI. A. (2014). Changepoint: AnRPackage for Changepoint Analysis. J. Stat. Soft.58 (3), 1–19. 10.18637/jss.v058.i03
48
KimG.ClarkeC. R.LaroseH.TranH. T.HaakD. C.ZhangL.et al (2017). Herbicide Injury Induces DNA Methylome Alterations in Arabidopsis. PeerJ5, e3560. 10.7717/peerj.3560
49
KrishnanG.BagathM.PragnaP.VidyaM. K.AleenaJ.ArchanaP. R.et al (2017). “Mitigation of the Heat Stress Impact in Livestock Reproduction,” In Payan Carreira, R. (IntechOpen).TheriogenologyAvailable at: https://www.intechopen.com/chapters/55491.
50
KwonS. G.HwangJ. H.ParkD. H.KimT. W.KangD. G.KangK. H.et al (2016). Identification of Differentially Expressed Genes Associated with Litter Size in Berkshire Pig Placenta. PLOS ONE11 (4), e0153311. 10.1371/journal.pone.0153311
51
LambethJ. D. (2004). NOX Enzymes and the Biology of Reactive Oxygen. Nat. Rev. Immunol.4 (3), 181–189. 10.1038/nri1312
52
LanX.CretneyE. C.KroppJ.KhateebK.BergM. A.PeñagaricanoF.et al (2013). Maternal Diet during Pregnancy Induces Gene Expression and DNA Methylation Changes in Fetal Tissues in Sheep. Front. Genet.4, 49. 10.3389/fgene.2013.00049
53
LaodimT.ElzoM. A.KoonawootrittrironS.SuwanasopeeT.JattawaD. (2017). Identification of SNP Markers Associated with Milk and Fat Yields in Multibreed Dairy Cattle Using Two Genetic Group Structures. Livestock Sci.206, 95–104. 10.1016/j.livsci.2017.10.015
54
LeaA. J.AltmannJ.AlbertsS. C.TungJ. (2016). Resource Base Influences Genome‐wide DNA Methylation Levels in Wild Baboons (Papio cynocephalus). Mol. Ecol.25 (8), 1681–1696. 10.1111/mec.13436
55
LeeB.-Y.LeeK.-N.LeeT.ParkJ.-H.KimS.-M.LeeH.-S.et al (2015). Bovine Genome-wide Association Study for Genetic Elements to Resist the Infection of Foot-And-Mouth Disease in the Field. Asian Australas. J. Anim. Sci.28 (2), 166–170. 10.5713/ajas.14.0383
56
LiH.DurbinR. (2009). Fast and Accurate Short Read Alignment with Burrows-Wheeler Transform. Bioinformatics25 (14), 1754–1760. 10.1093/bioinformatics/btp324
57
LienhardM.GrimmC.MorkelM.HerwigR.ChavezL. (2014). MEDIPS: Genome-wide Differential Coverage Analysis of Sequencing Data Derived from DNA Enrichment Experiments. Bioinformatics30 (2), 284–286. 10.1093/bioinformatics/btt650
58
LiuJ.ZhangL.WangD.ShenH.JiangM.MeiP.et al (2003). Congenital Diaphragmatic Hernia, Kidney Agenesis and Cardiac Defects Associated with Slit3-Deficiency in Mice. Mech. Dev.120 (9), 1059–1070. 10.1016/S0925-4773(03)00161-8
59
LokkK.ModhukurV.RajashekarB.MärtensK.MägiR.KoldeR.et al (2014). DNA Methylome Profiling of Human Tissues Identifies Global and Tissue-specific Methylation Patterns. Genome Biol.15 (4), 3248. 10.1186/gb-2014-15-4-r54
60
LuB.JiangY. J.ManM. Q.BrownB.EliasP. M.FeingoldK. R. (2005). Expression and Regulation of 1-Acyl-Sn-Glycerol- 3-phosphate Acyltransferases in the Epidermis. J. Lipid Res.46 (11), 2448–2457. 10.1194/jlr.M500258-JLR200
61
LukićB.FerenčakovićM.ŠalamonD.ČačićM.OrehovačkiV.IacolinaL.et al (2020). Conservation Genomic Analysis of the Croatian Indigenous Black Slavonian and Turopolje Pig Breeds. Front. Genet.11, 261. 10.3389/fgene.2020.00261
62
ManceauM.DominguesV. S.LinnenC. R.RosenblumE. B.HoekstraH. E. (2010). Convergence in Pigmentation at Multiple Levels: Mutations, Genes and Function. Phil. Trans. R. Soc. B365 (1552), 2439–2450. 10.1098/rstb.2010.0104
63
MaraiI. F. M.El-DarawanyA. A.FadielA.Abdel-HafezM. A. M. (2007). Physiological Traits as Affected by Heat Stress in Sheep-A Review. Small Ruminant Res.71 (1–3), 1–12. 10.1016/j.smallrumres.2006.10.003
64
MartinG. B.RodgerJ.BlacheD. (2004). Nutritional and Environmental Effects on Reproduction in Small Ruminants. Reprod. Fertil. Dev.16 (4), 491. 10.1071/RD04035
65
McDonaldK. W.McClureC. J. W.RolekB. W.HillG. E. (2012). Diversity of Birds in Eastern North America Shifts north with Global Warming. Ecol. Evol.2 (12), 3052–3060. 10.1002/ece3.410
66
MeyerhoefferD. C.WettemannR. P.ColemanS. W.WellsM. E. (1985). Reproductive Criteria of Beef Bulls during and after Exposure to Increased Ambient Temperature. J. Anim. Sci.60 (2), 352–357. 10.2527/jas1985.602352x
67
MiH.MuruganujanA.EbertD.HuangX.ThomasP. D. (2019). PANTHER Version 14: More Genomes, a New PANTHER GO-Slim and Improvements in Enrichment Analysis Tools. Nucleic Acids Res.47 (D1), D419–D426. 10.1093/nar/gky1038
68
MohammadiA.AlijaniS.RafatS. A.Abdollahi-ArpanahiR. (2020). Genome-Wide Association Study and Pathway Analysis for Female Fertility Traits in Iranian Holstein Cattle. Ann. Anim. Sci.20 (3), 825–851. 10.2478/aoas-2020-0031
69
MommersteegM. T. M.YehM. L.ParnavelasJ. G.AndrewsW. D. (2015). Disrupted Slit-Robo Signalling Results in Membranous Ventricular Septum Defects and Bicuspid Aortic Valves. Cardiovasc. Res.106 (1), 55–66. 10.1093/cvr/cvv040
70
MonteiroA. P. A.TaoS.ThompsonI. M. T.DahlG. E. (2016). In Utero heat Stress Decreases Calf Survival and Performance through the First Lactation. J. Dairy Sci.99 (10), 8443–8450. 10.3168/jds.2016-11072
71
MotaL. F. M.LopesF. B.Fernandes JúniorG. A.RosaG. J. M.MagalhãesA. F. B.CarvalheiroR.et al (2020). Genome-wide Scan Highlights the Role of Candidate Genes on Phenotypic Plasticity for Age at First Calving in Nellore Heifers. Sci. Rep.10 (1), 6481. 10.1038/s41598-020-63516-4
72
NamousH.PeñagaricanoF.Del CorvoM.CapraE.ThomasD. L.StellaA.et al (2018). Integrative Analysis of Methylomic and Transcriptomic Data in Fetal Sheep Muscle Tissues in Response to Maternal Diet during Pregnancy. BMC Genomics19 (1), 123. 10.1186/s12864-018-4509-0
73
OuelletV.LaportaJ.DahlG. E. (2020). Late Gestation Heat Stress in Dairy Cows: Effects on Dam and Daughter. Theriogenology150, 471–479. 10.1016/j.theriogenology.2020.03.011
74
PeñagaricanoF.WangX.RosaG. J.RadunzA. E.KhatibH. (2014). Maternal Nutrition Induces Gene Expression Changes in Fetal Muscle and Adipose Tissues in Sheep. BMC Genomics15 (1), 1034. 10.1186/1471-2164-15-1034
75
PetersonD. G.MatitashviliE. A.BaumanD. E. (2003). Diet-Induced Milk Fat Depression in Dairy Cows Results in Increased Trans-10, Cis-12 CLA in Milk Fat and Coordinate Suppression of mRNA Abundance for Mammary Enzymes Involved in Milk Fat Synthesis. J. Nutr.133 (10), 3098–3102. 10.1093/jn/133.10.3098
76
PierceC. F.SpeidelS. E.ColemanS. J.EnnsR. M.BaileyD. W.MedranoJ. F.et al (2020). Genome-wide Association Studies of Beef Cow Terrain-Use Traits Using Bayesian Multiple-SNP Regression. Livestock Sci.232, 103900. 10.1016/j.livsci.2019.103900
77
PimentelE. C. G.BauersachsS.TietzeM.SimianerH.TetensJ.ThallerG.et al (2011). Exploration of Relationships between Production and Fertility Traits in Dairy Cattle via Association Studies of SNPs within Candidate Genes Derived by Expression Profiling. Anim. Genet.42 (3), 251–262. 10.1111/j.1365-2052.2010.02148.x
78
PonsuksiliS.MuraniE.SchwerinM.SchellanderK.TesfayeD.WimmersK. (2012). Gene Expression and DNA-Methylation of Bovine Pretransfer Endometrium Depending on its Receptivity after In Vitro-Produced Embryo Transfer. PLOS ONE7 (8), e42402. 10.1371/journal.pone.0042402
79
PrasadS. S.GargA.AgarwalA. K. (2011). Enzymatic Activities of the Human AGPAT Isoform 3 and Isoform 5: Localization of AGPAT5 to Mitochondria. J. Lipid Res.52 (3), 451–462. 10.1194/jlr.M007575
80
Puig-OliverasA.BallesterM.CorominasJ.RevillaM.EstelléJ.FernándezA. I.et al (2014). A Co-association Network Analysis of the Genetic Determination of Pig Conformation, Growth and Fatness. PLOS ONE9 (12), e114862. 10.1371/journal.pone.0114862
81
QuinlanA. R.HallI. M. (2010). BEDTools: A Flexible Suite of Utilities for Comparing Genomic Features. Bioinformatics26 (6), 841–842. 10.1093/bioinformatics/btq033
82
R Core Team (2019). R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing. Available at: https://www.R-project.org/.
83
RaeymaekersJ. A. M.ChaturvediA.HablützelP. I.VerdonckI.HellemansB.MaesG. E.et al (2017). Adaptive and Non-adaptive Divergence in a Common Landscape. Nat. Commun.8 (1), 267. 10.1038/s41467-017-00256-6
84
RazzaghiA.ValizadehR.NaserianA. A.MesgaranM. D.CarpenterA. J.GhaffariM. H. (2016). Effect of Dietary Sugar Concentration and sunflower Seed Supplementation on Lactation Performance, Ruminal Fermentation, Milk Fatty Acid Profile, and Blood Metabolites of Dairy Cows. J. Dairy Sci.99 (5), 3539–3548. 10.3168/jds.2015-10565
85
RichardsE. J. (2011). Natural Epigenetic Variation in Plant Species: a View from the Field. Curr. Opin. Plant Biol.14 (2), 204–209. 10.1016/j.pbi.2011.03.009
86
RobinsonJ. J.AshworthC. J.RookeJ. A.MitchellL. M.McEvoyT. G. (2006). Nutrition and Fertility in Ruminant Livestock. Anim. Feed Sci. Tech.126 (3), 259–276. 10.1016/j.anifeedsci.2005.08.006
87
RobinsonM. D.McCarthyD. J.SmythG. K. (2010). edgeR: A Bioconductor Package for Differential Expression Analysis of Digital Gene Expression Data. Bioinformatics26 (1), 139–140. 10.1093/bioinformatics/btp616
88
RussoV.FontanesiL.DavoliR.Nanni CostaL.CagnazzoM.ButtazzoniL.et al (2002). Investigation of Candidate Genes for Meat Quality in Dry-Cured Ham Production: the Porcine Cathepsin B (CTSB ) and Cystatin B (CSTB ) Genes. Anim. Genet.33 (2), 123–131. 10.1046/j.1365-2052.2002.00835.x
89
SallehM. S.MazzoniG.HöglundJ. K.OlijhoekD. W.LundP.LøvendahlP.et al (2017). RNA-seq Transcriptomics and Pathway Analyses Reveal Potential Regulatory Genes and Molecular Mechanisms in High- and Low-Residual Feed Intake in Nordic Dairy Cattle. BMC Genomics18 (1), 258. 10.1186/s12864-017-3622-9
90
SavolainenO.LascouxM.MeriläJ. (2013). Ecological Genomics of Local Adaptation. Nat. Rev. Genet.14 (11), 807–820. 10.1038/nrg3522
91
SejianV.MauryaV. P.NaqviS. M. K. (2010). Adaptability and Growth of Malpura Ewes Subjected to thermal and Nutritional Stress. Trop. Anim. Health Prod.42 (8), 1763–1770. 10.1007/s11250-010-9633-z
92
SevaneN.MartínezR.BrufordM. W. (2018). Genome-wide Differential DNA Methylation in Tropically Adapted Creole Cattle and Their Iberian Ancestors. Anim. Genet.50 (1), 15–26. 10.1111/age.12731
93
ShinD.LeeC.ParkK.-D.KimH.ChoK.-h. (2017). Genome-association Analysis of Korean Holstein Milk Traits Using Genomic Estimated Breeding Value. Asian-australas J. Anim. Sci.30 (3), 309–319. 10.5713/ajas.15.0608
94
SkibielA. L.Dado-SennB.FabrisT. F.DahlG. E.LaportaJ. (2018). In Utero exposure to thermal Stress Has Long-Term Effects on Mammary Gland Microstructure and Function in Dairy Cattle. PLOS ONE13 (10), e0206046. 10.1371/journal.pone.0206046
95
SkibielA. L.PeñagaricanoF.AmorínR.AhmedB. M.DahlG. E.LaportaJ. (2018). In Utero Heat Stress Alters the Offspring Epigenome. Sci. Rep.8 (1), 14609. 10.1038/s41598-018-32975-1
96
StoreyJ. D.BassA. J.DabneyA.RobinsonD. (2015). Qvalue: Q-Value Estimation for False Discovery Rate Control. R package version 2.8.0. . Http://Github.Com/Jdstorey/Qvalue.
97
SugimotoM.FujikawaA.WomackJ. E.SugimotoY. (2006). Evidence that Bovine Forebrain Embryonic Zinc finger-like Gene Influences Immune Response Associated with Mastitis Resistance. Proc. Natl. Acad. Sci.103 (17), 6454–6459. 10.1073/pnas.0601015103
98
SuzukiM. M.BirdA. (2008). DNA Methylation Landscapes: Provocative Insights from Epigenomics. Nat. Rev. Genet.9 (6), 465–476. 10.1038/nrg2341
99
TaberletP.ValentiniA.RezaeiH. R.NaderiS.PompanonF.NegriniR.et al (2008). Are Cattle, Sheep, and Goats Endangered Species?Mol. Ecol.17 (1), 275–284. 10.1111/j.1365-294X.2007.03475.x
100
TakeuchiK.ReueK. (2009). Biochemistry, Physiology, and Genetics of GPAT, AGPAT, and Lipin Enzymes in Triglyceride Synthesis. Am. J. Physiology-Endocrinology Metab.296 (6), E1195–E1209. 10.1152/ajpendo.90958.2008
101
TaoS.BubolzJ. W.do AmaralB. C.ThompsonI. M.HayenM. J.JohnsonS. E.et al (2011). Effect of Heat Stress during the Dry Period on Mammary Gland Development. J. Dairy Sci.94 (12), 5976–5986. 10.3168/jds.2011-4329
102
ValléeM.BeaudryD.RobergeC.MatteJ. J.BlouinR.PalinM.-F. (2003). Isolation of Differentially Expressed Genes in Conceptuses and Endometrial Tissue of Sows in Early Gestation1. Biol. Reprod.69 (5), 1697–1706. 10.1095/biolreprod.103.019307
103
WaggottD.HaiderS.BoutrosP. C. (2017). Bedr: Genomic Region Processing Using Tools Such as “BEDTools”, “BEDOPS” and “Tabix”. R package version 1.0.4. Https://CRAN.R-Project.Org/Package=bedr.
104
WangK. H.BroseK.ArnottD.KiddT.GoodmanC. S.HenzelW.et al (1999). Biochemical Purification of a Mammalian Slit Protein as a Positive Regulator of Sensory Axon Elongation and Branching. Cell96 (6), 771–784. 10.1016/S0092-8674(00)80588-7
105
WangX.LanX.RadunzA. E.KhatibH. (2015). Maternal Nutrition during Pregnancy Is Associated with Differential Expression of Imprinted Genes and DNA Methyltranfereases in Muscle of Beef Cattle Offspring1. J. Anim. Sci.93 (1), 35–40. 10.2527/jas.2014-8148
106
WeirB. S.CockerhamC. C. (1984). Estimating F -Statistics for the Analysis of Population Structure. Evolution38 (6), 1358–1370. 10.2307/240864110.1111/j.1558-5646.1984.tb05657.x
107
XiangR.MacLeodI. M.BolormaaS.GoddardM. E. (2017). Genome-wide Comparative Analyses of Correlated and Uncorrelated Phenotypes Identify Major Pleiotropic Variants in Dairy Cattle. Sci. Rep.7 (1), 9248. 10.1038/s41598-017-09788-9
108
YangJ.JiangJ.LiuX.WangH.GuoG.ZhangQ.et al (2016). Differential Expression of Genes in Milk of Dairy Cattle during Lactation. Anim. Genet.47 (2), 174–180. 10.1111/age.12394
109
YangJ.LiW.-R.LvF.-H.HeS.-G.TianS.-L.PengW.-F.et al (2016). Whole-Genome Sequencing of Native Sheep Provides Insights into Rapid Adaptations to Extreme Environments. Mol. Biol. Evol.33 (10), 2576–2592. 10.1093/molbev/msw129
110
ZederM. A. (2008). Domestication and Early Agriculture in the Mediterranean Basin: Origins, Diffusion, and Impact. Proc. Natl. Acad. Sci.105 (33), 11597–11604. 10.1073/pnas.0801317105
111
ZhangX.ZhangS.MaL.JiangE.XuH.ChenR.et al (2017). Reduced Representation Bisulfite Sequencing (RRBS) of Dairy Goat Mammary Glands Reveals DNA Methylation Profiles of Integrated Genome-wide and Critical Milk-Related Genes. Oncotarget8 (70), 115326–115344. 10.18632/oncotarget.23260
112
ZhaoC.GuiL.LiY.PlathM.ZanL. (2015). Associations between Allelic Polymorphism of the BMP Binding Endothelial Regulator and Phenotypic Variation of Cattle. Mol. Cell Probes29 (6), 358–364. 10.1016/j.mcp.2015.09.007
113
ZhouD. a.WeiZ.DengS.WangT.ZaiM.WangH.et al (2013). SASH1 Regulates Melanocyte Transepithelial Migration through a Novel Gαs-SASH1-IQGAP1-E-Cadherin Dependent Pathway. Cell Signal.25 (6), 1526–1538. 10.1016/j.cellsig.2012.12.025
114
ZhouX.LindsayH.RobinsonM. D. (2014). Robustly Detecting Differential Expression in RNA Sequencing Data Using Observation Weights. Nucleic Acids Res.42 (11), e91. 10.1093/nar/gku310
115
ZoldošV.BirušI.MuratovićE.ŠatovićZ.VojtaA.RobinO.et al (2018). Epigenetic Differentiation of Natural Populations of Lilium Bosniacum Associated with Contrasting Habitat Conditions. Genome Biol. Evol.10 (1), 291–303. 10.1093/gbe/evy010
Summary
Keywords
adaptation, methylated DNA, whole genome, Morocco, small ruminant, acclimation
Citation
Denoyelle L, de Villemereuil P, Boyer F, Khelifi M, Gaffet C, Alberto F, Benjelloun B and Pompanon F (2021) Genetic Variations and Differential DNA Methylation to Face Contrasted Climates in Small Ruminants: An Analysis on Traditionally-Managed Sheep and Goats. Front. Genet. 12:745284. doi: 10.3389/fgene.2021.745284
Received
21 July 2021
Accepted
02 September 2021
Published
28 September 2021
Volume
12 - 2021
Edited by
Xingbo Zhao, China Agricultural University, China
Reviewed by
Yu Jiang, Northwest A and F University, China
Juliane Friedrich, University of Edinburgh, United Kingdom
Updates
Copyright
© 2021 Denoyelle, de Villemereuil, Boyer, Khelifi, Gaffet, Alberto, Benjelloun and Pompanon.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: François Pompanon, francois.pompanon@univ-grenoble-alpes.fr
This article was submitted to Evolutionary and Population Genetics, a section of the journal Frontiers in Genetics
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.