Genome-wide Methylation Patterns Under Caloric Restriction in Daphnia magna

The degradation of epigenetic control with age is associated with progressive diseases of ageing, including cancers, immunodeficiency and diabetes. Reduced caloric intake slows the effects of aging and age-related diseases, a process likely to be mediated by the impact of caloric restriction on epigenetic factors such as DNA methylation. We used whole genome bisulphite sequencing to study how DNA methylation patterns change with diet in a small invertebrate, the crustacean Daphnia magna. Daphnia show the classic response of longer life under CR, and they reproduce clonally, which permits the study of epigenetic changes in the absence of genetic variation. Global CpG methylation was 0.7-0.9%, and there was no difference in overall methylation levels between normal and calorie restricted replicates. However, 453 regions were differentially methylated (DMRs) between the normally fed and calorie restricted (CR) replicates. Of these 61% were hypomethylated in the CR group, and 39% were hypermethylated in the CR group. Gene Ontogeny (GO) term enrichment of hyper and hypo-methylated genes showed significant over- and under-representation in three molecular function terms and four biological process GO terms. Notable among these were kinase and phosphorylation activity, which have a well-known functional link to cancers.


Introduction 50
Epigenetic modifications play a key role in maintaining gene expression and organismal development. This is particularly evident when epigenetic control degrades, resulting in 52 progressive diseases in humans, including cancers, immunodeficiency and diabetes [1]. The degradation of epigenetic control with age is proposed to occur in a drift-like process. One 54 mechanism that may rescue age-related epigenetic dysregulation is caloric restriction (CR): reduced caloric intake without malnutrition or loss of nutrients. CR slows the effects of 56 aging and postpones the development of age-related diseases [2][3][4][5]. In rhesus monkeys and mice, CR of 30% and 40% respectively appears to reduce epigenetic drift in methylation and 58 increases lifespan, which in rodents can be an extension of up to 50% [6]. Similar results have been seen in yeast, spiders, worms, fish and non-human primates [4,7,8]. CR may also 60 delay a spectrum of diseases such as cancer, kidney disease, autoimmune disease and diabetes [9][10][11], as well as neurodegenerative diseases [12,13]. 62 DNA Methylation, a reversible covalent modification that regulates gene expression, is the best-studied epigenetic mechanism. DNA Methylation of cytosines occurs when DNA 64 methyltransferase enzymes (DNMTs) transfer a methyl group onto cytosine [14] to create a 5-methylcytosine. Most commonly at a cytosine immediately followed by guanine (CpG 66 site). There are three DNMT enzymes: DNMT3 establishes methylation de novo, DNMT1 maintains methylation, and DNMT2, which has no known role in DNA methylation. A 68 reduction in expression levels of DNMT enzymes is associated with ageing, leading to a global loss of genomic methylation [15]. In mammals around 70% of CpGs are methylated 70 [16], however in invertebrates the rate in species sampled to date is lower, from 0% in flies to 15% in the oyster Crassostrea gigas [17,18]. The model crustaceans Daphnia magna and 72 Daphnia pulex (Arthropoda: Crustacea) have genomic CpG methylation of 0.52% and 0.7% respectively [19]. 74 4 CpG methylation can increase or decrease gene expression dependent on the location of the methylation. In promoter regions, which can be rich in CpGs and are known as CpG 76 islands, it represses expression of the gene. Further to this, many CpG islands are also enriched for permissive chromatin modification, which condenses the structure of 78 chromatin and further prevents transcription. In contrast, methylation of gene bodies leads to an increase in expression of the effected gene. Invertebrates have few CpG Islands, and 80 methylation predominantly occurs in gene bodies, and is enriched in exonic sequence [20][21][22][23][24][25], where it may enhance transcription or mediate alternative splicing [25][26][27]. 82 Interestingly, in silkworms there is no correlation between methylation in promoters and gene expression [28], suggesting invertebrates and vertebrates differ in their usage of CpG 84 methylation.
The relationship between diet and CpG methylation, and subsequent impact on ageing and 86 health, is well established. Indeed, DNA methylation may be a predictor of biological age. CR in mice caused a two-year difference in biological (0.8) versus chronological (2.8) age, 88 while in rhesus monkeys CR resulted in a biological age of 20 years for monkeys of chronologically aged 27 years [29]. Specific examples of a diet by methylation interaction 90 include the expression of DNMTs which have elevated expression in response to CR in cancer cells, which counteracts the global hypomethylation [30] observed during ageing. CR 92 also causes a reduction in lipid metabolism gene expression by DNA methylation of gene bodies in mouse livers [31]. As a result, older mice undergoing CR were protected from 94 fatty degeneration, visceral fat accumulation, and hepatic insulin resistance compared to controls. In rats and monkeys short-term CR in older individuals ameliorates the effects of 96 ageing with respect to disease markers, oxidative stress and damage, and increases the expression of longevity related genes [32,33]. The reverse is seen in obesity-like phenotypes 98 in rodent models. For example, in Agouti mice, the agouti viable yellow metastable epiallele 5 (A vy ) interacts with an upstream retrotransposon intracisternal A particle (IAP) [34]. 100 Unmethylated IAP results in yellow mice and negative health effects associated with obesity, whereas methylation at IAP results in brown healthy mice [35]. Waterland et al (2003)[36] 102 demonstrated the that supplementing mothers with folic acid, vitamin B12, choline and betaine shifted the offspring of obesity phenotype mice to smaller, brown mice indicative of 104 increased methylation at IAP.
Our work aims to determine if an experimentally controlled nutritional environment directs 106 changes in methylation status in a small invertebrate, the crustacean Daphnia magna. We do this by whole genome bisulphite sequencing of CR and normally-fed (NF) replicates, 108 identifying and characterizing regions of differential methylation. Daphnia show the classic response of longer life under CR and strong maternal effects; the offspring of calorie-110 restricted mothers being larger and more resistant to pathogens than their counterparts from better fed mothers. Provisioning of offspring, e.g. with carbohydrates, protein or fats, 112 is one explanation for these maternal-effect phenotypes, and epigenetic processes, such as methylation, are also potentially key regulators in these plastic responses to fluctuating 114 environments.
Daphnia have many attributes that make them favourable for epigenetic study. First, they 116 reproduce clonally, which permits the study of epigenetic changes in the absence of genetic variation. This also allows powerful study of genetic variation: clonal replicates are 118 equivalent to identical twin studies, but with an experimentally chosen number of -uplets.

