# Genetic determinism exists for the global DNA methylation rate in sheep

Dominique Hazard, Florence Plisson-Petit, Carole C. Moreno, Stéphane

Fabre, Laurence Drouilhet

► To cite this version:

Dominique Hazard, Florence Plisson-Petit, Carole C. Moreno, Stéphane Fabre, Laurence Drouilhet. Genetic determinism exists for the global DNA methylation rate in sheep. Frontiers in Genetics, 2020, 11, 10.3389/fgene.2020.616960. [https://doi.org/10.3389/fgene.2020.616960](https://doi.org/10.3389/fgene.2020.616960). [arXiv:2006.01752](https://arxiv.org/abs/2006.01752)

HAL Id: hal-03137552

https://hal.inrae.fr/hal-03137552

Submitted on 19 Feb 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
<figure><img src="image_2.png" /></figure>

# Genetic Determinism Exists for the Global DNA Methylation Rate in Sheep

Dominique Hazard*, Florence Plisson-Petit, Carole Moreno-Romieux, Stéphane Fabre and Laurence Drouilhet

GenPhySE, Université de Toulouse, INRAE, ENVT, Castanet Tolosan, France

Recent studies showed that epigenetic marks, including DNA methylation, influence production and adaptive traits in plants and animals. So far, most studies dealing with genetics and epigenetics considered DNA methylation sites independently. However, the genetic basis of the global DNA methylation rate (GDMR) remains unknown. The main objective of the present study was to investigate genetic determinism of GDMR in sheep. The experiment was conducted on 1,047 Romane sheep allocated into 10 half-sib families. After weaning, all the lambs were phenotyped for global GDMR in blood as well as for production and adaptive traits. GDMR was measured by LUminometric Methylation Analysis (LUMA) using a pyrosequencing approach. Association analyses were conducted on some of the lambs (n = 775) genotyped by using the Illumina OvineSNP50 BeadChip. Blood GDMR varied among the animals (average 70.7 ± 6.0%). Female lambs had significantly higher GDMR than male lambs. Inter-individual variability of blood GDMR had an additive genetic component and heritability was moderate (h² = 0.20 ± 0.05). No significant genetic correlation was found between GDMR and growth or carcass traits, birthcoat, or social behaviors. Association analyses revealed 28 QTLs associated with blood GDMR. Seven genomic regions on chromosomes 1, 5, 11, 17, 24, and 26 were of most interest due to either high significant associations with GDMR or to the relevance of genes located close to the QTLs. QTL effects were moderate. Genomic regions associated with GDMR harbored several genes not yet described as being involved in DNA methylation, but some are already known to play an active role in gene expression. In addition, some candidate genes, *CHD1*, *NCOA3*, *KDM8*, *KAT7*, and *KAT6A* have previously been described to be involved in epigenetic modifications. In conclusion, the results of the present study indicate that blood GDMR in domestic sheep is under polygenic influence and provide new insights into DNA methylation genetic determinism.

**Keywords:** epigenetics, blood, heritability, genetic correlation, QTL

## INTRODUCTION

Genetic selection in several livestock species has recently progressed to include heritable DNA polymorphisms (i.e., genomic selection) for the improvement of production traits (Goddard et al., 2010). However, recent studies also showed that non-genetic information responsible for phenotypic differences between animals can also be inherited across generations (Johannes et al., 2008; Danchin et al., 2011). One such non-genetic inherited piece of information

**OPEN ACCESS**

**Edited by:**

Eveline M. Ibeagha-Awemu,
Agriculture and Agri-Food Canada
(AAFC), Canada

**Reviewed by:**

Hans D. Daetwyler,
La Trobe University, Australia
Matthew L. Spangler,
University of Nebraska-Lincoln,
United States

**\*Correspondence:**

Dominique Hazard
dominique.hazard@inrae.fr

**Specialty section:**

This article was submitted to
Livestock Genomics,
a section of the journal
Frontiers in Genetics

**Received:** 13 October 2020

**Accepted:** 07 December 2020

**Published:** 23 December 2020

**Citation:**

Hazard D, Plisson-Petit F,
Moreno-Romieux C, Fabre S and
Drouilhet L (2020) Genetic
Determinism Exists for the Global
DNA Methylation Rate in Sheep.
Front. Genet. 11:616960.
doi: 10.3389/fgene.2020.616960
concerns epigenetic marks. Epigenetic modifications include biochemical modifications of DNA (methylation) or proteins tied to DNA (methylation and acetylation of histone proteins). These epigenetics processes act in an interrelated way to influence gene expression and hence phenotype in response to environmental conditions, thereby demonstrating that the eukaryotic genome can respond dynamically to changes in the environment to which every individual is exposed (Geiman and Robertson, 2002).

One of the key mechanisms in the regulation of gene expression in mammals is cytosine methylation at CpG dinucleotides, the most common DNA modification in eukaryotes. The importance of epigenetics in mammals and plants has been demonstrated in several studies and epigenetics contribute to the inheritance of traits of interest. In plants, DNA methylation in several differentially methylated regions was shown to explain 60–90% of the heritability of complex traits (Cortijo et al., 2014). In humans, DNA methylation explained 7–11% of the body mass index while SNP explained 5–8% of the same trait (Dick et al., 2014). In livestock, current studies aim to identify epi-loci underlying traits of agronomical interest. For instance, in sheep, genome-wide analysis of DNA methylation profiles has already revealed numerous epi-loci associated with prolificacy or body size (Cao et al., 2015; Zhang et al., 2017). Epigenetic sources of inheritance are not currently taken into consideration in selection strategies in livestock production but these studies could help identify the origin of variability of production traits and later on, improve genetic progress.

Taking into account variations in epigenetic marks may enable better modeling of phenotypic variation and improving the precision of genetic evaluation of traits and breeds. In addition, taking epigenetic events into account may help not only understand but also improve adaptive processes in both plants and mammals. Particular epigenetic recombinant inbred lines of *Arabidopsis thaliana* showed highly variable DNA methylation, while being genetically very similar (Latzel et al., 2013). Using these plant model lines, it was demonstrated an important role of epigenetics in adaptive process facing saline conditions for example (Kooke et al., 2015). Moreover, evolutionary studies showed that epigenetic inheritance may be important for natural selection and transmission of advantageous traits, for a review (see Donohue, 2014). Recent advances in epigenetic studies in natural populations addressed how epigenetic inheritance may influence adaptive evolution by focusing on epigenetic stability and inheritance itself as a potentially evolving trait (Herman et al., 2014; Schlichting and Wund, 2014). These articles brought into focus the genetic and ecological basis of epigenetic stability and raised a number of questions including “is the stability of epigenetic marks heritable?”

Turck and Coupland (2014) showed that DNA methylation has the potential to be mitotically and meiotically stable, whereas histone modification is involved in environmentally induced epigenetic regulation that may be reversible. These authors hypothesized that both the degree of DNA methylation and the stability of environmentally induced changes in gene expression via histone modification are associated with changes in the DNA sequence that can cause heritable variation in epigenetic regulation. New scientific questions arise concerning

the relationship between genetics and epigenetics as part of the improvement of livestock species and their adaptation to changing production systems. In the present study, we used a complementary approach to those used in existing studies by considering the global methylation of genomic DNA as a novel quantitative phenotype to provide essential information on the processes involved in epigenetic inheritance. DNA methylation in differentially methylated regions representing a very small proportion of genomic DNA global methylation, we hypothesized that the stable fraction of DNA methylation is genetically determined and can be inherited across generations. To test this hypothesis, we quantified the global DNA methylation rate in blood samples collected from domestic sheep. Sheep not only has the advantage of being a species of agronomical interest, but also represents an interesting animal model to study adaptive processes. The aim of the present study was to estimate the genetic component of the global DNA methylation rate and to identify the genomic loci underlying the phenotype.

# MATERIALS AND METHODS

## Animals and Management

The experimental animals were Romane lambs, a fixed crossbreed between Romanov x Berrichon du Cher (Ricordeau et al., 1992). Data collected from a total of 1,047 male and female lambs over a period of 5 years were used in this study. The lambs originated from 10 rams, each ram was used for only 1 year. All the animals were born in April, identified at birth using electronic ear tags, and reared outdoors with their dams under extensive conditions in the experimental “la Fage” farm on 280 ha of rangeland (Causse du Larzac, Roquefort sur Soulzon, South of France). The farming system, management and environment characteristics are detailed in Gonzalez-Garcia et al. (2014). The lambs were weighed regularly from lambing to weaning to estimate average daily gain. All the lambs were weaned at 85 ± 4 days of age. Two weeks after weaning, experimental lambs were individually exposed to two behavioral tests. Tests of each animal were all performed on the same day, with a total of around 35 lambs tested per day on 10–11 consecutive test-days per year. After the behavioral tests (i.e., approximately 1 month after weaning), lambs were blood sampled for genotyping and determination of the global DNA methylation rate.

## Global DNA Methylation Rate

The global DNA methylation rate (GDMR) was measured from whole blood samples taken from sheep. Genomic DNA was extracted from the blood samples using a high-salt extraction method (Roussot et al., 2003). Methylation analyses were performed using the LUMA assay (Karimi et al., 2006a,b). Briefly, this method relies on the use of methylation-sensitive and –insensitive restriction endonucleases: *HpaII* and *MspI*, respectively. The target sequence for both enzymes is CCGG, *HpaII* is not able to cut if the internal cytosine is methylated (C<sup>m</sup>CGG), whereas *MspI* cuts the restriction site whatever the methylation status. Moreover, *EcoRI* (recognition site: GAATTC)
was used as an internal control for normalization of DNA amount. For each sample, DNA was independently digested by EcoRI + HpaII and EcoRI + MspI restriction enzymes (New England Biolabs) and then analyzed using a Q24 Pyromark sequencer and Q24 1.0.10 software (Qiagen). The dispensation order for nucleotides was GTGTCACATGTGTG. Methylation rate was calculated from peak heights as [1 - [(HpaII(G)/ EcoRI_Hpa(T))/(MspI(G)/ EcoRI_Msp(T))] × 100]. The same internal control DNA sample was used for each pyrosequencing run. Finally, considering complete digestions, GDMR integrates data from nearly six million CCGG sites detected in silico on sheep reference genome OAR v3.1.

## Behavioral Traits

Two behavioral tests were used to evaluate sociability for conspecifics and reaction to a human. A complete description of the devices and test procedure used for behavioral measurements are given in Hazard et al. (2016). Briefly, the arena test (AT) aimed to evaluate sociability for conspecifics in three successive phases: (1) attraction for conspecifics with visual contact (Arena test phase 1, AT1), (2) reactivity to social separation from conspecifics by preventing visual contact between the tested lamb and conspecifics (Arena test phase 2, AT2), and (3) attraction for conspecifics in presence of a motionless human standing in front of conspecifics (Arena test phase 3, AT3). The variables used were the frequency of low-pitched bleats (i.e., the lamb bleated with its mouth closed) recorded during AT1 (AT1-LBLEAT), and the frequency of high-pitched bleats (i.e., the lamb bleated with its mouth open) recorded during AT2 (AT2-HBLEAT). Locomotor activity was assessed by measuring the number of zones crossed during AT2 (AT2-LOCOM). Vigilance postures (i.e., animal motionless, head in an upright position, and ears perpendicular to the head) were measured during AT2 (AT2-VIGIL). The time spent in each virtual zone was recorded and the ewe's proximity score to conspecifics during AT3 (AT3-PROX) was calculated by weighting the time spent in each zone.

The corridor test (CT) aimed to evaluate the reaction to a human. The second phase of the CT was used to evaluate the lamb's reaction to a walking human. The mean distance between the tested lamb and the walking human was recorded (CT2-DIST).

## Zootechnical Traits

At birth, the type of lamb birthcoat was graded on a scale of one to nine based on composition and structure of the coat, i.e., grade 1, single wooly coat with no halo- or coarse hair, hair length < 10 mm; grade 9, double hairy coat (coarse hair mixed with fine down) > 25 mm in length. Further details on this trait are provided in Allain et al. (2014). Lambs were weighed at birth (birth weight, kg) and weaning (weaning weight, kg). Lambs were also weighed regularly from lambing to weaning to estimate average daily gain. Growth rates (i.e., average daily gain, ADG, in g) were estimated in all lambs from birth to 30, from 30 to 60 and from 60 to 90 days of age (called ADG0-30, ADG30-60, ADG60-90, respectively). In addition, some carcass traits were measured at slaughter only for males, including dressing yield (DY, percentage), conformation score (CONF, score 0-6),

compactness (COMP) (i.e., width/length ratio of the carcass, percentage), external fat score (FAT Score, score 0-9), and back fat depth at the 12th rib (FAT depth, mm).

## Statistical Analysis

Analyses of variance were performed to assess the significant fixed effects affecting each measured trait. The tested effects were sex (male or female), born and reared litter size, age of the dam, age and weight of the lambs and year of measurement. The litter size effect was classified according to the number of lambs born and suckled (class 1, ewes lambing and suckling singletons; class 2, ewes lambing twins but suckling only one lamb; class 3, ewes lambing and suckling twins; and class 4, ewes lambing and suckling more than two lambs). The age of the dam effect included ewes that were one, two, and three or more than 3 years old (classes 1, 2, and 3, respectively). The year of measurement effect covered the 5 years over which data were collected. The GLM procedures of the SAS® software (SAS Institute Inc, 1999) were applied to the variables to test the fixed effects and first order interactions. The factors of variation showing significant effect (i.e., $P < 0.05$) on the considered trait were included in subsequent genetic analyses. The fixed effects considered depended on the trait. In addition, for the GDMR trait, we investigated differences between lamb families by testing a family effect (10 classes), nested in the year effect, according to the GLM procedure described above. This was only used for GDMR because the genetic variability of the other traits included in the present study has already been described (Bibé et al., 2002; Allain et al., 2014; Hazard et al., 2014).

## Genetic Analyses

The (co)variance components of each trait were estimated by restricted maximum likelihood applied to an animal model using ASREML 3.0 software (Gilmour et al., 2009). For each trait, the model included the appropriate fixed effects (i.e., sex, born and reared litter size, age of the dam, year of measurement) and a direct additive genetic effect of the animal considered as random effect. For some traits analyzed, a litter or a dam permanent environmental effect was added as random effect in the model based on present results or the results of previous studies that described an appropriate model for the trait concerned (Bibé et al., 2002; Allain et al., 2014; Hazard et al., 2014). For GDMR, weight and growth traits, the dam was considered as a permanent environmental effect. For birthcoat, the animal's litter was considered as a permanent environmental effect. The following complete animal mixed model was fitted:

$$y = X\beta + Z_{a}a + W_{c}c + e \quad (1)$$

where $y$ is the vector of observations for the trait(s) being analyzed, $\beta$ is the vector of appropriate fixed effects (sex of lambs, litter size born and suckled, age of the dam and year of measurement), $a$ is the vector of random genetic effects and $c$ is the vector of permanent environmental effects, when appropriate, with incidence matrices $X$, $Z_a$, and $W_c$, respectively, and $e$ is the vector of residual effects. Univariate analyses were performed to estimate the variance of each trait. Bivariate analyses were performed to estimate the genetic correlations between the
GDMR and the other traits using the same model as that used in univariate analyses.

## Genomic Analyses

Genomic analyses were only done for the GDMR trait. Nine of the 10 families were used (i.e., family 10 was not genotyped). Among the 1,047 lambs phenotyped for GDMR, 800 lambs were genotyped (i.e., after filtering an individual call rate ≥ 98%) as well as their nine respective sires using the Illumina OvineSNP50 BeadChip (i.e., 54,241 SNPs). Outlier individuals were removed, and genomic analyses were performed using 775 phenotyped and genotyped lambs. SNP quality was checked as described by Hazard et al. (2014) (i.e., SNP call rates < 97%, a minor allele frequency < 1%, inconsistent for Hardy-Weinberg disequilibrium). Sex chromosomes were not included for analysis. Finally, 40,725 autosomal SNPs were retained for QTL analyses.

Genome wide association analysis was performed of the whole population genotyped using joint analysis considering simultaneously linkage association and linkage disequilibrium (LDLA) to take advantage of both family structure and linkage disequilibrium (Legarra and Fernando, 2009) using QTLMAP software (Gilbert et al., 2008). The LDLA method consisted in interval mapping and LDLA model considered the sire haplotypes effects of the LD model in addition to sire QTL effect. A haplotype size of four SNP was used, and when haplotype frequency was less than 1%, haplotypes were assigned to a rare haplotype group. The chromosome-wise *p*-values were estimated assuming that, conditional on the QTL position, the likelihood ratio test statistics followed a χ<fnref n="2" />-distribution with k degrees of freedom, k being the number of genetic effects (Piepho, 2001). In our study, k was equal to the number of haplotypes plus the number of families for LDLA. Genome-wise *p*-values were obtained by correcting for multiple testing assuming that the number of independent tests was equal to the number of chromosomes analyzed (i.e., 26 independent tests corresponding to the 26 autosomes) (Benjamini and Hochberg, 1995; Knott et al., 1998). Confidence intervals were determined using the “2-LOD drop-off” criterion and assuming 1 LOD = 4.61 LRT (Lynch and Walsh, 1998).

For each QTL, the confidence interval expanded by 2 Mb on each side was used to extract gene annotation by using the Biomart tool from Ensembl release 101 of the sheep reference genome OAR v3.1¹. Functional annotations for epigenetic modifiers were extracted from dbEM<fnref n="2" /> (Nanda et al., 2016) and from the review about epigenetic modifiers by Feinberg et al. (2016).

## RESULTS

### Phenotypes and Genetic Analyses

Descriptive statistics of GDMR are summarized in Table 1. The mean value of GDMR was 70.71% and the coefficient of

variation was 8.40%. Analysis of the gender effect indicated that female lambs had higher GDMR than male lambs. Marked differences were observed between extreme values within a gender. Phenotypic variability of male lambs was higher than female lambs. No significant effects of litter size, age of the dam, age of the lamb and weight (at birth or at weaning) were found, but a significant effect of family (nested in the year effect) was found on GDMR. One family had a low GDMR at 66.60% (family 3) contrasting with four families having a higher GDMR ranging from 71.95 to 74.11% (families 1, 2, 7, and 8). Four families had intermediate GDMR ranging from 68.86 to 70.85% (families 4, 5, 6, 9, and 10). The difference between extreme families reached 1.25 SD. Differences in phenotypic variability were observed between families.

The estimated variance components of GDMR are listed in Table 2. The estimated heritability for GDMR was moderate (0.20 ± 0.05) and the permanent environment effect of the dam was very low. The additive genetic coefficient of variation for GDMR was also low (3.5%).

Descriptive statistics and genetic parameters for additional traits recorded in the lambs used in the present study are reported in Table 3. Heritability estimates for behavioral, zootechnical and carcass traits were generally moderate to high (0.20 ± 0.05 to 0.74 ± 0.07). Only low-pitched bleats and conformation score traits showed low heritability. No genetic correlation was found between GDMR and behavioral, birthcoat or zootechnical traits.

### Genomic Analyses

The significant QTLs found using LDLA and reaching genome-wide (GW) or chromosome-wide (CW) significance thresholds are reported in Table 4 and Figure 1. Four haplotype-trait associations reached the GW significance thresholds and 16 haplotype-trait associations reached the 1% CW significance threshold. Eight significant associations reaching the 5% CW threshold were also detected. Associations found in LDLA mapped to 14 chromosomes. The four associations that reached the GW threshold of significance localized on chromosomes 1 (OAR1, 13.20 and 179.50 Mb) and 17 (OAR17, 40.14 and 55.54 Mb). Associations reaching the 1% CW significance threshold were mapped on chromosomes OAR3, 5, 13, 14, 16, 17, 18, 21, 23, 24, and 26. Effects of QTLs found in LDLA ranged from 0.68 to 3.20% of the total phenotypic variance of GDMR.

By crossing information from sheep reference genome annotation (OAR v3.1), the database of epigenetic modifiers (dbEM, (Nanda et al., 2016) and those reported by Feinberg et al. (2016), we highlighted *CHD1*, *NCO3A*, *KDM8*, *KAT7*, and *KATA6* as the most obvious candidate genes from each QTL (+2 Mb each side of the confidence interval), likely to play a role in the GDMR variation (Table 4).

## DISCUSSION

The global DNA methylation rate measured by LUMA reached 70% in sheep blood. This result is consistent with global DNA methylation percentages reported in human blood analyzed by LUMA (for a review, see Soriano-Tarraga et al., 2013). The

<footnote n="1">¹http://www.ensembl.org</footnote>

<footnote n="2">²http://crdd.osdd.net/raghava/dbem</footnote>
<table><caption>**TABLE 1** | Number of records, mean (± SD), minimum, maximum, quartiles, and least squares means for global DNA methylation rate in the whole population and for each sex and each family.</caption><thead><tr><th></th><th></th><th>Number of lambs</th><th>Mean (±SD)</th><th>Min</th><th>Max</th><th>Q1</th><th>Q3</th><th>P-value</th><th>LSMeans</th></tr></thead><tbody><tr><td>Population</td><td></td><td>1047</td><td>70.71 (5.97)</td><td>23.08</td><td>87.94</td><td>68.72</td><td>74.24</td><td></td><td></td></tr><tr><td rowspan="2">Sex</td><td>Male</td><td>426</td><td>68.75 (7.21)</td><td>23.08</td><td>80.63</td><td>66.54</td><td>72.60</td><td>***</td><td>69.12</td></tr><tr><td>Female</td><td>621</td><td>72.05 (4.47)</td><td>40.61</td><td>87.94</td><td>70.01</td><td>74.85</td><td></td><td>72.17</td></tr><tr><td rowspan="10">Family (year)</td><td>1</td><td>93</td><td>73.40 (4.37)</td><td>55.82</td><td>81.28</td><td>71.17</td><td>76.81</td><td>***</td><td>73.05</td></tr><tr><td>2</td><td>77</td><td>74.73 (3.55)</td><td>60.56</td><td>80.79</td><td>72.82</td><td>77.21</td><td></td><td>74.11</td></tr><tr><td>3</td><td>116</td><td>66.94 (8.15)</td><td>27.04</td><td>75.34</td><td>65.34</td><td>71.46</td><td></td><td>66.63</td></tr><tr><td>4</td><td>91</td><td>71.27 (6.7)</td><td>23.08</td><td>79.85</td><td>69.83</td><td>75.05</td><td></td><td>70.85</td></tr><tr><td>5</td><td>115</td><td>68.95 (4.95)</td><td>37.96</td><td>77.06</td><td>66.24</td><td>72.05</td><td></td><td>68.86</td></tr><tr><td>6</td><td>117</td><td>69.40 (5.62)</td><td>43.39</td><td>87.94</td><td>66.82</td><td>72.18</td><td></td><td>69.29</td></tr><tr><td>7</td><td>101</td><td>73.00 (3.89)</td><td>55.01</td><td>77.89</td><td>71.30</td><td>75.87</td><td></td><td>72.57</td></tr><tr><td>8</td><td>108</td><td>72.20 (3.31)</td><td>61.02</td><td>78.83</td><td>70.52</td><td>74.43</td><td></td><td>71.95</td></tr><tr><td>9</td><td>123</td><td>69.92 (5.60)</td><td>41.62</td><td>79.56</td><td>68.65</td><td>73.03</td><td></td><td>69.81</td></tr><tr><td>10</td><td>106</td><td>69.64 (6.73)</td><td>28.03</td><td>77.22</td><td>68.19</td><td>73.32</td><td></td><td>69.33</td></tr></tbody></table>

<table><caption>**TABLE 2** | Estimates of heritability, repeatability, permanent and residual effects (±S.E.) for the global DNA methylation rate.</caption><thead><tr><th rowspan="2"></th><th rowspan="2">n</th><th colspan="3">Component</th><th rowspan="2">Repeatability (r<sup>2</sup>)</th><th rowspan="2">Total σ<sub>p</sub><sup>2</sup></th></tr><tr><th>Genetic (h<sup>2</sup>)</th><th>Permanent (m<sup>2</sup>)</th><th>Residual (e<sup>2</sup>)</th></tr></thead><tbody><tr><td>GDMR</td><td>1,047</td><td>0.20 ± 0.05</td><td>0.02 ± 0.03</td><td>0.78 ± 0.05</td><td>0.22 ± 0.05</td><td>31.3 ± 1.5</td></tr></tbody></table>

<table><caption>**TABLE 3** | Number of records, mean (±SD), estimates of heritability, and permanent environmental effects for behavioral, zootechnical and carcass traits, and genetic correlations between these traits and the global DNA methylation rate.</caption><thead><tr><th></th><th>n</th><th>Mean (±SD)</th><th>h<sup>2</sup></th><th>p<sup>2</sup></th><th>r<sub>g</sub></th></tr></thead><tbody><tr><td>AT1-LBLEAT</td><td>1,047</td><td>0.28 (0.28)</td><td>0.17 (0.07)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>AT2-HBLEAT</td><td>1,047</td><td>2.91 (1.39)</td><td>0.46 (0.07)</td><td>ne</td><td>0.05 (0.03)</td></tr><tr><td>AT2-LOCOM</td><td>1,017</td><td>19.44 (8.5)</td><td>0.26 (0.07)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>AT2-VIGIL</td><td>1,017</td><td>20.68 (9.81)</td><td>0.32 (0.07)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>AT3-PROX</td><td>990</td><td>26.26 (16.89)</td><td>0.27 (0.06)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>CT2-DIST</td><td>1,047</td><td>5.37 (1.14)</td><td>0.20 (0.06)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>Birthcoat (score 0-9)</td><td>1,032</td><td>6.06 (3.38)</td><td>0.74 (0.07)</td><td>0.15 (0.04)</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>Birth weight (kg)</td><td>1,047</td><td>3.9 (0.8)</td><td>0.26 (0.09)</td><td>0.27 (0.05)</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>Weaning weight (kg)</td><td>1,043</td><td>23.3 (4.2)</td><td>0.36 (0.09)</td><td>0.20 (0.04)</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>ADG0-30 (g)</td><td>1,047</td><td>360.3 (73.04)</td><td>0.30 (0.08)</td><td>0.24 (0.04)</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>ADG30-60 (g)</td><td>1,047</td><td>227.4 (50.5)</td><td>0.34 (0.08)</td><td>0.06 (0.04)</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>ADG60-90 (g)</td><td>1,046</td><td>216.0 (58.0)</td><td>0.32 (0.07)</td><td>0.10 (0.04)</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>DY (%)</td><td>386</td><td>42.07 (2.19)</td><td>0.41 (0.15)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>CONF (score 0-6)</td><td>390</td><td>3.44 (0.8)</td><td>0.09 (0.08)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>COMP (%)</td><td>388</td><td>32.26 (1.46)</td><td>0.58 (0.20)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>FAT score (score 0-9)</td><td>390</td><td>5.07 (1.39)</td><td>0.22 (0.12)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr><tr><td>FAT depth (mm)</td><td>389</td><td>1.87 (1.36)</td><td>0.36 (0.13)</td><td>ne</td><td>&lt;0.05 (&lt;0.01)</td></tr></tbody></table>

value of GDMR in the present study is also in the range of calculated DNA methylation either using the global CpG methylation rate (50–55%) measured by reduced representation

bisulfite sequencing (RRBS) in sheep muscle (Couldrey et al., 2014) or the total methyl content measured by high performance liquid chromatography reported in somatic tissues (80%) in

<footnote>***P-value < 0.001. Min, minimum; max, maximum; Q1, first quartile; Q3, third quartile; LSMeans, least squares means; Family effect was nested in the year effect.</footnote>

<footnote>h², m², e², proportion of total phenotypic variance attributed to additive genetic, permanent, and residual effects, respectively, Total σ<sub>p</sub><sup>2</sup>, total phenotypic variance; r² = h² + m²; n, number of animals.</footnote>

<footnote>AT1/2/3, arena test phase 1/2/3; CT2, corridor test phase2; LBLEAT, low bleats; HBLEAT, high bleats; LOCOM, locomotion; VIGIL, vigilance behavior; PROX, proximity score; DIST, mean distance between tested animal and a human; ADG, average daily gain; DY, dressing yield; CONF, conformation score; COMP, compactness. h², p², proportion of total phenotypic variance attributed to additive genetic or permanent environmental effects, respectively, r<sub>g</sub>, genetic correlations; standard errors are in parentheses; ne, not estimated.</footnote>
<table><caption>**TABLE 4** | List of QTLs detected in joint analysis (linkage disequilibrium linkage analysis) associated with the global DNA methylation rate.</caption><thead><tr><th>Chr.</th><th>Significance</th><th>level<sup>a</sup></th><th>QTL effect (%)<sup>b</sup></th><th>Confidence interval (Mb)</th><th>Flanking markers<sup>c</sup></th><th>Annotated protein coding genes within confidence interval</th><th>Epigenetic modifier candidate</th></tr></thead><tbody><tr><td>1</td><td><strong>GW</strong></td><td>*</td><td>2.48</td><td>13.1-13.4</td><td>OAR1_13033760.1</td><td>s62110.1</td><td>RRAGC, GJA9, RHDB</td></tr><tr><td>1</td><td>CW</td><td>*</td><td>2.75</td><td>69.3-69.5</td><td>OAR1_74015424.1</td><td>OAR1_74065628.1</td><td>MTF2</td></tr><tr><td>1</td><td>CW</td><td>*</td><td>3.0</td><td>138.0-138.3</td><td>OAR1_149183890.1</td><td>OAR1_149346875.1</td><td>-</td></tr><tr><td>1</td><td><strong>GW</strong></td><td>*</td><td>2.4</td><td>179.4-179.6</td><td>OAR1_193623565.1</td><td>OAR1_193664187.1</td><td>-</td></tr><tr><td>2</td><td>CW</td><td>*</td><td>1.44</td><td>199.7-199.9</td><td>OAR2_211501156.1</td><td>OAR2_211573792.1</td><td>-</td></tr><tr><td>3</td><td>CW</td><td>**</td><td>1.0</td><td>114.5-114.7</td><td>OAR3_122117752.1</td><td>OAR3_122207054.1</td><td>-</td></tr><tr><td>5</td><td>CW</td><td>**</td><td>0.83</td><td>96.6-96.9</td><td>OAR5_105186382.1</td><td>OAR5_105374983.1</td><td>-</td></tr><tr><td>11</td><td>CW</td><td>*</td><td>2.95</td><td>10.5-10.9</td><td>s74938.1</td><td>OAR11_10495137.1</td><td>MAD13, INTS2, BRIP1</td></tr><tr><td>11</td><td>CW</td><td>*</td><td>0.85</td><td>37.4-37.6</td><td>OAR11_39854191.1</td><td>OAR11_39948154.1</td><td>HOXB1/2/3/4/5/6/7/9</td></tr><tr><td>13</td><td>CW</td><td>*</td><td>2.9</td><td>17.8-18.5</td><td>OAR13_20537279.1</td><td>OAR13_20578061.1</td><td>PARD3</td></tr><tr><td>13</td><td>CW</td><td>*</td><td>2.7</td><td>32.1-32.3</td><td>OAR13_35499997.1</td><td>OAR13_35648372.1</td><td>NSUN6, EPC1</td></tr><tr><td>13</td><td>CW</td><td>**</td><td>1.17</td><td>74.5-74.7</td><td>s65442.1</td><td>s62908.1</td><td>CDH22, SLC35C2, ELMO2</td></tr><tr><td>14</td><td>CW</td><td>**</td><td>2.45</td><td>7.7-7.9</td><td>s16689.1</td><td>OAR14_8081413.1</td><td>CMIP, PLCG2, SNORA70</td></tr><tr><td>14</td><td>CW</td><td>**</td><td>2.06</td><td>10.7-10.9</td><td>s26487.1</td><td>s43948.1</td><td>-</td></tr><tr><td>14</td><td>CW</td><td>**</td><td>2.77</td><td>14.8-15.2</td><td>OAR14_15219778.1</td><td>OAR14_15268863.1</td><td>GTP2, DNAJA2, NETO2, ITFG1</td></tr><tr><td>14</td><td>CW</td><td>**</td><td>2.85</td><td>17.5-18.3</td><td>s22016.1</td><td>s31525.1</td><td>ZNF423, CNEP1R1, HEATR3, TENT4B, ADCY7, BRD7</td></tr><tr><td>16</td><td>CW</td><td>**</td><td>2.38</td><td>4.5-4.8</td><td>OAR16_4698169.1</td><td>OAR16_4737110_X.1</td><td>DUSP1, ERGIC1, RPL26L1, ATP6V0E1</td></tr><tr><td>17</td><td>CW</td><td>**</td><td>2.7</td><td>26.7-26.9</td><td>OAR17_29456634.1</td><td>s52728.1</td><td>-</td></tr><tr><td>17</td><td><strong>GW</strong></td><td>***</td><td>2.71</td><td>40.0-40.2</td><td>OAR17_43310387.1</td><td>OAR17_43356226.1</td><td>PPID, ETFDH, C4orf46, RXFP1</td></tr><tr><td>17</td><td><strong>GW</strong></td><td>*</td><td>2.8</td><td>55.4-55.6</td><td>OAR17_60522910.1</td><td>s57641.1</td><td>CCDC60</td></tr><tr><td>18</td><td>CW</td><td>**</td><td>1.68</td><td>6.6-6.8</td><td>OAR18_6256173.1</td><td>s37275.1</td><td>TTC23, SYNM</td></tr><tr><td>21</td><td>CW</td><td>**</td><td>0.80</td><td>3.7-4.1</td><td>OAR21_4678890.1</td><td>s05427.1</td><td>-</td></tr><tr><td>23</td><td>CW</td><td>**</td><td>0.68</td><td>8.3-8.5</td><td>OAR23_9159961.1</td><td>OAR23_9204273.1</td><td>-</td></tr><tr><td>23</td><td>CW</td><td>**</td><td>2.7</td><td>34.4-34.6</td><td>OAR23_36433223.1</td><td>s07796.1</td><td>GATA6</td></tr><tr><td>24</td><td>CW</td><td>*</td><td>2.5</td><td>7.4-7.7</td><td>OAR24_8719574.1</td><td>s03253.1</td><td>TMEM114, METTL22, ABAT, TMEM186, PMM2, CARHSP1</td></tr><tr><td>24</td><td>CW</td><td>**</td><td>3.2</td><td>23.6-24.0</td><td>s12637.1</td><td>OAR24_26210500.1</td><td>-</td></tr><tr><td>24</td><td>CW</td><td>**</td><td>1.08</td><td>32.6-33.0</td><td>s64882.1</td><td>s17643.1</td><td>RCC1L, NCF1, GTF2IRD1, CLIP2</td></tr><tr><td>26</td><td>CW</td><td>**</td><td>2.0</td><td>33.4-33.6</td><td>OAR26_38152631.1</td><td>OAR26_38251830.1</td><td>ADAM18, IDO1, IDO2</td></tr></tbody></table>

Significant QTLs reaching the chromosome-wide or the genome-wide thresholds are listed. <sup>a</sup>*p* < 5%; <sup>**</sup>*p* < 1%; <sup>***</sup>*p* < 0.1%. <sup>b</sup>Percentage of phenotypic variance explained by the haplotype with the most effect (absolute value). <sup>c</sup>SNPs flanking the QTL with significant association. CW, chromosome wide; GW, genome wide; Position and genes are based on the reference genome Oar_v3.1.
<figure><img src="image_3.png" /><figcaption>**FIGURE 1** | Chromosome plots of the likelihood ratio test values obtained for global DNA methylation rate by LDLA analysis. The likelihood ratio test (LRT) is plotted against SNP haplotype positions (four consecutive SNPs) for chromosomes 1, 2, 3, 5, 11, 13, 14, 16, 17, 18, 21, 23, 24, and 26. The horizontal lines indicate the 5, 1, and 0.1% genome wide thresholds. One percent chromosome wide thresholds were OAR3, 47.6; OAR5, 37.4; OAR13, 42.6; OAR14, 35.9; OAR16, 38.7; OAR17, 42.3; OAR18, 40.6; OAR21, 37.3; OAR23, 45.7; OAR24, 42.8, and OAR26, 47.1. 5% chromosome wide thresholds were OAR1, 45.5; OAR2, 45.0; OAR11, 41.3; OAR13, 38.3; OAR24, 38.4.</figcaption></figure>

humans (Ehrlich et al., 1982). The large range observed for DNA methylation rate across the genome is likely due to the methods used for quantification of DNA methylation and/or the tissues studied. Indeed, in their review, Soriano-Tarraga et al. (2013) reported that variation in global DNA methylation in human blood analyzed by LUMA ranged from 52 to 78% depending on the DNA isolation method used. In addition, Doherty and Couldrey (2014) reported that RRBS libraries are biased because they contain promoter regions which have high CG content known to be largely devoid of DNA methylation. Consequently, RRBS libraries are expected to display lower methylation on average across genome than unbiased libraries.

Variability was observed for the GDMR phenotype recorded in this study. Part of this variability originated from differences between male and female lambs. Sex differences in the methylome have already been reported in humans. In human

liver, females displayed higher average DNA methylation in the X-chromosome whereas males presented higher DNA methylation in autosomes in this tissue (García-Calzón et al., 2018). Genome-wide sex differences in locus-specific DNA methylation across autosomes have also been reported in other human tissues such as blood (Singmann et al., 2015). In the present study, the positions of DNA methylation across genome were unknown, so we are unable to conclude whether gender difference in DNA methylation found in lambs was of similar magnitude on the X-chromosome and autosomes. In addition, it should be kept in mind that we investigated DNA methylation in young lambs (4 months of age), not in adult sheep. A study investigating gender differences in DNA methylation in human blood at birth reported that CpG located on autosomes, as well as DMR, were hypermethylated in girls compared to boys (Yousefi et al., 2015). Because gender differences in DNA methylation
across autosomes in human blood varied between newborns and adults, we cannot exclude that a sex by age interaction could affect DNA methylation in sheep blood.

We hypothesized that GDMR has a genetic component in sheep. The variance component decomposition analysis of GDMR revealed moderate heritability (20%). To our knowledge, this is the first heritability estimate for GDMR in livestock species. Estimation of heritability of DNA methylation levels in human whole blood at thousands of sites has shown that heritability estimates vary across the genome but is on average 20%, consistently with our estimate in sheep blood (van Dongen et al., 2016). Additional studies in humans also indicated that the average twin based heritability of DNA methylation across genome-wide CpGs varied between 5 and 19% depending on the tissue (i.e., $h^2$ = 5, 7, and 12% in placenta, human umbilical vascular endothelial cells and cord blood mononuclear cells from newborns twins, respectively, $h^2$ = 18% in peripheral blood lymphocytes from adolescent twins [Gordon et al., 2012; McRae et al., 2014]). Interestingly, Shah et al. (2014) reported that DNA methylation in whole blood remained stable over the human lifetime, as indicated by the high correlation (0.68) between DNA methylation repeatability in older people (median age 70–90) and heritability estimated in teenagers.

Genetic relationships between GDMR and production or adaptive traits were also investigated in the present study. We did not find any relationship between GDMR and production or adaptive traits in domestic sheep. Our results contrasted with the negative correlation between total DNA methylation level and daily growth rate in fish embryos (Ou et al., 2019). To our knowledge, no published works provide evidence of a relationship between GDMR in blood and the specific adaptive traits used in the present study (social behaviors and birthcoat). Some relationships would be expected, since several studies suggested epigenetic influence on adaptive processes in mammals (Massicotte et al., 2011; Lea et al., 2016). Nevertheless, only a few adaptive traits were explored in this study and we cannot exclude the possibility that such relationships exist for other adaptive traits or in older animals. Additionally, regardless of the absence of relationships between GDMR and adaptive traits in our study, we cannot exclude, by considering individually DNA methylation sites, to find epi-loci associated with production or adaptive traits investigated in the present study. Indeed, epi-loci have previously been reported to be associated with animal behaviors (Champagne and Curley, 2009; Roth, 2013; Verhulst et al., 2016; Seebacher and Krause, 2019).

Even if the greater attention paid to epigenetics in farm animals has led to an increasing number of GWA epigenetic studies, genomic studies of DNA methylation rate had remained to be performed. In our study, we were particularly interested in detecting QTLs for GDMR in domestic sheep. The QTL analysis for GDMR resulted in mapping 28 QTLs on 14 chromosomes, with four regions (on OAR1 and OAR17) reaching the genome wide level of significance. Each detected QTL explained less than 3% of the phenotypic variance. Our experimental design enabled to reach a power of 90% for the detection of a QTL explaining 8% of the phenotypic variance while the power

was 20% for a QTL that explains 1% of phenotypic variance (Luo, 1998). Consequently, the undetected QTLs likely had a limited effect on GDMR. Thus, the low proportion of variance explained by each QTL supports the hypothesis that GDMR in domestic sheep is under polygenic influence and unlikely under the control of a major gene. Nevertheless, several QTLs mapped in the present study could probably act in combination to account for substantial genetic variation of GDMR.

To better understand the underlying molecular mechanisms associated with these QTLs, we looked for possible overlapping between the location of the identified QTLs (confidence interval enlarged by 2 Mb each side) and the location of annotated genes relevant to epigenetics. Despite the large number of genes located close to the QTL regions identified, only few genes seem relevant as epigenetic modifiers through their action on methylation or acetylation processes on DNA or histones. Indeed, by crossing the list of positional genes in enlarged QTL regions with that of known epigenetic modifiers either from the dbEM database (Nanda et al., 2016) or those listed in Feinberg et al. (2016), only *CHD1*, *KDM8*, *NCOA3*, *KAT7*, and *KATA6* appeared as the most obvious candidate genes. *CHD1* gene encodes the chromodomain-helicase-DNA-binding protein 1 (OAR5: 95.4 Mb), that functions as substrate recognition component of the transcription regulatory histone acetylation complex SAGA. Mutations in *CHD1* are frequently associated with prostate cancers (Berger et al., 2011). Functional analyses have shown that *CHD1* was involved in transient DNA methylation and several loci were hypermethylated in *CHD1* deleted strains of Neurospora crassa (Belden et al., 2011). *KDM8* (Lysin demethylase 8 or Jumonji C domain-containing demethylase 5, JDMD5; OAR24: 24.8 Mb) is a histone H3 demethylase with specificity for Lys-36. This demethylase promotes homologous recombination (Amendola et al., 2017), cell proliferation (Hsia et al., 2010) and protects from nerve demyelination (Fuhrmann et al., 2018). Interestingly, a second gene encoding an important histone lysine demethylase, *KDM2B* (Lysine-specific demethylase 2B, OAR17: 53.23 Mb), was found close to one of the enlarged QTL regions (OAR17: 55.4–55.6 Mb). In this region, we also identified *SETD1B* gene (Histone-lysine N-methyltransferase SETD1B, 53.06 Mb) which encodes a component of a histone methyltransferase complex that specifically methylates Lys-4 of histone H3 and is responsible for the epigenetic control of chromatin structure and gene expression. A specific hypermethylation signature was associated with loss of function mutations in the *SETD1B* gene (Krzyzewska et al., 2019). The *NCOA3* gene (nuclear receptor coactivator 3; OAR13: 32.0 Mb) encodes a nuclear receptor coactivator with histone acetyltransferase activity. It particularly interacts with other transcriptional activators such p300/CBP-associated factor and CREB binding protein (CREBBP) as part of a multi-subunit coactivation complex. The *CREBBP* gene, also located near a significant QTL position (OAR 24: 3.1 Mb) has been recently reported to be involved in the crosstalk between DNA methylation and histone acetylation
(Zhang et al., 2020). Finally, two lysine acetyltransferase coding genes of the MYST family, *KAT6A* (lysine acetyltransferase 6A, also named *MOZ*, OAR26: 35.3 Mb) and *KAT7* (lysine acetyltransferase 7, also named *HBO1*, OAR11: 36.4 Mb) has been reported to play a key role in acetylation of lysine residues in histone H3 and/or H4 (Kitabayashi et al., 2001). Interestingly, KDM8, NCOA3 and KAT7 act as androgen or estrogen receptor co-regulators (Sharma et al., 2000; Vija et al., 2013; Wagner et al., 2013; Wang et al., 2019). Polymorphisms that alter the function of these genes could explain that female lambs had higher GDMR than male lambs.

The level and state of histone acetylation regulated by histone acetyl transferase enzymes are widely reported to be an epigenetic regulation of gene expression (Kouzarides, 2007). Therefore, we cannot exclude the possibility that global DNA methylation in sheep is associated with histone acetylation, since relationships between DNA methylation and histone de-acetylation have already been described (Geiman and Robertson, 2002; Vaissière et al., 2008). Several other coding genes, located close to QTLs we identified, have been reported to be involved in the regulation of transcription or regulation of chromatin architecture. Although no scientific evidence has suggested their involvement in DNA methylation to date, these genes could mediate this process through their active role in the coregulatory complexes that associate epigenetic modifiers to modulate transcription.

## CONCLUSION

Inter-individual variability in global DNA methylation rate in blood has an additive genetic component in sheep and heritability was moderate (h² = 0.20 ± 0.05). This opens the way for a possible genetic selection of this trait and create experimental divergent lines to investigate further the relationship between GDMR and adaptive traits as already evidenced in plants (Johannes et al., 2009; Kooke and Keurentjes, 2015). Moreover, this work reported the first SNP-based QTL detection study for GDMR as a quantitative trait in livestock species. The evidence of 28 QTLs associated with blood GDMR, each explaining a small proportion of the phenotypic variance (between 1 and 3%), most likely indicated a polygenic determinism of this trait. The further identification of genes and their polymorphisms underlying the global DNA methylation rate assessed in our study will paves the way for a deeper understanding the genetic component of such a trait.

## REFERENCES

* Allain, D., Foulquié, D., Autran, P., Francois, D., and Bouix, J. (2014). Importance of birthcoat for lamb survival and growth in the romane sheep breed extensively managed on rangelands. *J. Anim. Sci.* 92, 54–63. doi: 10.2527/jas.20 13-6660

* Amendola, P. G., Zaghet, N., Ramalho, J. J., Johansen, J. V., Boxem, M., and Salcini, A. E. (2017). JMJD-5/KDM8 regulates H3K36me2 and is required for late steps of homologous recombination and genome integrity. *PLoS Genet.* 13:e1006632. doi: 10.1371/journal.pgen.1006632

# DATA AVAILABILITY STATEMENT

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

# ETHICS STATEMENT

The animal study was reviewed and approved by the Local ethics committee SCIENCE ET SANTE ANIMALES N°115 (Ecole Nationale Vétérinaire de Toulouse).

# AUTHOR CONTRIBUTIONS

DH, LD, CM-R, and SF conceived and designed the genetic experimental design and were involved in planning the study. FP-P was in charge of the phenotyping of GDMR. DH performed statistical and genetic analyses. DH, LD, SF, FP-P, and CM-R analyzed and interpreted the results. DH drafted the manuscript. LD, SF, and CM-R helped in writing, reviewing, and editing the final manuscript. All authors have read and approved the final manuscript.

# FUNDING

This study was granted by the INRAE metaprogram SELGEN (GENEPIGEN project), particularly for the phenotyping of GDMR. Genotyping was funded by ANR and APIS-GENE organizations in the frame of SheepSNPQTL project (ANR-08-GENM-039). The INRAE Division of Animal Genetics supported computational facilities and experimental costs for the production of the animals and the collection of behavioral, zootechnical and carcass traits.

# ACKNOWLEDGMENTS

We are indebted to Paul Autran for the management of the experimental farm, Didier Foulquié and Fabien Carrière for managing the experimental flock, for animal care as well as for their active role in collecting experimental data. We also thank Alain Boissy and Eric Delval (INRAE UMR Herbivores) for their involvement in behavioral experiments, and Dominique François for his involvement in the management of animal breeding.

* Belden, W. J., Lewis, Z. A., Selker, E. U., Loros, J. J., and Dunlap, J. C. (2011). CHD1 remodels chromatin and influences transient DNA methylation at the clock gene frequency. *PLoS Genet.* 7:e1002166. doi: 10.1371/journal.pgen.1002166

* Benjamini, Y., and Hochberg, Y. (1995). Controlling the false discovery rate - a practical and powerful approach to multiple testing. *J. R. Stat. Soc. Ser. B Methodol.* 57, 289–300. doi: 10.1111/j.2517-6161.1995.tb02031.x

* Berger, M. F., Lawrence, M. S., Demichelis, F., Drier, Y., Cibulskis, K., Sivachenko, A. Y., et al. (2011). The genomic complexity of primary human prostate cancer. *Nature* 470, 214–220. doi: 10.1038/nature09744
Bibé, B., Brunel, J., Bourdillon, Y., Laradoux, D., Gordy, M., Weisbecker, J., et al. (2002). "Genetic parameters of growth and carcass quality of lambs at the French progeny-test station Berrytest," in *Proceedings of the 7th World Congress on Genetics Applied to Livestock Production*, (Montpellier).

Cao, J. X., Wei, C. H., Liu, D. M., Wang, H. H., Wu, M. M., Xie, Z. Y., et al. (2015). DNA methylation Landscape of body size variation in sheep. *Sci. Rep.* 5:13950. doi: 10.1038/srep13950

Champagne, F. A., and Curley, J. P. (2009). Epigenetic mechanisms mediating the long-term effects of maternal care on development. *Neurosci. Biobehav. Rev.* 33, 593–600. doi: 10.1016/j.neubiorev.2007.10.009

Cortijo, S., Wardenaar, R., Colomé-Tatché, M., Gilly, A., Etcheverry, M., Labadie, K., et al. (2014). Mapping the epigenetic basis of complex traits. *Science* 343, 1145–1148. doi: 10.1126/science.1248127

Couldrey, C., Brauning, R., Bracegirdle, J., Maclean, P., Henderson, H. V., and McEwan, J. C. (2014). Genome-wide DNA methylation patterns and transcription analysis in sheep muscle. *PLoS One* 9:101853. doi: 10.1371/journal.pone.0101853

Danchin, E., Charmantier, A., Champagne, F. A., Mesoudi, A., Pujol, B., and Blanchet, S. (2011). Beyond DNA: integrating inclusive inheritance into an extended theory of evolution. *Nat. Rev. Genet.* 12, 475–486. doi: 10.1038/nrg3028

Dick, K. J., Nelson, C. P., Tsaprouni, L., Sandling, J. K., Aïssi, D., Wahl, S., et al. (2014). DNA methylation and body-mass index: a genome-wide analysis. *Lancet* 383, 1990–1998. doi: 10.1016/s0140-6736(13)62674-62674

Doherty, R., and Couldrey, C. (2014). Exploring genome wide bisulfite sequencing for DNA methylation analysis in livestock: a technical assessment. *Front. Genet.* 5:126. doi: 10.3389/fgene.2014.00126

Donohue, K. (2014). THE epigenetics of adaptation: focusing on epigenetic stability as an evolving trait. *Evolution* 68, 617–619. doi: 10.1111/evo.12347

Ehrlich, M., Gama-Sosa, M. A., Huang, L.-H., Midgett, R. M., Kuo, K. C., McCune, R. A., et al. (1982). Amount and distribution of 5-methylcytosine in human DNA from different types of tissues or cells. *Nucleic Acids Res.* 10, 2709–2721. doi: 10.1093/nar/10.8.2709

Feinberg, A. P., Koldobskiy, M. A., and Gondor, A. (2016). Epigenetic modulators, modifiers and mediators in cancer aetiology and progression. *Nat. Rev. Genet.* 17, 284–299. doi: 10.1038/nrg.2016.13

Fuhrmann, D., Mernberger, M., Nist, A., Stiewe, T., and Elsasser, H. P. (2018). Miz1 controls schwann cell proliferation via H3K36(me2) demethylase Kdm8 to prevent peripheral nerve demyelination. *J. Neurosci.* 38, 858–877. doi: 10.1523/jneurosci.0843-17.2017

García-Calzón, S., Perfilyev, A., de Mello, V. D., Pihlajamäki, J., and Ling, C. (2018). Sex differences in the methylome and transcriptome of the human liver and circulating HDL-cholesterol levels. *J. Clin. Endocrinol. Metab.* 103, 4395–4408. doi: 10.1210/jc.2018-2423

Geiman, T. M., and Robertson, K. D. (2002). Chromatin remodeling, histone modifications, and DNA methylation—how does it all fit together? *J. Cell. Biochem.* 87, 117–125. doi: 10.1002/jcb.10286

Gilbert, H., Le Roy, P., Moreno, C., Robelin, D., and Elsen, J. M. (2008). QTLMAP, a software for QTL detection in outbred populations. *Ann. Hum. Genet.* 72, 694–694.

Gilmour, A. R., Gogel, B. J., Cullis, B. R., and Thompson, R. (2009). *ASReml User Guide Release 3.0*. Hemel Hempstead: VSN International Ltd.

Goddard, M. E., Hayes, B. J., and Meuwissen, T. H. E. (2010). Genomic selection in livestock populations. *Genet. Res.* 92, 413–421. doi: 10.1017/S0016672310000613

Gonzalez-Garcia, E., de Figuereido, V. G., Foulquie, D., Jousserand, E., Autran, P., Camous, S., et al. (2014). Circannual body reserve dynamics and metabolic profile changes in Romane ewes grazing on rangelands. *Domest. Anim. Endocrinol.* 46, 37–48. doi: 10.1016/j.domaniend.2013.10.002

Gordon, L., Joo, J. E., Powel, J. E., Ollikainen, M., Novakovic, B., Li, X., et al. (2012). Neonatal DNA methylation profile in human twins is specified by a complex interplay between intrauterine environmental and genetic factors, subject to tissue-specific influence. *Genome Res.* 22, 1395–1406. doi: 10.1101/gr.136598.111

Hazard, D., Bouix, J., Chassier, M., Delval, E., Foulquie, D., Fassier, T., et al. (2016). Genotype by environment interactions for behavioral reactivity in sheep. *J. Anim. Sci.* 94, 1459–1471. doi: 10.2527/jas2015-0277

Hazard, D., Moreno, C., Foulquié, D., Delval, E., François, D., Bouix, J., et al. (2014). Identification of QTLs for behavioral reactivity to social separation and humans in sheep using the OvineSNP50 BeadChip. *BMC Genom.* 15:778. doi: 10.1186/1471-2164-15-778

Herman, J. H., Spencer, H. G., Donohue, K., and Sultan, S. E. (2014). How stable 'should' epigenetic modifications be? Insights from adaptive plasticity and bet hedging. *Evolution* 68, 632–643. doi: 10.1111/evo.12324

Hsia, D. A., Tepper, C. G., Pochampalli, M. R., Hsia, E. Y. C., Izumiya, C., Huerta, S. B., et al. (2010). KDM8, a H3K36me2 histone demethylase that acts in the cyclin A1 coding region to regulate cancer cell proliferation. *Proc. Natl. Acad. Sci. U. S. A.* 107, 9671–9676. doi: 10.1073/pnas.1000401107

Johannes, F., Colot, V., and Jansen, R. C. (2008). Epigenome dynamics: a quantitative genetics perspective. *Nat. Rev. Genet.* 9, 883–890. doi: 10.1038/nrg2467

Johannes, F., Porcher, E., Teixeira, F. K., Saliba-Colombani, V., Simon, M., Agier, N., et al. (2009). Assessing the impact of transgenerational epigenetic variation on complex traits. *PLoS Genet.* 5:e1000530. doi: 10.1371/journal.pgen.1000530

Karimi, M., Johansson, S., and Ekstrom, T. J. (2006a). Using LUMA A luminometric-based assay for global DNA-methylation. *Genetics* 1, 45–48.

Karimi, M., Johansson, S., Stach, D., Corcoran, M., Grander, D., Schalling, M., et al. (2006b). LUMA (LUMinometric methylation assay) - a high throughput method to the analysis of genomic DNA methylation. *Exp. Cell Res.* 312, 1989–1995. doi: 10.1016/j.yexcr.2006.03.006

Kitabayashi, I., Aikawa, Y., Nguyen, L. A., Yokoyama, A., and Ohki, M. (2001). Activation of AML1-mediated transcription by MOZ and inhibition by the MOZ-CBP fusion protein. *EMBO J.* 20, 7184–7196. doi: 10.1093/emboj/20.24.7184

Knott, S. A., Marklund, L., Haley, C. S., Andersson, K., Davies, W., Dellegren, H., et al. (1998). Multiple marker mapping of quantitative trait loci in a cross between outbred wild boar and large white pigs. *Genetics* 149, 1069–1080.

Kooke, R., Johannes, F., Wardenaar, R., Becker, F., Etcheverry, M., Colot, V., et al. (2015). Epigenetic basis of morphological variation and phenotypic plasticity in *Arabidopsis thaliana*. *Plant Cell* 27, 337–348. doi: 10.1105/tpc.114.133025

Kooke, R., and Keurentjes, J. J. B. (2015). Epigenetic variation contributes to environmental adaptation of *Arabidopsis thaliana*. *Plant Signal. Behav.* 10:e1057368. doi: 10.1080/15592324.2015.1057368

Kouzarides, T. (2007). Chromatin modifications and their function. *Cell* 128, 693–705. doi: 10.1016/j.cell.2007.02.005

Krzyżewska, I. M., Maas, S. M., Henneman, P., Lip, K. V. D., Venema, A., Baranano, K., et al. (2019). A genome-wide DNA methylation signature for SETD1B-related syndrome. *Clin. Epigenet.* 11:156. doi: 10.1186/s13148-019-0749-743

Latzel, V., Allan, E., Silveira, A. B., Colot, V., Fischer, M., and Bossdorf, O. (2013). Epigenetic diversity increases the productivity and stability of plant populations. *Nat. Commun.* 4:2875. doi: 10.1038/ncomms3875

Lea, A. J., Altmann, J., Alberts, S. C., and Tung, J. (2016). Resource base influences genome-wide DNA methylation levels in wild baboons (*Papio cynocephalus*). *Mol. Ecol.* 25, 1681–1696. doi: 10.1111/mec.13436

Legarra, A., and Fernando, R. L. (2009). Linear models for joint association and linkage QTL mapping. *Genet. Sel. Evol.* 41:43. doi: 10.1186/1297-9686-41-43

Luo, Z. W. (1998). Detecting linkage disequilibrium between a polymorphic marker locus and a trait locus in natural populations. *Heredity* 80, 198–208. doi: 10.1046/j.1365-2540.1998.00275.x

Lynch, M., and Walsh, B. (1998). *Genetics and Analysis of Quantitative Traits*. Sunderland, MA: Sinauer Associates.

Massicotte, R., Whitelaw, E., and Angers, B. (2011). DNA methylation: a source of random variation in natural populations. *Epigenetics* 6, 421–427. doi: 10.4161/epi.6.4.14532

McRae, A. F., Powell, J. E., Henders, A. K., Bowdler, L., Hemani, G., Shah, S., et al. (2014). Contribution of genetic variation to transgenerational inheritance of DNA methylation. *Genome Biol.* 15:R73. doi: 10.1186/gb-2014-15-5-r73

Nanda, J. S., Kumar, R., and Raghava, G. P. S. (2016). DbEM: a database of epigenetic modifiers curated from cancerous and normal genomes. *Sci. Rep.* 6:19340. doi: 10.1038/srep19340

Ou, M., Mao, H. H., Luo, Q., Zhao, J., Liu, H. Y., Zhu, X. P., et al. (2019). The DNA methylation level is associated with the superior growth of the hybrid fry in snakehead fish (*Channa argus* x *Channa maculata*). *Gene* 703, 125–133. doi: 10.1016/j.gene.2019.03.072
Piepho, H. P. (2001). A quick method for computing approximate thresholds for quantitative trait loci detection. *Genetics* 157, 425-432.

Ricordeau, G., Tchamitchian, L., Brunel, J., Nguyen, T., and François, D. (1992). La race ovine INRA 401: un exemple de souche synthétique. *INRA Produ. Anim.* 5, 255-262. doi: 10.20870/productions-animales.1992.5.hs.4300

Roth, T. L. (2013). Epigenetic mechanisms in the development of behavior: advances, challenges, and future promises of a new field. *Dev. Psychopathol.* 25(4 Pt 2), 1279-1291. doi: 10.1017/S0954579413000618

Roussot, O., Feve, K., Plisson-Petit, F., Pitel, F., Faure, J. M., Beaumont, C., et al. (2003). AFLP linkage map of the Japanese quail *Coturnix japonica*. *Genet. Sel. Evol.* 35, 559-572. doi: 10.1186/1297-9686-35-6-559

SAS Institute Inc (1999). *SAS/STAT® User's Guide, Version 8*. Cary, NC: SAS Institute Inc.

Schlichting, C. D., and Wund, M. A. (2014). Phenotypic plasticity and epigenetic marking: an assessment of evidence for genetic accommodation. *Evolution* 68, 656-672. doi: 10.1111/evo.12348

Seebacher, F., and Krause, J. (2019). Epigenetics of social behaviour. *Trends Ecol. Evol.* 34, 818-830. doi: 10.1016/j.tree.2019.04.017

Shah, S., McRae, A. F., Marioni, R. E., Harris, S. E., Gibson, J., Henders, A. K., et al. (2014). Genetic and environmental exposures constrain epigenetic drift over the human life course. *Genome Res.* 24, 1725-1733. doi: 10.1101/gr.176933.114

Sharma, M., Zarnegar, M., Li, X. Y., Lim, B., and Sun, Z. J. (2000). Androgen receptor interacts with a novel MYST protein, HBO1. *J. Biol. Chem.* 275, 35200-35208. doi: 10.1074/jbc.M004838200

Singmann, P., Shem-Tov, D., Wahl, S., Grallert, H., Fiorito, G., Shin, S.-Y., et al. (2015). Characterization of whole-genome autosomal differences of DNA methylation between men and women. *Epigenetics Chromatin* 8:43. doi: 10.1186/s13072-015-0035-33

Soriano-Tarraga, C., Jimenez-Conde, J., Giralt-Steinhauer, E., Ois, A., Rodriguez-Campello, A., Cuadrado-Godia, E., et al. (2013). DNA isolation method is a source of global DNA methylation variability measured with LUMA. experimental analysis and a systematic review. *PLoS One* 8:e60750. doi: 10.1371/journal.pone.0060750

Turck, F., and Coupland, G. (2014). Natural variation in epigenetic gene regulation and its effects on plant developmental traits. *Evolution* 68, 620-631. doi: 10.1111/evo.12286

Vaissière, T., Sawan, C., and Herceg, Z. (2008). Epigenetic interplay between histone modifications and DNA methylation in gene silencing. *Mutat. Res. Rev. Mut. Res.* 659, 40-48. doi: 10.1016/j.mrrev.2008.02.004

van Dongen, J., Nivard, M. G., Willemse, G., Hottenga, J. J., Helmer, Q., Dolan, C. V., et al. (2016). Genetic and environmental influences interact with age and

sex in shaping the human methylome. *Nat. Commun.* 7:11115. doi: 10.1038/ncomms11115

Verhulst, E. C., Mateman, A. C., Zwier, M. V., Caro, S. P., Verhoeven, K. J. F., and van Oers, K. (2016). Evidence from pyrosequencing indicates that natural variation in animal personality is associated with DRD4 DNA methylation. *Mol. Ecol.* 25, 1801-1811. doi: 10.1111/mec.13519

Vija, L., Meduri, G., Comperat, E., Vasilii, V., Izard, V., Ferlicot, S., et al. (2013). Expression and characterization of androgen receptor coregulators, SRC-2 and HBO1, during human testis ontogenesis and in androgen signaling deficient patients. *Mol. Cell. Endocrinol.* 375, 140-148. doi: 10.1016/j.mce.2013.05.004

Wagner, M., Koslowski, M., Paret, C., Schmidt, M., Tureci, O., and Sahin, U. (2013). NCOA3 is a selective co-activator of estrogen receptor alpha-mediated transactivation of PLAC1 in MCF-7 breast cancer cells. *BMC Cancer* 13:570. doi: 10.1186/1471-2407-13-570

Wang, H. J., Pochampalli, M., Wang, L. Y., Zou, J. X., Li, P. S., Hsu, S. C., et al. (2019). KDM8/JMJD5 as a dual coactivator of AR and PKM2 integrates AR/EZH2 network and tumor metabolism in CRPC. *Oncogene* 38, 17-32. doi: 10.1038/s41388-018-0414-x

Yousefi, P., Huen, K., Davé, V., Barcellos, L., Eskenazi, B., and Holland, N. (2015). Sex differences in DNA methylation assessed by 450 K BeadChip in newborns. *BMC Genom.* 16:911. doi: 10.1186/s12864-015-2034-y

Zhang, B. L., Gu, X. H., Han, X., Gao, Q., Liu, J., Guo, T. W., et al. (2020). Crosstalk between DNA methylation and histone acetylation triggers GDNF high transcription in glioblastoma cells. *Clin. Epigenet.* 12:47. doi: 10.1186/s13148-020-00835-833

Zhang, Y. L., Li, F. Z., Feng, X., Yang, H., Zhu, A. X., Pang, J., et al. (2017). Genome-wide analysis of DNA Methylation profiles on sheep ovaries associated with prolificacy using whole-genome Bisulfite sequencing. *BMC Genom.* 18:759. doi: 10.1186/s12864-017-4068-4069

**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.

Copyright © 2020 Hazard, Plisson-Petit, Moreno-Romieux, Fabre and Drouilhet. 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.