Population differences in Chinook salmon (Oncorhynchus tshawytscha) DNA methylation: Genetic drift and environmental factors

Abstract Local adaptation and phenotypic differences among populations have been reported in many species, though most studies focus on either neutral or adaptive genetic differentiation. With the discovery of DNA methylation, questions have arisen about its contribution to individual variation in and among natural populations. Previous studies have identified differences in methylation among populations of organisms, although most to date have been in plants and model animal species. Here we obtained eyed eggs from eight populations of Chinook salmon (Oncorhynchus tshawytscha) and assayed DNA methylation at 23 genes involved in development, immune function, stress response, and metabolism using a gene‐targeted PCR‐based assay for next‐generation sequencing. Evidence for population differences in methylation was found at eight out of 23 gene loci after controlling for developmental timing in each individual. However, we found no correlation between freshwater environmental parameters and methylation variation among populations at those eight genes. A weak correlation was identified between pairwise DNA methylation dissimilarity among populations and pairwise F ST based on 15 microsatellite loci, indicating weak effects of genetic drift or geographic distance on methylation. The weak correlation was primarily driven by two genes, GTIIBS and Nkef. However, single‐gene Mantel tests comparing methylation and pairwise F ST were not significant after Bonferroni correction. Thus, population differences in DNA methylation are more likely related to unmeasured oceanic environmental conditions, local adaptation, and/or genetic drift. DNA methylation is an additional mechanism that contributes to among population variation, with potential influences on organism phenotype, adaptive potential, and population resilience.