Methylated sites prediction
Trimmed bisulphite-converted reads aligned to the genome using Bismark [40] exhibited 126 lower mapping efficiencies than standard short-read alignments typical of WGBS [41], with 20-32% or reads not aligning to the reference genome (read filtering and coverages per 128 replicate Table 1, Bismark report outputs Supplementary File 1). Of the aligned reads, 29-38% of reads were discarded as PCR duplicates, and 6-10% of the remainder contained 130 predicted CHH or CHG methylated sites which were also removed from analyses. This resulted in replicate average read coverages of 8-12-fold (read filtering and coverages per 132 replicate Table 1, Bismark report outputs Supplementary File 1). Global CpG methylation was 0.7-0.9% in all samples, and no difference in methylation levels is observed between 134 normal and calorie restricted replicates. Removal of polymorphic CpG sites, where a polymorphic C/T can be miscalled as methylated C, using variants predicted from the 136 bisulphite-unconverted data had little effect on the total number of sites (table 2). After filtering 99.2% of sites were retained per replicate on average per replicate, with average 138 total sites across replicates going from 6.9 million to 6.85 million. Hierarchical clustering of replicates by methylation status in methylKit [42] demonstrated that mother has a stronger 140 effect on global methylation status than nutritional treatment ( Figure 1).

