Analysis of DNA methylation at birth and in childhood reveals changes associated with season of birth and latitude

Background Seasonal variations in environmental exposures at birth or during gestation are associated with numerous adult traits and health outcomes later in life. Whether DNA methylation (DNAm) plays a role in the molecular mechanisms underlying the associations between birth season and lifelong phenotypes remains unclear. Methods We carried out epigenome-wide meta-analyses within the Pregnancy And Childhood Epigenetic Consortium to identify associations of DNAm with birth season, both at differentially methylated probes (DMPs) and regions (DMRs). Associations were examined at two time points: at birth (21 cohorts, N = 9358) and in children aged 1–11 years (12 cohorts, N = 3610). We conducted meta-analyses to assess the impact of latitude on birth season-specific associations at both time points. Results We identified associations between birth season and DNAm (False Discovery Rate-adjusted p values < 0.05) at two CpGs at birth (winter-born) and four in the childhood (summer-born) analyses when compared to children born in autumn. Furthermore, we identified twenty-six differentially methylated regions (DMR) at birth (winter-born: 8, spring-born: 15, summer-born: 3) and thirty-two in childhood (winter-born: 12, spring and summer: 10 each) meta-analyses with few overlapping DMRs between the birth seasons or the two time points. The DMRs were associated with genes of known functions in tumorigenesis, psychiatric/neurological disorders, inflammation, or immunity, amongst others. Latitude-stratified meta-analyses [higher (≥ 50°N), lower (< 50°N, northern hemisphere only)] revealed differences in associations between birth season and DNAm by birth latitude. DMR analysis implicated genes with previously reported links to schizophrenia (LAX1), skin disorders (PSORS1C, LTB4R), and airway inflammation including asthma (LTB4R), present only at birth in the higher latitudes (≥ 50°N). Conclusions In this large epigenome-wide meta-analysis study, we provide evidence for (i) associations between DNAm and season of birth that are unique for the seasons of the year (temporal effect) and (ii) latitude-dependent variations in the seasonal associations (spatial effect). DNAm could play a role in the molecular mechanisms underlying the effect of birth season on adult health outcomes. Supplementary Information The online version contains supplementary material available at 10.1186/s13148-023-01542-5.

may also act as molecular mediators between seasonlinked environmental factors and health outcomes.In plants and animals, the role of DNAm in the regulation of several seasonally associated genes is well documented [18][19][20].In humans, the seasonal periodicity of gene expression profiles has been demonstrated for genes involved in immunity and physiology [21].In addition, there is evidence to suggest that environmental factors alter DNAm levels in a season-specific manner in healthy individuals [22].It is proposed that these changes in DNAm profiles may influence health outcomes later in life.Taken together, a growing body of evidence suggests a role for differential DNAm in the exposure-induced changes in gene expression and disease [23].However, studies on the seasonality of DNA methylation at birth in humans have not assessed variations across a wide range of geographies on a large scale.
Early life environmental exposures result in specific changes in DNAm that may persist throughout the life course.For example, maternal smoking in pregnancy results in a specific DNA methylation signature in offspring [24] that persists well into adulthood [25].Similarly, it is plausible that DNAm changes associated with birth season, a proxy for in utero exposures, could be a mechanistic pathway underlying the effect of season of birth on later life health and disease.Lockett et al. tested the effect of season of birth on the risk of allergic outcomes in adulthood as well as the association between

Background
Plants and animals adapt to altered seasonal cues such as temperature and light by regulating their transcriptional programmes.It has been well established in humans that many traits, including physiological (e.g.blood pressure and cholesterol), behavioural traits (e.g.conception, suicidal tendencies), and life span, display distinct seasonal patterns [1][2][3][4][5].The incidence of complex diseases such as cardiovascular and autoimmune disease, and psychiatric disorders shows seasonal fluctuations [6][7][8][9].Several studies have demonstrated an association between season of birth and development of disease both in childhood and later in life.These adverse events include asthma and allergy-related diseases [10][11][12][13], neonatal immune development [14], multiple sclerosis [15] and schizophrenia [9,16], suggesting long-lasting effects on a broad range of characteristics.Exposures that vary seasonally such as outdoor temperature, humidity, ultraviolet (UV) light, pollen levels, Vitamin D, melatonin levels, air pollution and availability of nutrients have been hypothesised to drive the season-of-birth-linked disease phenotypes (reviewed in [17]).Uncovering the molecular mechanisms underlying the effect of season of birth on lifelong phenotypes can help elucidate the biological basis for the observed associations.
Epigenetic modifications, including alterations in DNA methylation (DNAm), play a key role in the gene regulation of a variety of cellular processes.Such modifications season of birth and differential DNAm at birth, and at age 18 [26].They demonstrated a season-associated pattern of DNAm at some CpG sites associated with apoptotic genes as well as season-dependent risk patterns of allergic diseases such as eczema and rhinitis at 18 years of age.The study was limited in sample size (N = 175) but was validated in a small independent cohort.
Here, we report the results of a meta-analysis of epigenome-wide association studies (EWAS) from 27 independent cohorts (N = 12927 participants) from the Pregnancy And Childhood Epigenetics (PACE [27]) consortium.Blood DNAm data at birth and from children aged between 1 and 11 years were used to investigate the association between season of birth and differential blood DNA methylation at these two time points.In addition, we examined if the latitude where the children were born influenced the season-of-birth-associated DNAm profiles as exposures such as sunlight, relative humidity, pollen count, infectious agents and air pollution can vary with latitude in each season.

Methods
This meta-analysis study intended to investigate the epigenome-wide associations of season of birth and DNAm at two time points, at-birth, and childhood, on a large-scale using data from multiple cohorts across differing geographies.