| INTRODUC TI ON
Local adaptation occurs when organisms evolve in response to selective pressures in their immediate environment, resulting in increased individual fitness within their native habitat relative to non-native habitats (García de Leániz et al., 2007;Kawecki & Ebert, 2004;Savolainen et al., 2013). Traditionally, the main mechanism considered to be underlying local adaptation has been genetic adaptation: selection acts upon the phenotypes produced by standing genetic variation, resulting in increased frequency of beneficial alleles and thus evolution of populations over multiple generations (Bernatchez, 2016). Additional mechanisms are now also accepted as contributing to local adaptation. Chromosomal translocations can result in co-adapted gene complexes resistant to crossing-over (Barth et al., 2019;Kess et al., 2020;Kirkpatrick & Barton, 2006;Lehnert et al., 2019). Maternal effects can influence offspring phenotype to prepare offspring for a predicted environment based on the experiences of the mother Galloway, 2005;Galloway & Etterson, 2007). Phenotypic plasticity, whereby an organismal phenotype, is shifted toward an "ideal" phenotype based on the environment without underlying genetic changes (Hutchings, 2011;Pfennig et al., 2010;Torres-Dowdall et al., 2012), which can be mediated by differences in gene expression (Fangue et al., 2006;Whitehead & Crawford, 2006;Wellband & Heath, 2013). These mechanisms can all occur simultaneously in an organism, leading to a wide variety of mechanistic contributions toward local adaptation.
Adaptive population differences in gene expression have been reported in a broad variety of taxa. Differences in gene expression occur among populations of killifish (Fundulus heteroclitus) across a natural thermal cline (Fangue et al., 2006), among rainbow trout (O. mykiss) from different tributaries subjected to stress challenges (Wellband & Heath, 2013), between populations of the copepod Tigriopus californicus residing in different thermal regimes (Schoville et al., 2012), among populations of Drosophila subobscura across latitudinal and thermal clines in Europe (Porcelli et al., 2016), and both within and among populations of teleost fish from the genus Fundulus (Oleksiak et al., 2002). Further, patterns in gene expression variation may also reflect parallel evolution due to similar environmental conditions (reviewed in Fraser et al., 2011). While local adaptation through gene expression variation has been frequently reported, the mechanisms underlying these differences in gene expression are poorly characterized, though environmental, genetic, and epigenetic variation could contribute to locally adapted gene expression profiles.
DNA methylation is one potential mechanism underlying transcriptional differences observed among populations in the context of local adaptation. DNA methylation is the addition of a methyl group to a cytosine (C) base that precedes a guanine (G) in the DNA sequence, known as a CpG site (Head, 2014). Numerous studies have shown that DNA methylation is highly sensitive to environmental signals (Barfield et al., 2014;Bossdorf et al., 2008;Foust et al., 2016;Herrera & Bazaga, 2010;Richards et al., 2010) and is involved in acclimation to environmental stress (Metzger & Schulte, 2017; Morán et al., 2013). Due to the potential to modify methylation in response to environmental cues, methylation presents an important mechanistic intersection between acclimation and adaptation, particularly with extensive evidence for rapid (or "contemporary") evolution over short time scales (Stockwell et al., 2003). Methylation has been shown to be a highly targeted process (Venney et al., 2016(Venney et al., , 2020. Therefore, short-term changes in methylation can occur that allow an organism to cope with its environment, without the lag times associated with selection on standing genetic variation (Bossdorf et al., 2008;Hu & Barrett, 2017;Richards et al., 2010), consistent with rapid evolution. Due to the sensitivity of methylation to environmental cues, it presents a novel mechanism for organisms to adapt to their environment and adds an additional level of complexity in organismal phenotypic variation and evolution (Bossdorf et al., 2008). Furthermore, methylation may respond to environmental stress, allowing for targeted short-term responses to environmental changes, which cannot occur through genetic adaptation (Hu & Barrett, 2017). If methylation results in phenotypic plasticity, it may act in lieu of genetic adaptation, since the detrimental phenotype is no longer present to be selected against, or it may prolong the persistence of organisms in stressful environments until selection and genetic adaptation can occur (Crispo, 2008).
Population-level variation in methylation has been reported in a variety of species and appears to have an underlying genetic basis.
While methylation can be taxon-specific, studies in several taxa have identified a link between genetic and epigenetic variation (Fraser et al., 2012;Herrera & Bazaga, 2010;Liu et al., 2012). For example, a study in Spanish violets (Viola cazorlensis) across an elevation gradient identified a strong correlation between methylation and genetic variation using pairwise distance-based AFLP analyses (Herrera & Bazaga, 2010). Similar results were found using restriction enzymebased methods for whole-genome DNA methylation estimation and sequence polymorphism in female great roundleaf bat (Hipposideros armiger) populations (Liu et al., 2012), when comparing CpG-specific methylation and sequence variation in oak (Quercus lobata Née) populations (Platt et al., 2015), and for correlations between methylation differences and allele frequencies among human ethnicities (Fraser et al., 2012). However, a study in salt marsh perennials (Spartina alterniflora) was unable to link genetic differences with variation in methylation through AFLP-based approaches and instead found a strong correlation with environmental variation (Foust et al., 2016).
Thus, the relationship among epigenetic variation, genetic variation, and environmental heterogeneity is unclear, yet characterizing the interactions between these three drivers of population-level phenotypic variation is important in determining the role DNA methylation may play in driving local adaptation and thus population and species resiliency. While many studies have shown methylation differences among populations, most studies have focused on agriculturally important laboratory-reared species, while studies of natural populations are limited (Richards et al., 2010), making the role of DNA methylation in population differentiation unclear.
Chinook salmon (Oncorhynchus tshawytscha) are a culturally, ecologically, and economically important species of Pacific salmon.
There is ample evidence for local adaptation based on functional differences among populations of Chinook salmon resulting in increased fitness in their native environments (Fraser et al., 2011).
Adaptive genetic variation occurs at selected immune and growthrelated candidate loci indicating genetic adaptation to their environment, while divergence at neutral (microsatellite) loci is related to isolation and genetic drift (Heath et al., 2006). Adaptation can occur within Chinook salmon stocks, for example, as evidenced by intrapopulation genetic differences in circadian clock genes based on migration timing, in the absence of neutral genetic variation (O'Malley et al., 2013). Variants impacting life-history traits associated with environmental differences have also been reported in recently colonized Chinook salmon populations (Unwin et al., 2000), as well as differences in genetic variance components and fitnessrelated traits . Thus, there is abundant evidence for adaptive differences among populations of Chinook salmon, though most studies focus on genetic differences. While there have been studies documenting neutral and functional genetic variation among populations of Chinook salmon, it is unclear how rapid adaptation occurs when local conditions change or salmon colonize new habitats. However, studies have shown evidence for rapid adaptation to hatchery rearing, resulting in differences in gene expression (Christie et al., 2016), reproductive success (Christie et al., 2012), and DNA methylation (Gavery et al., 2018;Le Luyer et al., 2017). Due to the role of DNA methylation in rapid evolution of salmonids, it is possible that DNA methylation is important for responding to environmental changes, as well as maintaining standing genetic variation in salmon.
The goal of this study is to determine the potential role of DNA methylation in maintaining differences (adaptive or drift-related) among populations and to assess genetic and environmental drivers of population-level differences in methylation. We characterize locus-specific population differences in DNA methylation in Chinook salmon while testing for the effects of developmental timing and sampling year. We associate population differences in methylation with freshwater environmental differences and genetic drift on levels of methylation at selected genes. We obtained eyed eggs from eight populations of Chinook salmon and measured DNA methylation using a gene-targeted PCR-based DNA methylation assay for next-generation sequencing. We expected that populations would exhibit different levels of DNA methylation at specific functional loci. Such patterns of methylation differences among populations could be due to environmental acclimation (Foust et al., 2016), underlying adaptive genetic variation (Fraser et al., 2012;Herrera & Bazaga, 2010;Liu et al., 2012), or maternal effects at the eyed egg stage (Venney et al., 2020). These differences would likely show developmental and interannual variation due to strict developmental control of DNA methylation and interannual environmental variation. We hypothesized that population differences in methylation should occur in specific genes in response to unique environmental conditions and/or selective pressures among natural environments.
We tested for correlations between locus-specific methylation and freshwater environmental variables from the native rivers of each population to determine whether local environmental factors influence gene-specific DNA methylation differences. We also tested for a correlation between genetic drift (variation at neutral marker loci) and methylation differences among populations to determine whether methylation differences could be explained by population divergence due to genetic drift (and/or geographic isolation) distance. DNA methylation presents an additional, potentially powerful and rapid mechanism through which populations can respond to their environments, cope with environmental stress, and evolve. Embryos were dissected from 48 eyed eggs per population (n = 10). The yolk was removed and the embryos were digested whole in 10 µl of 20 mg/ml proteinase K and 1,000 µl of digestion buffer (100 mM NaCl, 50 mM Tris-HCl pH 8.0, 10 mM EDTA, 0.5% SDS) at 37°C for 24 hr. We used 150 µl of the digested product for DNA extraction in a high-throughput automated plate-based DNA extraction protocol (Venney et al., 2020). Methylation analysis was performed with 21 published bisulfite sequencing primers (Venney et al., 2016) and two novel bisulfite sequencing primer sets for growth hormone 2 (GH2, forward primer 5′-TTATTAAACCTTTCTAAAAACACAC-3′, reverse primer 5′-ATTTAAATTTTAATTTTTTATAGGG-3′, 241 bp fragment excluding primer sequences) and heat shock factor 1b (hsf1b, forward primer 5′-AGGATTAGGATTTTGAAGAGGATTT-3′, reverse primer 5′-AATTAATTTTTCATCATCTACACATTAACA-3′, 132 bp fragment excluding primer sequences). All primers were designed for gene regions with little to no sequence variation to minimize the effects of genetic variation on the interpretation of DNA methylation data. Assayed genes were selected for their roles in early development, stress and immune function, metabolism, early growth, and differentiation. Amplicons ranged from 79 to 249 bp, with a total of 4,111 bp sequenced excluding primer sequences (Table 1).

| Bisulfite conversion and sequencing
PCRs were performed using a two-stage PCR approach (Venney et al., 2016) where the first stage amplified the targeted gene region, and the second stage ligated sample barcode and adaptor sequences to the amplicon. Barcode sequences are 10-12 bp unique sequences that allow for the identification of individual samples in massively parallel (next-generation) sequencing. Samples were split among three sequencing runs and sequenced with an Ion 318™ Chip using an Ion PGM™ Sequencing 400 bp kit on the Ion Torrent Personal Genome Machine® (PGM™) with an expected 500 reads per gene per sample.

| Bisulfite sequencing data processing
Sequence data files were demultiplexed using mothur (Schloss et al., 2009), primer sequences were trimmed, and one fastq sequence file was created per individual. Bisulfite sequence data were aligned with known genomic sequences using bwa-meth (Pedersen et al., 2014) with a maximum of two mismatches per sequence to ensure sequences represented the target genes. Tabulated methylation data from bwameth were imported into R (R Development Core Team, 2021) for quality filtering to ensure the same CpG sites were compared across all samples: CpG sites were excluded from the analysis if they were sequenced (a) with fewer than five reads per gene per sample (though generally greatly exceeded this requirement) and (b) in less than 70% of individuals. Rosner's test for extreme outliers was used to exclude significant outlier data points, which were likely reflections of low sequence depth rather than biologically meaningful variation. The final processed data provided average percent methylation for each indi-    These datasets were analyzed using custom R scripts (R Development Core Team, 2021). In brief, datasets were loaded into R using adegenet v.2.1.1 (Jombart, 2008), and dendrograms were constructed using the aboot function of poppr v.2.8.3 (Kamvar et al., 2014) with the edwards.
Pairwise distance matrices for microsatellite and SNP data were compared to methylation matrices to determine whether populationlevel differences in methylation corresponded with expected divergence due to isolation and genetic drift. A Euclidean distance matrix for population-level methylation variation was generated in the R package ade4 (Dray & Dufour, 2007)

| Population differences in methylation
Population and the population by gene interaction significantly affected methylation levels across all genes combined (both p < .001, R 2 = 0.10), indicating that while populations differ in overall methylation levels, they also differ in extents of gene-specific methylation.
Direct between-gene differences in methylation were not quantifiable, as gene methylation values were standardized and centered around zero by using the ATU model residuals (p = 1.0).

| Principal component regressions for environmental effects on methylation
Six principal components explained 98.9% of variation in the environmental dataset. These PCs were retained in the analysis based on PC eigenvalue >1 and the relative values of PCs in the Scree plot (

| D ISCUSS I ON
DNA methylation presents an evolutionary mechanism for individuals to rapidly respond to environmental changes and improve their survival in natural systems. In contrast, novel beneficial genetic mutations and natural selection acting upon existing variation can be slow processes that take place over generations (Bossdorf et al., 2008 TA B L E 4 p-values and R 2 values from ANOVAs and Mantel tests for population effects on DNA methylation in Chinook salmon have primarily focused on sources of individual variation, rather than population-level differences in methylation (Hu & Barrett, 2017), yet population-level differences in methylation could explain heritable variation among populations which cannot be explained solely by genetic variation (Bossdorf et al., 2008). We observed significant population differences in methylation across all genes combined, as well as a significant population by gene interaction, indicating that populations differ in overall methylation, and also that the extent of methylation differences between populations varies among individual genes.
Methylation differences among populations have been reported in several other studies (Barfield et al., 2014;Foust et al., 2016;Fraser et al., 2012;Herrera & Bazaga, 2010;Liu et al., 2012;Platt et al., 2015;Richards et al., 2010) with the potential to contribute to rapid acclimation and/or adaptation to stressors (Bossdorf et al., 2008;Hu & Barrett, 2017;Richards et al., 2010). The population-level differences in methylation we report may represent an important evolutionary mechanism that could contribute to the extensive adaptive variation observed in natural populations of Chinook salmon (Fraser et al., 2011). However, the patterns of considerable population-level variation in DNA methylation reported here are of broad interest when considering potential mechanisms of phenotypic differentiation in natural populations in general. Population-level differences in methylation could reflect shortterm acclimation to the local environment, or local adaptation due to environmental selection on heritable phenotypes. While several studies have identified population differences in methylation, most focus on methylation at the genome-wide or whole-genome level rather than using a candidate gene approach. We observed population-level differences in methylation at specific genes in Chinook salmon eyed eggs: four heat shock protein genes (hsc71, hsp47, hsp70, and hsp90), three immune genes (p53, RAG1, and or population-level differences in heat shock protein expression (Sørensen et al., 2003;Tine et al., 2010). Previous studies in teleost fish have identified differences in immune response among populations (Evans et al., 1997(Evans et al., , 2010Fraser et al., 2011), as well as differences in hormone concentrations and endocrine function (Carr & Patiño, 2011;Sopinka et al., 2017). Differences in gene methylation could reflect acclimation or adaptation to local environments, though further research is required to determine whether population-level differences in gene-specific methylation result from acclimation or adaptation. However, significant differences in methylation between BQ and Harr with no significant temporal effects suggest local adaptation. Future research measuring methylation in reciprocal transplants or in common garden experiments with natural populations could determine whether population-level variation in methylation is retained and therefore whether it likely represents acclimation or adaptation. Regardless of the underlying process, the genes showing significant population effects are logical targets for differential DNA methylation due to differences in environmental context and stressors among populations.
DNA methylation is often influenced by environmental context (Barfield et al., 2014;Bossdorf et al., 2008;Foust et al., 2016;Herrera & Bazaga, 2010;Richards et al., 2010). We used principal component analysis and regression to test for environmental effects on DNA methylation among populations using environmental data from the natal streams of the studied Chinook salmon populations. We found no significant effects after correcting for multiple comparisons, which was unexpected, as many studies have reported environmental effects on methylation (Angers et al., 2010;Dimond & Roberts, 2016;Foust et al., 2016;Le Luyer et al., 2017;Morán et al., 2013). The lack of significant environmental correlates is likely due to our use of Chinook salmon eggs. At the egg stage, the embryo is isolated and protected from the environment, which may reduce its response to environmental variation, though it is still possible that eggs respond to local environmental conditions through changes in methylation.
If methylation in eggs is genetically encoded rather than environmentally determined, these differences in methylation may reflect genetic adaptation rather than acclimation. Alternatively, Chinook salmon exhibit strong maternal effects on DNA methylation at the eyed egg stage (Venney et al., 2020) which may increase variation within a population and reduce correlations between gene-specific DNA methylation and environmental variables. Parents experience the freshwater environment prior to spawning and therefore could alter egg methylation signals in response to the offspring's predicted environment. Thus, the population-level differences in methylation observed in Chinook salmon may be due to acclimation or adaptation to freshwater environmental signals from their parents, eyed egg acclimation to the environment, or due to genetic differences among populations (Fraser et al., 2012;Liu et al., 2012).
Population epigenetic studies in different taxa vary in their conclusions as to the link between epigenetic differences among populations and genetic divergence. A study in salt marsh perennial plants found no link between genetic and epigenetic differences across environmental gradients, but a strong correlation with environmental conditions (Foust et al., 2016). However, another study linked DNA methylation differences in Spanish violets to genetic differences identified by AFLP in response to elevation (Herrera & Bazaga, 2010).
A significant correlation between genetic and epigenetic variation was also reported among female great roundleaf bat populations (Liu et al., 2012) and due to differences in allele frequency among human ethnic groups (Fraser et al., 2012). Here we compared epigenetic differences among populations to neutral genetic variation at microsatellite loci to determine whether differences in DNA methylation among populations align with genetic drift. The correlation between microsatellite F ST and Euclidean pairwise dissimilarity in methylation among populations (p = .02, i 2 = 0.19) was likely primarily driven by higher correlations between epigenetic differences at GTIIBS and Nkef and neutral genetic divergence. However, there was no significant correlation between SNP divergence and methylation pairwise dissimilarity across all eight genes (p = .12, R 2 = 0.064), likely due to weaker single-gene correlations between GTIIBS and Nkef methylation and SNP divergence. While divergence in methylation among populations may be attributed in part to genetic drift, neutral genetic divergence in Chinook salmon is affected by geographic distance (Beacham et al., 2006;Heath et al., 2006). Given that geographic distance is expected to be related to ecosystem dissimilarity, it is possible that weak signals of drift may simply reflect environmental similarities among proximate populations. The weak correlation between neutral genetic markers and differences in methylation among populations suggests that while drift acts on methylation, mechanisms other than drift (such as selective mechanisms) likely also contribute to differences in methylation among populations. It is also possible that population differences in methylation are due to genetic control of methylation processes, that is, different genotypes result in different methylation patterns (Liu et al., 2012). We show that differences in methylation among populations are not well explained by genetic drift alone, suggesting that methylation is likely also influenced by a combination of genomic differences among populations, environmental acclimation, and local adaptation.
We found that ATUs (a measure of developmental timing in salmon), and the interaction between population and sampling year influenced DNA methylation. DNA methylation patterns have been shown to change through development in fish (Fang et al., 2013;Fellous et al., 2018;Venney et al., 2020), and thus, we expected differences in methylation levels in the eyed eggs as they developed.
In mangrove rivulus (Kryptolebias marmoratus), changes in methylation occurred during development throughout organogenesis leading up to hatch (Fellous et al., 2018). However, while developmental changes in methylation are well-characterized, interannual changes in methylation are not. We found a significant population by sampling year interaction on one gene after correcting for multiple comparisons when controlling for ATU in Harr and BQ 2015 and 2017 samples. The significant population by year effect suggests that there is some interannual variation in methylation within populations which is likely due to acclimation, though population-level differences persist across years. These differences could be due to changes in freshwater and marine environments experienced by the parents and offspring from year to year. This raises the question of whether the egg's freshwater environment, or the parental marine and/or freshwater environments are influencing offspring methylation patterns. Since the population of origin (Harr vs. BQ) significantly affected methylation of four genes after FDR correction, and sampling year explained very little phenotypic variation in methylation (Table 3), population clearly has a greater effect on methylation state than sampling year. Our results reinforce the importance of controlling for potential confounding variables such as organism age/developmental stage and year of sampling, since methylation is a highly sensitive and dynamic mechanism for controlling gene expression.
Population epigenetic status is an important new consideration in evolutionary and ecological studies (Bossdorf et al., 2008) since DNA methylation could act as a highly dynamic evolutionary mechanism upon which selection could act (Bossdorf et al., 2008;Hu & Barrett, 2017). Unlike genetic adaptation, which requires standing variation and selection, methylation changes are rapid and dynamic, adding an additional layer of complexity and specificity for organisms to acclimate and adapt to their environment (Bossdorf et al., 2008;Hu & Barrett, 2017). In this study, we provide evidence for differences in methylation among populations of Chinook salmon, consistent with previous population epigenetic studies in other taxa (Barfield et al., 2014;Foust et al., 2016;Fraser et al., 2012;Herrera & Bazaga, 2010;Liu et al., 2012;Platt et al., 2015;Richards et al., 2010). Despite reported strong environmental effects on DNA methylation (Barfield et al., 2014;Bossdorf et al., 2008;Foust et al., 2016;Herrera & Bazaga, 2010;Richards et al., 2010), we found no link between freshwater environmental parameters and population differences in methylation. This may be due to (a) methylation corresponding to the marine environment experienced by the parents rather than freshwater variables considered here; (b) strong maternal effects on methylation at the eyed egg stage in Chinook salmon (Venney et al., 2020), which could decrease DNA methylation-environment correlations due to varying environmental experiences of individual mothers; or (c) key environmental variables that affect methylation but were not included in our analysis. We identified weak correlations between genetic drift and DNA methylation, indicating that while some changes in methylation state among populations are likely due to drift, other differences could be the result of selection (Bossdorf et al., 2008) or are linked to underlying functional genetic differences (Fraser et al., 2012).
Characterizing sources of phenotypic variation among natural populations are critical to understanding individual variation and the adaptive potential and resiliency of natural populations. DNA methylation is an important source of phenotypic variation, and a substrate for characterizing adaptive response in nature since an organism's environment and experiences can influence methylation levels (Bossdorf et al., 2008;Burggren, 2014). Furthermore, methylation signals can be passed on to offspring generations and beyond (Kamstra et al., 2018;Santangeli et al., 2019), resulting in rapid adaptation and evolutionary change in response to changing environments.