GO term enrichment in methylated genes and DMRs
GO term enrichment was explored using the "weight01" and "classic" algorithms in topGO 156 [44] and molecular function (MF) and biological process (BP) terms are reported. The enrichment analysis (DMRs) using the weight01 algorithm showed significant over and under 158 representation in three molecular function (MF) GO terms and four biological process (BP) terms. The more permissive classic algorithm revealed significant enrichment in twenty-five 160 MF GO terms and twenty BP. These results have been condensed into their most-specific terms and direction of methylation in Table 3 (

Global methylation and differentially methylated regions 168 170
The global CpG methylome of ~0.7%, consistent across replicates, is possibly higher than for the previously sequenced D. magna CpG methylome of 0.5% [19], though this earlier 8 study was performed on a different strain and used different data filtering methods [19]. In line with previously observed methylation patterns in arthropods [20,23,24,47], the majority 174 (86%) of DMRs are found in gene bodies. Furthermore, most are present in exonic regions, although more so for hypo-(99%) than hypermethylated (80%) regions. This suggests DNA 176 methylation is regulating expression of targeted genes in Daphnia as for other invertebrates [20,47]. Although this requires confirmation by gene expression data, the expectation is that In what follows, we discuss the genes that are associated with the functional enrichment of differentially methylated regions (GO terms, Table 3), giving particular attention to genes 182 whose expression is known to respond to CR, or are linked to progressive disease of ageing or cancers. 184

Hypermethylation under CR
Perhaps the most prominent difference between CR and control Daphnia is the 186 hypermethylation of protein phosphorylation and tyrosine kinase activity (Table 3,  (Akt) [48], which also has roles in cell-cycle progression. In humans, aberrant expression of CAM-KK is a known factor in several cancers, and is considered a therapeutic target for 192 prostate and stomach cancers [49,50]. CAM-KK responses to CR are less well understood, but can protect against atherosclerosis by activation of AMP-activated protein kinase and AMPK, sirtuins, mTORC1 [54] and sirtuin triggered pathways.

200
The next major hypermethylated GO term, the molecular function term for ATP binding (GO:0005524), is associated with a gene group that includes kinases, ligases, and eleven  [55]. This suggests upregulation of energy production in our CR lines (though we also 206 find evidence of down-regulation of lipid metabolism, discussed below). Mutations in this gene in humans causes adrenomyeloneuropathy, characterised by an accumulation of 208 unbranched saturated fatty acids [56]. In ABCD1 knockout mice cholesterol levels are higher than in wildtypes, and are unaffected by cholesterol feeding [57]. Most pertinent to 210 CR, however, is the defective antioxidant response correlated with ABCD1 dysfunction [58], because reducing oxidative stress is a proposed mechanism by which CR increases 212 longevity [6,59,60]. GMP synthase expression gradually decreases with age resulting in lower cognitive performance [61]. The RNA helicase DDX39A has no connection to CR, 214 but its overexpression is associated with poor cancer prognosis [62][63][64][65][66].

232
The tRNA aminoacylation ligase (GO:0043039) group (Table 3) is associated with four tRNA ligases (for glutamate, proline, histidine, and phenylalanine). There is no research 234 directly linking these ligase genes to CR, but increased expression via methylation may be a response to low abundance of these amino acids, which indicates that future studies of diet 236 should vary protein availability. Indeed, it may be that protein restriction is more important than overall calorie restriction for longevity (cites). Additionally, fragments of tRNAs called 238 5' tRNA halves are a class of signalling molecules that are modulated by CR and ageing in mice [74,75], in which CR 'rescues' older 5' tRNA halves in line with other CR phenotypes. 240 An as yet undiscovered mechanism of diet and ageing could involve tRNA ligases regulating levels of 5' tRNA halves in response to CR. 242 Respiratory electron transport chain (GO:0022904; Table 3) contains two proteins: 244 cytochrome b-c1 complex subunit and NADH (Nicotinamide adenine dinucleotide reduced form) dehydrogenase 1 alpha subcomplex subunit. Their methylation may relate to more 246 efficient respiration because of CR to extract maximum energy from food. Interestingly, in yeast, CR is associated with increased longevity due to a reduction in NADH levels because 248 of NADH dehydrogenase activity [76] to create NAD+ (oxidised form). NADH is a competitive inhibitor of yeast sirtuin, leading to its activation on decreased NADH levels 250 [76,77]. This could also be occurring in Daphnia under CR if methylation of the NADH gene results in the expected increase in gene expression. The RNA methylation group for 252 hypermethylated genes (GO:0031167) contains two methyltransferase protein 20s (not DNA methyltransferases) which do not have a clear link to CR in the literature. 254

Hypomethylation under CR 256
Phospholipid/glycerolipid metabolism is reduced under CR (Table 3, GO:0006644 and GO:0046486), and both processes are associated with the same genes. These genes are 258 GPI inositol-deacylase, cardiolipin synthase, phosphatidylinositol-glycan biosynthesis class W protein, and phosphatidylserine synthase. Assuming that decreased methylation lowers gene 260 expression, this result is in keeping with previous work on effects of CR on phospholipids.
In mice myocardium, phospholipids undergo a reduction in mass and are remodelled when 262 facing CR [78], which is speculated to maximise energy efficiency. The same drop in phospholipids was observed in humans undergoing acute CR [79], and more generally 264 reduces the rick of atherosclerosis and heart disease [80,81]. This is potentially a further common mechanism of response to CR in which DNA methylation is an important reflect the broad-range of substrates ATP binding protein are able to efflux.

Conclusion 274
We have shown that caloric restriction effects the methylation status of a subset of genes, despite the low overall CpG methylation found in Daphnia. There is a strong concordance 276 between these results and CR experiments in humans, mice and yeast among other species.
We show that hypermethylated genes and processes are in line with upregulation in 278 previous CR and hypomethylated genes with downregulation. Although we have focused on the effect of caloric restriction on DNA methylation status, there are alternative potential 280 epigenetic responses to CR, including small RNAs (sRNAs) and histone modifications.
Previously, we established that CR induces differential miRNA expression in D. magna under 282 an equivalent experimental design [82], but other sRNAs, for example piRNAs and tsRNAs, could also have a role in CR-dependent gene regulation [83][84][85]. Histone modifications in 284 response to CR or protein restriction (PR) are known from work on humans, rats and mice [86][87][88][89] and are proposed to increase longevity [87] by delaying and repressing ageing-286 related processes and diseases. Future studies could vary a range of dietary components (overall calories, proteins or fatty acids) and examine the joint effects of a range of 288 epigenetic mechanisms.

Daphnia preparation and experiment 292
Six replicates of control (i.e. well-fed) Daphnia magna were compared to six replicates of 294 caloric restricted Daphnia to identify differentially methylated regions. We used a single clone (known to us as Clone 32) of D. magna, collected from the Kaimes pond near Leitholm 296 in the Scottish Borders [90]. Maternal lines were first acclimatized for three generations.
For this, individuals were kept in artificial pond medium at 20°C and on a 12h:12h light:dark 298 cycle and fed 2.5 x 10 6 cells of the single-celled green algae Chlorella vulgaris daily. Following three generations of acclimatisation (detailed in [82]), 40 offspring from each mother were 300 isolated and split to form a replicate. Twenty were fed a normal diet of 5x10 6 algal cells/day and the remaining twenty that were fed a calorie restricted diet of 1x10 6 algal cells/day. Each 302 replicate was split and reared in four sub-replicate jars of five Daphnia which were subsequently pooled at DNA extraction. Hence, normal food and calorie restricted 304 replicates and were paired by mother and each consisted of twenty Daphnia total. The experiment was ended after the birth of 2 nd clutch (approximately day 12 of the treatment 306 generation). Daphnia were ground by motorized pestle in Digsol and proteinase K and incubated overnight at 37°C and stored at -70°C until DNA extraction.  [82]. This was to increase mapping efficiency, and accuracy of the analysis, by reducing polymorphism 330 between the reference (assembled from a different clone) and our data. inspection of Bismark m-bias plots. TrimGalore! was also used to remove base calls with a Phred score of 20 or lower, adapter sequences, and sequences shorter than 20 bp. FastQC 336 0.11.4 [92] was used to before and after quality control to inspect the data. Bisulphite calls were made with Bismark 0. 16.3 [40]. Bismark alignments were performed with options "-338 score_min L,0,-0.6". PCR duplicates were removed using deduplicate_bismark script.
Bismark reports indicated that libraries were not fully bisulphite converted and raw 340 methylation calls were approximately 3% for CpG, CHH and CHG sites. Previous research has shown that CHH and CHG methylation is negligible in the D. magna genome [19]. As a 342 result, we used the filter_non_conversion script to remove reads that contain either CHH and CHG methylation sites as diagnostic of a non-bisulphite converted read. Finally, 344 methylated sites were identified using bismark_methylation_extractor and reports created with bismark2report. 346 Further variants were predicted per replicate using the unconverted library reads by 348 following the GATK pipeline [93,94] and converted strain 32 reference sequence. One round of base quality score recalibration was sufficient using variants previously identified 350 from strain 32. Sites with single nucleotide polymorphisms at methylated positions were removed from the analysis using BEDtools [95]. 352

Differential Methylation Analysis 354
All analyses were performed using methylation calls from the bisulphite converted libraries only. Hierarchical sample clustering of genome-wide methylation patterns across replicates 356 was generated using methylKit [42]. The six-normal food-and six calorie restricted replicates were then compared using bsseq [43] Bioconductor package in R to identify 358 regions of differential methylation. We ran bsseq with a paired t-test to control for the batch effect of mother on methylation. DMRs were selected using t-statistic cutoff of -4.6 360 and 4.6, a greater than 0.1 average difference in methylation between groups, and at least three methylated CpGs. Genes overlapping DMRs were identified using the D.