Cohort participants
Cohorts participating in the PACE consortium with data on season of birth and DNAm at birth and/or later in childhood were invited to participate in a meta-analysis of associations between season of birth and DNA methylation.Twenty-one independent cohorts (N = 9418) contributed to the at-birth and 12 cohorts (N = 3610) to the childhood meta-analyses.Most cohorts contributed data either for at-birth or childhood analysis (21/27), while six cohorts (~ 18% of the at-birth data and ~ 45% of the childhood data) contributed to both time points.Less than 20% of the study population had data at two time points from the same child.Data for 60 individuals (EARLI, non-European) were excluded from the at-birth analysis due to insufficient numbers (as per the PACE recommendations of ≥ 15/group) of births in the spring season (n (spring) = 4, Additional file 1: Table S1B) resulting in a total of 9358 for the at-birth analysis (Fig. 1).All study participants originated from countries in the northern hemisphere (32.7-71.2°N).Additional file 1 (Table S1A-C) lists all cohorts along with the distributions of covariates used in the cohort-specific EWASs.Additional file 2 provides cohort-specific definition of covariates, DNAm data collection and pre-processing methods.EWASs and meta-analysis were carried out on at-birth and childhood samples separately.All participating cohorts were from the northern hemisphere (range 37.2-71.2°N).Latitude ≥ 50°N represents a subset of cohorts from 50 to 71.2°N (referred to as "higher latitude subgroup" in the text).Latitude < 50°N is a subset of cohorts from 32.7 to 50°N (referred to as "lower latitude subgroup" in the text.a 27 independent cohorts, but a total of 33 datasets as 6 cohorts contributed to both at-birth and childhood analyses.b At-birth blood samples (cord and heel prick).c Whole blood samples (age 1-11 years).d This analysis was not done due to its small sample size

Primary analysis
The outcome was DNAm at birth and childhood.For new-born samples, all but two cohorts used cord blood, while CBC and CHS used dried blood spots collected post-delivery (Additional file 1: Table S1).The primary exposure was season of birth as a categorical variable with autumn as the reference season.Seasons were defined as winter (December-February), spring (March-May), summer (June-August) and autumn (September-November) since all participants were from the northern hemisphere.

Secondary analysis
Associations of season of birth and DNAm at-birth were also investigated in samples stratified by the cohort latitude to explore its influence on DNAm signals.Cohorts were dichotomized into a "lower latitude subgroup" (32.7-50°N with N = 3838) and "higher latitude subgroup" (≥ 50°N with N = 5580).The choice of 50° as the cut off was arbitrary, although it provided ~ 8 h of sunlight/day at the winter equinox.We focussed on comparing (i) the impact of latitudes on the pattern of association between season of birth and methylation signals at-birth and (ii) the pattern of season of birth specific associations across the two time points (at birth and in childhood) in each latitude.However, due to the small sample size contributing to childhood DNAm analysis in the lower latitude (N = 936), the DNAm results of the child cohorts in the higher latitude subgroup alone were used for comparison with the at-birth DNAm associations.

Measurement and pre-processing of DNA methylation
Details of cohort-specific methylation assays, data handling, quality control and normalisation are given in Additional file 2. Briefly, DNAm assays were carried out on bisulfite-converted genomic DNA extracted from either cord blood/heel prick samples collected at delivery or shortly after delivery and from whole blood samples from children aged 1-11 years.Methylation data for these samples were generated using illumina Infinium ® Human-Methylation450 (most cohorts) or EPIC BeadChip assay (EAGeR, IoWF2, BiB of HELIX, (Additional file 1: Table S1A).All cohorts estimated methylation levels as beta values (β) using the formula β = M M+U +100 where M and U are methylated and unmethylated signal intensities, respectively.Normalised beta values were used in preference to M values M = log 2 β 1−β to afford direct comparison of our findings with similar findings in the literature where most epigenetic association studies employ beta values.While all cohorts followed a pre-specified analysis plan including recommendations for bioinformatics pipelines (Additional file 3), the individual cohorts conducted probe filtering, normalisation of methylation data and correction for batch effects using their own preferred methods.Probes on sex chromosomes and control probes were removed prior to meta-analysis.Cross-reactive probes and polymorphic probes were dealt with post-meta-analysis.Each cohort also provided information on genomic inflation in their cohort-specific EWAS models.Leave-one-out analyses were carried out to check if inflation in any one cohort unduly influenced the overall meta-analysis results.

Covariates
All cohort-specific models, at birth or childhood, included the following covariates: maternal smoking status during pregnancy, maternal socio-economic status (SES), maternal age and gestational age at delivery, new-born's sex, and estimated cell type proportions.These covariates have been reported to be associated with either season of birth or DNA methylation making them potential confounders of the association between season of birth and DNA methylation or both [26,[28][29][30].Cell type proportions were estimated by applying the Houseman method [31] to methylation data using cord blood [32] or peripheral blood reference panels [33] for at-birth and childhood analyses, respectively, as these reference-based methods were reported to be the best available methods for studies with large sample sizes [34].In addition, the childhood models included children's age and season of sample collection if data were available.Most cohorts used three categories for maternal smoking status during pregnancy (no smoking, quit early during pregnancy and smoked throughout pregnancy), but cohorts lacking this level of detail or with too few who smoked during pregnancy used any versus no smoking.The definition of maternal SES was cohort specific and was mainly based on maternal education, occupation and/or income.Ninety-six percent of the study population were of European ancestry despite being drawn from multiple countries.Ethnicity was, therefore, not controlled for in the models.Cohort-specific covariate data collection methods are described in Additional file 2. In general, most cohorts collected data through questionnaires/interviews and/or from medical records.Individual cohorts used their preferred batch correction methods, and they typically comprised adjustments for methylation profiling-specific variables such as plates, bead chip or other relevant technical covariates.Studies which used a sampling scheme included the sampling variable as a covariate.

Cohort-specific epigenome-wide association studies (EWAS)
Each cohort performed independent EWAS analysis following a pre-specified analysis plan and an R script (Additional file 3).The primary exposure variable, season-of-birth, was consistently coded by all cohorts with autumn as the reference season.Data on one child per family were included (randomly selected) in the cohortspecific EWAS if multiple new-borns/siblings from the same family were enrolled in the cohort.Only samples with complete data for all covariates were included in EWAS analyses.
Cohort-specific EWASs modelled the association between season of birth (exposure) and DNAm (outcome) at a CpG site of the at-birth or childhood blood DNAm data using analysis of covariance via robust regression applying Huber's weighting, implemented in the R function rlm().Birth seasons were defined as winter (December-February), spring (March-May), summer (June-August) and autumn (September-November). Briefly, cohorts generated a regression model comparing children born in spring, summer, and winter against autumn as the reference for each of the two time points (at-birth or childhood for which they had methylation data).Models were adjusted for cohort-specific covariates as described above and in Additional file 2.

Meta-analyses
Meta-analyses were carried out using fixed effect inverse variance weighted option of METAL [35].Meta-analyses for at-birth and childhood data were performed separately.All analyses were repeated by one other member of the team independently.Probes that were common to both EPIC and 450 k platforms alone were included in the meta-analysis [36].A total of 470870 probes were retained after the meta-analyses of EWAS summary results.Each meta-analysis result was checked visually using quantile-quantile (QQ) plots and quantitatively by estimating lambda values, an indicator of inflation.Bias and inflation were estimated following the BACON method of van Iterson et al. [37].Meta-analysed probes were annotated against hg19 reference genome and mapped to the nearest genes and genomic locations (CpG islands, shores, shelves) using the Illumina 450k manifest or ACME R package [38] when not available in the manifest.The annotated output from meta-analysis holds information on probes which overlap SNPs or are within 10 bases from the CpG sites (polymorphic probes).These were disregarded while interpreting the significantly differentially methylated CpGs.Associations were considered significant for CpGs that reached a statistical significance based on a false discovery rate (FDR) below 5% [39].The FDR-corrected probes were also checked for a more stringent Bonferroni criterion of satisfying a p value threshold of 1.06 × 10 −7 .

Differentially methylated regions
Differentially methylated region (DMR) analyses were carried out using comb-p in R [40] on the meta-analysed single-CpG EWAS results separately for at-birth and childhood samples.CpGs within 1000 base pairs were combined to define a region [41].The significance threshold was set to the Šidák-corrected p value of 0.05 in all DMR analyses.

Trait enrichment analysis
Trait enrichment analysis of top CpGs (raw p values for association in meta-analyses < 1.0 × 10 −5 ) from the primary and latitude analyses was carried out using EWAS Atlas [42].Traits were considered significant if the p values for their odds ratios were less than 0.05.

Metastable epialleles
Fisher's exact test was carried out in R (3.5.1) to compare the proportion of metastable epialleles amongst significant loci (p value < 0.05) in each of the birth seasons of new-born infants against the proportion in array background.The comparator was 2408 metastable epialleles listed in three separate studies [43][44][45].

Study population
Twenty-seven studies from the PACE consortium participated in the meta-analysis study to investigate possible associations between season of birth and DNAm at birth and in childhood [at-birth: 21 studies, children: 12 studies of which six cohorts contributed to both (Fig. 1, Additional file 1: Table S1A)].There was a total of 9358 and 3610 participants for the at-birth and childhood analyses, respectively (Fig. 1).A subgroup analysis, carried out on participants living in regions of ≥ 50°N, included 5580 new-borns (9 studies) and 2674 children aged 1-11 years (8 studies) [see Additional file 1: Table S1 The births were evenly distributed across the seasons (22-27.5%,Additional file 1: Table S1B-C).Approximately 50% of the babies and children included in the analyses were girls.The median maternal age at delivery was 30 years [range (min-max): 24-34.2]for the participants of the at-birth analyses and 30 years (range 26.7-32.1)for the childhood analyses.The median gestational age was 39.5 weeks (range 38.9-40.2) and 39.6 (range 36.7-40.2) for the at-birth and childhood participants, respectively.The proportion of mothers who smoked during pregnancy ranged from 0 to 15% in most cohorts in both the at-birth and childhood samples except for the mothers of IoW F2 (36.7%)Table 1 Differentially methylated CpGs (FDR < 0.05) in the birth seasons in at-birth and childhood analyses CpGs shown passed the FDR (False discovery rate) p value threshold of 5% in meta-analyses of EWAS summary results carried out with 21 and 12 cohorts for at-birth and childhood blood samples, respectively.The p values for the CpG sites that also passed the Bonferroni-corrected significance threshold of 1.06 and NEST (23.2% and 20% for the two ethnic groups) for the at-birth analyses.

Association between season of birth and DNA methylation at birth and in childhood
For each time point, at-birth and childhood, the individual cohorts generated an EWAS model (adjusted for covariates) comparing the three seasons (winter, spring, and summer) with autumn as a reference.The lambdas for the six models varied between 1.02 and 1.16 (Additional file 4: Table S2).Meta-analyses were then carried out on the summary results from the cohort-specific EWAS models for each of the time points (i) to identify at-birth DNAm signals associated with season of birth (winter, spring, or summer with autumn as the reference season) and (ii) to investigate whether such signals persist into childhood.Leave-one-out analyses showed that none of the individual cohorts unduly influenced the meta-analysis results.
The six CpGs that passed the pre-specified 5% FDR threshold for winter-, spring-and summer-born children for the at-birth and childhood analyses are presented in Table 1 and Figs. 2 and 3.The two FDR significant CpGs identified in the new-borns, cg26416241 (winter-born, mapped to NTM) and cg18978324 (spring-born, mapped to RABGGTA ), had statistical significance that passed the more stringent Bonferronicorrected threshold of 1.06 × 10 −7 (CpGs shown in bold, Table 1).The methylation levels of the at-birth CpGs were higher in the winter-and spring-born offspring) when compared to those born in the autumn.This trend was consistent across most of the cohorts (winter: 18; spring: 17 out of a total of 24 cohorts).There was no statistical evidence to suggest that the two CpGs (cg26416241, cg18978324) identified in the at-birth analysis (preliminary analysis) persisted into childhood (Additional file 4: Table S3).In the individual cohort samples of the relevant seasons, these two CpGs mostly had higher methylation levels when compared to the autumn levels in both at-birth samples (75-77% of the cohorts) and childhood samples (58% for both CpGs).
Mean differences in methylation between autumn and each of the other two seasons were small in absolute terms [mean change in DNAm winter vs. autumn (cg26416241): 0.0034 (0.3%) and spring vs. autumn (cg18978324): 0.006 (0.6%)].
In the analysis of childhood samples (aged 1-11 years), differential DNAm was observed at four CpG sites (cg19416462, cg01656588, cg03263237 and cg15437053) in the summer-born children (FDR < 0.05, Table 1) when compared to the autumn-born children.One of the four CpGs, cg03263237, was a cross-reactive probe [46] and was excluded from further biological interpretations.The observed between-study heterogeneity was low to moderate (I 2 < 50) for the FDR-significant CpGs.One of the four CpGs, cg19416462 (mapped to CD83), also passed the Bonferroni-corrected threshold p value.The remaining CpGs, after excluding the cross-reactive probe of the four childhood CpGs, two showed decreased methylation levels in children born in the summer compared to those born in the autumn in 92-100% of the cohorts, while cg15437053 was more methylated than the same CpG in autumn born children in 57% of the cohorts.The corresponding three CpGs in the at-birth samples of summer-born infants lacked even nominal significance but showed an overall tendency towards hypomethylation compared to the autumn-born infants (57-61% of the contributing cohorts).The mean differences in methylation (childhood analysis: summer vs. autumn) were small in absolute terms [mean change in DNAm (cg19416462): − 0.4%, (cg1656588): − 0.7%, (cg15437053): 0.4%].

Differentially methylated regions (DMRs) at birth and in childhood associated with season of birth
Several DMRs, containing at least two CpG sites, were associated with season of birth both in new-borns and in children aged 1-11 years (Table 2).There were 24 DMRs mapping to 27 known genes in the at-birth analysis.CpGs of the DMRs were mostly hyper-methylated, and their direction in the DMRs was the same in most cohorts (Additional file 5: Table S4).The DMRs identified (See figure on next page.)Fig. 2 Manhattan plots of the meta-analysed data showing the association of season of birth and DNA methylation at birth using neonatal blood data (A) and during childhood using whole blood data (B).Birth seasons with FDR-significant CpGs when compared to autumn as the reference season alone are shown in the figure.The blue and red lines indicate the threshold p values for false discovery rate (FDR) and the Bonferroni adjustments, respectively.The observed p values on the Y-axis are from models adjusted for covariates and cell types for the seasons indicated in the figures when compared to autumn as the reference season.Genes associated with the CpGs (circled) are indicated.All cohort-specific EWAS analyses were adjusted for gender of the child, gestational age at delivery, maternal age at delivery, maternal smoking during pregnancy, maternal socio-economic status, batch, child's age at the time of sample collection (in the case of childhood samples and if data were available) and estimated cell proportions.*This CpG (cg01801443, location: intergenic) has a SNP within 10 base pairs and is not included in  and the genes mapped to these differential methylation signals were unique to the seasons apart from two DMRs mapped to HLA-F (in both winter and spring) and to AZU1 (in spring and summer, Table 2).The childhood analysis identified 32 DMRs mapped to 35 genes.CpGs of the DMRs were mostly hypo-methylated, and their direction in the DMRs was the same in most cohorts (Additional file 5: Table S4).Like in the at-birth analyses, the season of birth associated DMRs of childhood samples was unique to the season of birth except for two mapped genes, LOC441666 in the winter-and summer-born and HTR2A in winter-and spring-born children.There was no overlap between the at-birth and childhood DMRs.

DNA methylation at birth and in childhood by latitude
We performed additional at-birth and childhood metaanalyses in subgroups stratified by latitude to check if the observed temporal associations varied with latitude.Cohorts, all from the northern hemisphere, were divided into two latitude groups-higher latitude (≥ 50°N, 9 cohorts, N = 5580, 59.3%) and the lower latitude [37.2-50°N (12 cohorts, N = 3823, 40.7%)].
In the at-birth analyses, we identified a single CpG, cg06251958, from the spring-born babies of the higher latitude subgroup (FDR-adjusted p value = 0.001, Table 3) and none from the winter-or summer-borns (all against autumn).This CpG (mapped to PSMC2) also passed the Bonferroni p value threshold (1.06 × 10 −7 ).There were seven significant (FDR p value < 0.05) at-birth CpGs, six in the spring-born and one in the winter-born infants in the lower latitude subgroup analyses.One of these seven CpGs, cg23369114, was a cross-reactive probe [46].There was no overlap in the FDR-adjusted differentially methylated CpGs between the seasons in each latitude subgroup nor between the two latitude subgroups (higher vs. lower) in any given birth season.The characteristics of cohorts from the higher and lower latitudes were similar apart from the latitude of the place of birth (Additional file 1: Table S1B).These results, therefore, suggest that the temporal associations between season of birth and DNAm were also latitude specific.
The childhood analyses of the higher latitude cohorts revealed two FDR-significant (also Bonferroni-significant) CpGs in the winter-born children and none in the spring-and summer-born.Furthermore, there was no overlap in CpG methylation signals between the at-birth and childhood samples.We did not perform childhood DNAm analysis for the lower latitude subgroup due to its small sample size (N = 936, Fig. 1).

Differentially methylated regions (DMRs) at birth and in childhood by latitude
Interrogation of at-birth DMRs of the higher latitude subgroup (≥ 50°N) identified 18 DMRs mapped to 20 unique genes in the winter-born and 25 DMRs (23 mapped genes) in the summer-born infants (Table 4).In contrast, there were relatively fewer DMRs (n = 7) in the spring-born infants of the higher latitude cohorts.Most of the at-birth DMRs identified in the new-borns were unique to the season of birth (winter, spring, or summer) except for DMRs which mapped to ELFN1, SEMA5B and LOC154449 (identified in winter-and summer-born children) and AZU1 (in spring-and summer-born).
Examination of the biological functions of the genes mapped to the at-birth higher latitude DMRs revealed that most of the genes have roles, amongst others, in the functions of the central nervous system (CNS) including psychological disorders, cancer or inflammation and immunity (see Additional file 6: Table S5, for a more detailed, but not exhaustive, list of genes associated with at-birth DMRs and their functions).The DMRs identified in the new-borns from the lower latitude (< 50°N) cohorts were also unique to seasons of birth (Additional file 7: Table S6).None of the DMRs identified in the higher latitude subgroup were present in the lower latitude samples (compare Table 4 and Additional 7: Table S6).
Interrogation of childhood DMRs of the higher latitude subgroup (≥ 50°N) identified 47 DMRs which were mapped to 54 genes (winter: 21 DMRs and 24 genes, spring: 16 DMRs and 18 genes, summer: 10 DMRs and 12 genes, Additional file 7: Table S7).Most of these DMRs were unique to the season of birth (winter, spring, or summer) except for a DMR that mapped to S100A13 identified in both winter-and spring-born children and a second DMR that mapped to MIR7159 identified in both winter-and summer-born children.Furthermore, there was little overlap between the DMRs identified in the atbirth (Table 4) and childhood analyses of the higher latitude subgroup (Additional file 7: Table S7).
Four at-birth DMRs of the summer-born babies from the higher latitude subset mapped to imprinted genes, Fig. 3 Quantile-quantile (Q-Q) plots for post-meta-analysis models for association between seasons of birth and DNA methylation (primary analyses) for at-birth (A and B) and childhood (C and D) samples.The Q-Q plots were generated by plotting observed p values (y-axis) against the expected uniform distribution of p values under the null hypothesis of no association (x-axis).Lambda and bias were estimated using BACON method [37].All cohort-specific EWAS analyses were adjusted for gender of the child, gestational age at delivery, maternal age at delivery, maternal smoking during pregnancy, maternal socio-economic status, batch, child's age at the time of sample collection and estimated cell proportions  GNAS/GNAS-AS1, MEG3 and MEST, with two of them mapping to the GNAS locus (GNAS and GNAS-AS1).These imprinted genes were specific for the summerborn babies and were not identified in the new-borns of spring or winter (Table 4).Whilst DMRs associated with GNAS/GNAS-AS1 and MEG3 were present only in the at-birth samples from the higher latitude subgroup, a DMR associated with the MEST gene (mesoderm-specific transcript) was present in both at-birth and in early childhood from the higher latitude subgroup (compare Table 4 with Additional file 7: Table S7.

Trait enrichment analysis
Trait Enrichment analysis was carried out on CpGs with unadjusted meta-analysed association p value < 1.0 × 10 5 using EWAS Atlas [42].This analysis revealed that, apart from enrichment for traits related to cancers or CNS-associated disorders, the enriched traits are mutually exclusive between at-birth and childhood analyses (Additional file 9).The traits with higher odds ratios for enrichment in the at-birth analysis were vitamin B12 supplements, polycystic ovary syndrome, maternal stress, systemic lupus erythematosus, and Behcet's disease.These traits were enriched in the at-birth but not in the childhood samples and were mostly in new-borns of mothers with the early months of their pregnancy in the months of autumn to winter, which may indicate influence of exposures in the earlier parts of foetal development.
Differentially methylated regions (DMRs) were identified using comb-p.The inputs for comb-p were the annotated meta-analysed outputs of at-birth and childhood blood EWASs     A similar trend of mutual exclusivity of enriched traits was seen when the new-borns of higher (≥ 50°N) and lower (< 50°N) latitudes were compared.While the higher latitude new-borns alone had traits enriched in exposure to pollutants (particulate matter, nitric oxide, and maternal arsenic exposure) and folic acid supplement, the lower latitude offspring had unique enriched traits associated with response to different treatments, mixed connective tissue disease, asthma, puberty, and household socio-economic status.Parental arsenic exposure and autistic spectrum disorder were common in the newborns and young children of higher latitude, although weakly enriched in the childhood analysis.

Metastable epialleles
We tested CpGs with unadjusted p value < 0.05 in metaanalyses for enrichment of metastable epialleles (loci whose epigenetic modifications are established during early embryonic development [47]).Our analyses found no evidence of enrichment of metastable epialleles.

Discussion
This study investigated the association of season of birth and DNAm at two time points (at-birth and in childhood), and in latitude subgroups, in multiple cohorts recruited from regions across the northern hemisphere.Season of birth can be a proxy for seasonal variations in the length of the day, temperature, availability of sunlight, exposure to UV light, pollen, nutrition, seasonal infectious agents, air pollution, C-section and many more.Here, we show evidence for the existence of season of birth specific associations with DNAm at-birth as well as in childhood.The differential methylation patterns in at-birth samples followed a season-dependent annual

TNXB intron
Differentially methylated regions (DMRs) were identified using comb-p.The inputs for comb-p were the annotated meta-analysed outputs of at-birth EWASs CHR chromosomes, UTR untranslated region, TSS transcription start site, cds coding sequence a Reference season in EWAS of individual cohorts: autumn b Regions with Šidák p value < 0.05 c P value for the most significant CpG in a DMR identified by comb-p periodicity (temporal effect).In addition, we found suggestive evidence for latitude-specific fluctuations in DMRs (a spatial effect).To the best of our knowledge, our study is a first to demonstrate such season of birth effects on DNAm as well as the latitude specificity in these effects in a large study population with participants from diverse geographical locations.

Season of birth specific differential DNA methylation
Our primary analyses found epigenome-wide significant season of birth associations at two and four CpGs in the at-birth and childhood blood samples, respectively.One of the four CpGs in the childhood samples was found to be a cross-reactive probe [46] and should be interpreted with caution.The remaining differentially methylated CpGs were unique for the birth seasons in the at-birth and childhood samples.Furthermore, there was no overlap in the CpG methylation signals observed across the two time points investigated (at-birth and childhood).Several DMRs associated with the birth seasons were identified in the at-birth and childhood analyses.CpG sites that are close to each other are known to be comethylated, and together these sites can function as a regional unit [48].At the individual CpG level, the exposure-induced changes in methylation of these CpGs are often small and with weak statistical evidence.However, analysis of pre-defined regions as functional units, each containing several co-methylated CpGs, may provide more statistical power to detect the associations with methylation signals.Half of the DMRs identified in this study were in gene locations such as coding sequences, transcriptional start sites and/or untranslated 5′ or 3′ ends.The DMRs identified were in the vicinity of genes with known functions, amongst others, in tumorigenesis, psychiatric/neurological disorders or inflammation and immunity.Most DMRs and their mapped genes in the atbirth (24 out of 26) and childhood (30 out of 32) DMRs were specific to the birth seasons.The seasonal specificity of DMRs was consistent with the observed temporal effect of DNAm signals at the CpG level both at birth and in childhood.Like in the case of differentially methylated season-of-birth-associated CpGs, there was little overlap between most of the at-birth and childhood DMRs identified which may indicate that the at-birth methylation signals associated with season of birth did not persist over time.

Association between season of birth and DNA methylation is latitude dependent
Latitude-stratified analyses (higher latitude: ≥ 50°N and lower latitude: < 50°N) were carried out to explore the impact of latitude on seasonal differences in at-birth and childhood DNAm.One of the significant CpGs of the at-birth lower latitude analysis, cg23369114 (Table 3) was a cross-reactive probe identified by Chen et al. [46] which warrants caution while interpreting the results.Like in the primary analyses, temporal associations (season-specific) were observed in the two latitude subgroups for the at-birth and the higher latitude childhood analyses both at the CpG (6/7 CpGs) and DMR levels.Furthermore, the DMRs of at-birth and childhood analyses of the higher latitude subgroup were mutually exclusive except for two DMRs which were present at both time points.This is consistent with the findings in the primary analysis that majority of the at-birth associations between season of birth and methylation signals did not persist during childhood.Like in the primary analyses, the observed associations between season of birth and DNAm signals (CpG/DMR) were unique for the latitudes in each of the seasons.These findings point to a latitude-dependent spatial effect of the association between DNAm and season of birth.The absence of sustained signals from birth to childhood in this study does not necessarily imply absence of persistence in general.There could be several reasons for this.The smaller sample size of the childhood analysis compared to that of the at-birth analysis might be one of the reasons.Choice of cord blood for the at-birth studies could have masked the birth to childhood persistence of DNAm signals in this study.For example, Paquette and Marsit highlighted the importance of placenta samples for studies on the influence of intrauterine HTR2A gene expression on early childhood and/or lifelong health outcomes [49].The current study identified a DMR on chromosome 13 (47472025-47472454) in the winter-born children of the primary and the higher latitude subgroup childhood analyses.The same DMR, mapped to HTR2A (codes for serotonin receptor), was present in the springborn children of the primary childhood analysis, but not in any of the at-birth analyses that were based mostly on cord blood data.Yet another factor could be the absence of paired comparisons in our study.The at-birth and childhood comparisons were done with > 80% of unpaired data as most cohorts contributed to either at-birth or childhood data.This meant that DNAm measurements at the two time points did not come from the same child.This might have made the true persistent signals, if there were any, undetectable.An analysis with paired data from the six cohorts (data available from the same child at both the time points) would have had limited power to detect methylation signals of small DMR sizes as such data were available for < 20% of the study population.The heterogeneity of the age group in the childhood analysis may also have contributed to the non-persistence of the at-birth DNAm signals.Age effects may be strong, and the age range was wide (1-11 years).It is likely that there may be some remaining effect despite controlling for child's age in the models.
It may be noted that persistence of DNAm signals from birth to childhood is not necessarily needed for differential methylation to have an effect.Differential methylation at birth, which affects developmental programmes at a critical period, may have long-lasting effects, even if the differential methylation does not persist into early childhood.
The season-related DNAm signals (CpGs and DMRs) identified in the primary and latitude analyses were non-overlapping.The primary analysis investigated the seasonal variations disregarding the information on the latitude of birth of the children.The results make use of the distribution of methylation signals across the seasons only.In contrast, the latitude analysis is a stratified analysis which examines the above seasonal variations in two separate latitude regions (lower and higher latitudes in the northern hemisphere).The results of such an analysis look for differences in the distribution of methylation signals across the seasons in two separate latitude subgroups which are likely to be different from that of the primary analysis.However, the DMR mapped to HTR2A was present both in the winter-and spring-born children of the preliminary childhood analysis as well as the winter-born children of the higher latitude childhood analysis.Alternatively, the non-overlapping methylation signals observed in the latitude subgroups, despite smaller sample sizes and stringent FDR cut offs, are likely to be the enriched latitude-specific signals which were not visible in the primary analysis.

Genes implicated in the season of birth association of DNA methylation
Several of the differentially methylated probes (DMPs) and DMRs identified in our analyses mapped to genes with well characterised functions.We highlight below some of the genes mapped to DMPs or DMRs identified in the at-birth and childhood analyses of the higher latitude subgroup which are likely to be associated with outcomes relevant to season of birth.
The most studied link between season of birth and a disease outcome is that of schizophrenia (SZ).In the general population in the northern hemisphere, individuals born in late winter and early spring have an increased risk of developing SZ later in life [9,50,51].In addition, the prevalence of SZ has also been shown to increase with increasing latitude with the highest prevalence rates at the poles [52,53].Correlation patterns of risk factors and SZ identified pre-natal vitamin D deficiency and infections like influenza and toxoplasmosis which are more prevalent in higher latitude and cold climates, as major risk factors of SZ.Interestingly, DNAm of cg12022621, mapped to the LAX1 gene (Lymphocyte Transmembrane Adaptor 1), has been shown to be associated with severity of certain symptoms of schizophrenia (SZ) in a case-control study [54].LAX1 is also the mapped gene for one of the top-ranked at-birth DMRs (Šidák-corrected p value: 7.2 × 10 −8 ) spanning six CpGs including cg12022621 (crude p value: 3 × 10 −4 , Table 4) and was identified only in the winter-born group of the higher latitude subset in the at-birth samples.While LAX1-associated SZ was absent from the childhood DMR analyses (Table 4), published literature reveals SZ to be one of the most highlighted outcomes of the major depressive disorder (MDD) family of mental health outcomes in winter-and spring-born children (Additional file 6: Table S5 and references therein).The corresponding associated genes in the childhood analyses were LY6G5C, HTR2A, SHANK 1 and NKAPL, in the winter-born children and, SNTG2 in the spring-borns [55][56][57][58][59].
Of note, SLC17A9, mapped by a DMR unique to the at-birth analyses of infants of the higher latitude subgroup born in winter (Additional file 6: Table S5), is associated with a skin-specific autoinflammatory disease, disseminated superficial actinic porokeratosis (DSAP) [60].Induction and exacerbation of DSAP are known to be linked to exposure to sunlight or artificial ultraviolet radiation [61].Similarly, other genes linked to skin disorders include PSORS1C3 (spring-born), and LTB4R (winter-born), both of which are known to be linked to the autoimmune disease psoriasis [62][63][64].Interestingly, prevalence of psoriasis is higher in populations of higher latitude regions of northern Europe [65].The DMR-associated genes responsible for skin disorders in the winterborn children of higher latitude subgroup were different from those in the at-birth samples.The childhood genes were MIR-1914, MIR-647 and MAP3K8 (common wart) and TNF (psoriasis) (Additional file 8: Table S8) [66][67][68].
A DMR on chromosome 10, identified in the summerborn children of the higher latitude subgroup childhood analysis, was mapped to FGFR2 gene.Interestingly, alterations in the FGFR2 gene seen in several cancer types made it a target for the development of personspecific treatments [69][70][71][72].DNAm status in the proximity of the FGFR2 gene, such as the one observed in this study, may further influence the FGFR2 alterations in the cancer patients.Therefore, knowledge of methylation status around FGFR2 of the cancer patients will help fine tune the design of the person-specific therapies further.
Several genes mapped to DMRs of at-birth blood samples of the winter-born of the higher latitude, not found in the lower latitude, had links to immunity and inflammation (LAX1: T-/B-/NK-cell activation; TLR5: pathogen recognition; LTB4R: pathogenesis of inflammatory diseases; PSORS1C3: psoriasis; TCL1A: modulation of immune responses; ADAM5: autoimmune diseases) [54,[62][63][64][73][74][75][76].LTB4R, leukotriene B4 receptor, appears to be the most epigenetically divergent gene in the peripheral blood of humans.LTB4R has been linked to a variety of inflammatory diseases such as asthma [77], allergic airway inflammation [78], inflammatory arthritis [79], atherosclerosis [80], inflammatory bowel disease [81] and psoriasis [82].This is significant as the pregnant mothers of the winter-born babies would have been exposed to various seasonal factors that are responsible for allergy/asthma and other causal exposures of inflammation during the first and the second trimesters of their pregnancy.Lockett et al. hypothesised that the season-associated DNAm of allergic diseases such as eczema most probably arose postnatally since no such association was observed at birth (in cord blood) in their study [26].However, the lack of evidence for a prenatal association at an epigenome-wide level could also have been due to the modest sample size of their study (N = 175).
A CpG, cg003488551, mapped to C7orf50 and DMR containing this CpG were reported to be associated with prenatal exposure to particulate air pollution in a meta-analysis by Gruzieva et al. [83].Major air pollutants show seasonal patterns with highest concentrations in the indoor heating seasons of November to February (winter months) in the northern hemisphere [84].It is possible the summer-born babies in the current study were exposed to indoor pollutants during most part of their in utero life.However, our study did not identify the same methylation signals in babies born in the higher latitude subgroup with extended and extreme winter conditions.This study also did not identify any known vitamin D-associated CpGs or genes mapped to DMRs in cohorts from the higher latitude (≥ 50°N).

Imprinted genes and metastable epialleles
Our study also identified DMRs associated with three imprinted genes in the at-birth analyses of the summer-born infants of the higher latitude, GNAS/GNAS-AS1, MEG3 and MEST (PEG1).Almost all the CpGs found in the DMRs of these imprinted genes were hyper-methylated.In addition to being season of birth and latitude specific (only in the summer-born babies of the higher latitude), they were not found in the older children from the higher latitude subgroup.However, a DMR associated with the MEST gene (mesoderm-specific transcript) was present both at birth (summer-born) and in early childhood (spring-born) from the higher latitude subgroup.Many imprinted genes function as regulators of embryonic or neonatal growth and may therefore influence a spectrum of heritable outcomes later in life.The majority of imprinted genes are expressed in the brain, and methylation of these genes in their imprinting control regions (ICRs) has been implicated in neuropsychiatric disorders (reviewed in [85]).For example, GNAS is a complex imprinted locus with five gene products and multiple DMRs in four of these genes.These DMRs are shown to be associated with a genetic disorder known as pseudohypoparathyroidism type-Ib (PHP-Ib) (reviewed in [86]).Furthermore, a UK Biobank GWAS study of 113,000 individuals with insomnia identified GNAS as a potential gene candidate in females [87].A previous study reported that the retention of a sex-specific association between a hypermethylated DMR associated with MEST and weight status from birth to early childhood [88].
Metastable epialleles are variably methylated loci with cross-tissue methylation signatures indicative of establishment in the early embryo [44].They, therefore, provide a useful tool for examining the timing of exposure driven DNAm changes in easily accessible tissues such as blood that may serve as a proxy for patterns of systemic methylation.Silver et al. demonstrated elevated DNAm levels at putative metastable epialleles in rural Gambian children who were conceived during the rainy season compared to those conceived in the dry season [47].We checked for the enrichment of metastable epialleles in each of the seasons of birth, but none were found, providing no evidence for a season of conception effect at these loci.

Implications of temporal and latitude-dependent associations between season of birth and differential DNA methylation
In the at-birth meta-analyses, there were only a few robust but weak genome-wide associations between CpG methylation and season of birth.DMR analysis is believed to be statistically more powerful than the analysis of individual probes as it combines methylation signals from nearby CpGs to give more reliable signals for associations between DNAm and exposures.It is possible that the seasonal variations in exposures influence DNAm via groups of CpG sites over an extended region as in the DMRs where contiguous differential methylation may be maintained.Our findings on the birth season-dependent variations of regional epigenetic signals concur with the findings of a study by Dopico et al. which demonstrated the existence of seasonal variations in gene expression profiles over a year in ethnically and geographically diverse populations [21].They attributed these seasonal variations in gene expression profiles to the seasonal changes in the cellular composition of blood.However, unlike in the study by Dopico et al., our models for association between birth season and DNAm were adjusted for cellular heterogeneity and therefore, the findings are unlikely to be the result of seasonality of blood com-

Strengths and Limitations
Our study investigated associations between season of birth and DNAm on a large scale with 9358 at-birth and 3610 childhood samples.Analyses with such a large sample size, especially for the at-birth samples, make it possible to interpret results with more confidence.The cohort-specific EWAS analyses were robust and adjusted for several pre-specified potential confounders, including cell type proportions.The cohorts originated from geographically diverse locations in terms of latitude in the northern hemisphere: 32.7-71.2°Nand 36.7-58.7°Nfor the at-birth and childhood samples.This enabled us to carry out latitude stratifications to investigate season of birth associations of DNAm.Of the 27 cohorts analysed in this study, more than 95% of the participants were of European ancestry.It remains to be seen whether our results are generalizable to other populations.
This study has limitations.This study was designed only to examine the associations between DNAm and season of birth as an exposure and not the factors for which season of birth is a proxy (see above).Heterogeneity of the age range in the childhood population in this study and its smaller sample size (childhood: N = 3569 vs. at-birth: N = 9358) made it harder to interpret with certainty the absence of at-birth DNAm signals in childhood data.Another limitation of our study is that the at-birth and childhood data were not from the same children for most of the cohorts (21/27 studies).A longitudinal study which follows up the same children at different latitudes in sufficient numbers is necessary to detect the persistence of at-birth DNAm signals from birth to childhood.The latitude analyses of this study have not been adjusted for longitude.The same latitude can have geographic locations that have very different climates, e.g.New York and Madrid (40°N).Furthermore, the choice of autumn as a reference season may have masked some of the other significant associations between season of birth and DNAm.This study was initiated as a follow-on from an earlier study on seasonality effects on asthma and allergy outcomes which showed the strongest effects on allergy phenotypes with autumn as the reference season [26].However, we expect to observe the circannual oscillations in DNAm even if a season other than autumn was used as the reference season.The cohorts which contributed the EWAS summary results pre-processed and analysed their data using their preferred pipelines and this may have influenced our results.However, Joubert et al. found that their results were robust to different normalization methods used across studies and cell type adjustment [24].Furthermore, Lussier et al. demonstrated that while different pipelines give different EWAS associations at a set significance threshold their magnitude and directions were consistent [91].

Conclusions
In this large epigenome-wide meta-analysis study, we provide evidence for an association between season of birth and differential DNAm that is unique for the seasons of the year (temporal effect) and suggestive evidence of a latitude (spatial) effect.Findings in this study add to the understanding of a potential epigenetic role in the seasonality of human disease.Our study suggests the existence of a circannual periodicity in DNAm patterns, much like the seasonal periodicity observed in gene expression profiles.

Fig. 1
Fig. 1 Schematic diagram for the season of birth-DNA methylation (DNAm) analyses.Numbers included in the analyses are indicated.PACE Pregnancy And Childhood Epigenetics, EWAS Epigenome Wide Association Study, DMP Differentially Methylated Probe, FDR False Positive Rate, DMR Differentially Methylated Region.EWASs and meta-analysis were carried out on at-birth and childhood samples separately.All participating cohorts were from the northern hemisphere (range 37.2-71.2°N).Latitude ≥ 50°N represents a subset of cohorts from 50 to 71.2°N (referred to as "higher latitude subgroup" in the text).Latitude < 50°N is a subset of cohorts from 32.7 to 50°N (referred to as "lower latitude subgroup" in the text.a 27 independent cohorts, but a total of 33 datasets as 6 cohorts contributed to both at-birth and childhood analyses.b At-birth blood samples (cord and heel prick).c Whole blood samples (age 1-11 years).d This analysis was not done due to its small sample size

10 − 7
are shown in bold.All cohort-specific EWAS analyses were adjusted for sex of the child, gestational age at delivery, maternal age at delivery, maternal smoking during pregnancy, maternal socio-economic status, batch, child's age at the time of sample collection (in the childhood analyses) and estimated cell proportions Coeff: regression coefficient (change in mean methylation compared to autumn reference); SE: standard error of the coefficient; Lambda: genomic inflation factor; CHR: chromosome; I 2 : a measure of heterogeneity between studies a Reference season for the EWAS analyses of individual cohorts: autumn b Estimated using the BACON method (van Iterson et al. [37]) c Hyper-methylation is indicated by "+" and hypo-methylation by "−" (minus sign) d One CpG (in winter-born; cg01801443; p value: 7.77 × 10 −10 ) with a SNP within 10 base pairs of the probe is not shown in the table eCross-reactive probe (Chen et al.[46])

CHR
chromosomes, UTR untranslated region, TSS transcription start site, cds coding sequence a Reference season in EWAS of individual cohorts: autumn b Regions with spatially adjusted and Šidák multiple testing corrected p value < 0.05 c P value for the most significant CpG in a DMR identified by comb-p the BACON method of van Iterson et al.Ref [37]: c + : hyper-methylation; -: hypo-methylation d Cross-reactive probe (Chen et al., Ref: [46]) position.Epigenetic marks including DNAm signals, whether acting at the individual CpG site or DMR level, are dynamic and may not originate stochastically.Oh & Petronis and Oh et al. proposed that DNAm variability occurring in a periodic fashion over twelve months (circannual oscillations) could contribute to variations in the severity of seasonal diseases [89, 90].We conclude that our findings on seasonal variations in DNAm signals including DMRs are in line with the circannual oscillations proposed by Oh & Petronis and are likely to reflect the influence of seasonal exposures on later life health events.

Table 2
Differentially methylated regions (DMRs) at birth and in childhood in different seasons

Table 3
Association between season of birth and DNA methylation in at-birth and childhood stratified by latitude ) are shown in bold.Cohort-specific EWAS models were adjusted for sex of the child, gestational age and maternal age at delivery, maternal smoking during pregnancy, maternal socio-economic status, batch, and cell type proportions Coeff: regression coefficient (change in mean methylation compared to autumn reference); SE: standard error of the coefficient; Lambda: genomic inflation factor; CHR: chromosome; I 2

Table 4
Season of birth associated differentially methylated regions (DMRs) at birth in cohorts from latitudes ≥ 50°N