DNA methylation in relation to gestational age and brain dysmaturation in preterm infants

Abstract Preterm birth is associated with dysconnectivity of structural brain networks and is a leading cause of neurocognitive impairment in childhood. Variation in DNA methylation is associated with early exposure to extrauterine life but there has been little research exploring its relationship with brain development. Using genome-wide DNA methylation data from the saliva of 258 neonates, we investigated the impact of gestational age on the methylome and performed functional analysis to identify enriched gene sets from probes that contributed to differentially methylated probes or regions. We tested the hypothesis that variation in DNA methylation could underpin the association between low gestational age at birth and atypical brain development by linking differentially methylated probes with measures of white matter connectivity derived from diffusion MRI metrics: peak width skeletonized mean diffusivity, peak width skeletonized fractional anisotropy and peak width skeletonized neurite density index. Gestational age at birth was associated with widespread differential methylation at term equivalent age, with genome-wide significant associations observed for 8870 CpG probes (P < 3.6 × 10−8) and 1767 differentially methylated regions. Functional analysis identified 14 enriched gene ontology terms pertaining to cell–cell contacts and cell–extracellular matrix contacts. Principal component analysis of probes with genome-wide significance revealed a first principal component that explained 23.5% of the variance in DNA methylation, and this was negatively associated with gestational age at birth. The first principal component was associated with peak width of skeletonized mean diffusivity (β = 0.349, P = 8.37 × 10−10) and peak width skeletonized neurite density index (β = 0.364, P = 4.15 × 10−5), but not with peak width skeletonized fraction anisotropy (β = −0.035, P = 0.510); these relationships mirrored the imaging metrics’ associations with gestational age at birth. Low gestational age at birth has a profound and widely distributed effect on the neonatal saliva methylome that is apparent at term equivalent age. Enriched gene ontology terms related to cell–cell contacts reveal pathways that could mediate the effect of early life environmental exposures on development. Finally, associations between differential DNA methylation and image markers of white matter tract microstructure suggest that variation in DNA methylation may provide a link between preterm birth and the dysconnectivity of developing brain networks that characterizes atypical brain development in preterm infants.


Introduction
Preterm birth, defined as birth at ,37 weeks of gestation, affects around 11% of births worldwide 1 and is a leading cause of neurodevelopmental and cognitive problems that extend across the life course. These include autism spectrum disorder, social difficulties, language impairment, attentiondeficit hyperactivity disorder, reduced intelligence quotient, educational underachievement and psychiatric diagnoses. [2][3][4][5][6][7][8][9] The neural phenotypes that underlie long-term functional impairment include diffuse white matter injury and subsequent dysmaturation of white matter and grey matter neuroaxonal structures, collectively termed the 'encephalopathy of prematurity'. 10 A consequence of the encephalopathy is generalized dysconnectivity of developing structural networks, which can be inferred from diffusion MRI (dMRI) and neurite orientation dispersion and density imaging (NODDI) during the neonatal period. [11][12][13][14][15] Specifically, normal maturation is characterized by a reduction in mean diffusivity (MD) and increases in both fractional anisotropy (FA) and neurite density index (NDI) in white matter; but MD is increased, and FA and NDI are decreased in preterm infants at term equivalent age, compared with control infants born at term. 16 These changes reflect an increase in water content and a decrease in white matter organization in preterm infants. The peak width of skeletonized mean diffusivity (PSMD) is a method for histogram-based calculation of MD distribution across the entire white matter skeleton; it provides a measure of generalized white matter microstructure, is robust to scanner variation and is predictive of cognition in later life. [17][18][19] In previous work, we extended the histogram model to neonatal data and included other dMRI and NODDI metrics. We found that PSMD and peak width skeletonized neurite density index (PSNDI) are altered in preterm infants at term equivalent age and that histogram-based measures have utility for investigating upstream determinants of dysmaturation such as systemic inflammation. 20,21 The mechanisms that link the environmental stress of preterm birth with atypical brain development are uncertain. Variation in DNA methylation (DNAm) is a possible mechanism; DNAm is involved in the regulation of gene expression and cell fate during foetal brain development. 22 Alterations in DNAm contribute to the pathogenesis of neurodevelopmental disorders [Rett syndrome, 23 Immunodeficiency, Centromeric region instability, Facial anomalies (ICE) syndrome 24 and Angleman and Prader-Willi syndromes]. 25 There is growing evidence that differential DNAm can mediate the effect of environmental pressures on brain structure and function across the life course. 26,27 The neonatal methylome is sensitive to prenatal factors such as maternal smoking, 28 maternal body mass index, 29 as well as birth weight. 30 It is altered in association with comorbidities of preterm birth, 26,31,32 and there is some evidence for legacy differences in the methylome two decades after preterm birth. 33 A meta-analysis investigating DNAm from umbilical cord blood identified widespread differential methylation associated with GA at birth (across the range 27−42 weeks) affecting 2% of sites and involving both hypo-and hypermethylation, as measured on the Illumina 450k array. 34 However, due to differences in the cellular composition of samples, epigenetic signatures observed in different tissues are likely to be distinct. 35 The main cellular component of saliva, buccal epithelium, may be more representative of the brain than umbilical cord blood because of ectodermal origin. [36][37][38] Here, our first aim was to determine whether low gestational age at birth was associated with variation in the salivary methylome at term equivalent age and to characterize the biological pathways implicated in the DNAm response to preterm birth. Our second aim was to investigate whether the DNAm signal of gestational age explains variance in measures of white matter microstructure at term equivalent age. We tested the hypotheses that low gestational age at birth is associated with widespread differential methylation apparent at the end of neonatal intensive care and that DNAm contributes to variance in peak width skeletonized metrics of white matter connectivity.

Participants
All participants were born at the Royal Infirmary of Edinburgh, UK. Preterm infants, defined as GA of ,33 weeks of gestation and term infants, defined as GA . 37 weeks based on first-trimester ultrasound scan dating, were recruited to the Theirworld Edinburgh Birth Cohort. This is a longitudinal study designed to investigate the effect of preterm birth on brain development. 39

DNA methylation preprocessing
Raw intensity (.idat) files were read into the R environment (version 3.4.4) using minfi. wateRmelon and minfi were used for preprocessing, quality control and normalization. 40,41 The pfilter function in wateRmelon was used to exclude samples with 1% of sites with a detection P-value . 0.05; sites with beadcount ,3 in 5% of samples and sites with 1% of samples with detection P-value . 0.05. Cross hybridizing probes and probes targeting single nucleotide polymorphisms with overall minor allele frequency ≥0.05 were also removed. 42 Control probes were also removed. Samples were removed if there was a mismatch between predicted sex (minfi) and recorded sex (n = 3). Data were danet normalized which includes background correction and dye bias correction. 41 Saliva contains different cells types, including buccal epithelial cells and leucocytes. Epithelial cell proportions were estimated with epigenetic dissection of intra-sample heterogeneity with the reduced partial correlation method implemented in the R package EpiDISH. 43 Probes located on sex chromosomes were removed before analysis. Data from one of each twin pair were removed randomly (n = 20).

MRI acquisition
MRI was obtained at the same appointment as saliva sample collection for DNAm analysis. Structural and dMRI were performed on 93 infants using a MAGNETOM Verio 3T clinical MRI scanner (Siemens Healthcare GmbH, Erlangen, Germany) and 12-channel phased-array head coil, which were used to acquire dMRI using a protocol consisting of 11 T 2 -and 64 diffusion-weighted (b = 750 s/mm 2 ) singleshot spin-echo echo-planar imaging (EPI) volumes acquired with 2 mm isotropic voxels (echo time (TE) = 106 ms and repetition time (TR) = 7300 ms).
One hundred and twenty-one infants were scanned using a MAGNETOM Prisma 3T clinical MRI scanner (Siemens Healthcare GmbH) and 16-channel phased-array paediatric head and neck coil. dMRI was acquired in two separate acquisitions: the first consisted of eight baseline volumes [b = 0 s/mm 2 (b0)] and 64 volumes with b = 750 s/mm 2 ; the second consisted of 8 b0, 3 volumes with b = 200 s/mm 2 , 6 volumes with b = 500 s/mm 2 and 64 volumes with b = 2500 s/ mm 2 . An optimal angular coverage for the sampling scheme was applied. 44 In addition, an acquisition of three b0 volumes with an inverse phase encoding direction was performed. All dMRI volumes were acquired using single-shot spin-echo EPI with 2-fold simultaneous multislice and 2-fold in-plane parallel imaging acceleration and 2 mm isotropic voxels; all three diffusion acquisitions had the same parameters (TR/TE 3500/78.0 ms). Images affected by motion artefact were re-acquired multiple times as required; dMRI acquisitions were repeated if the signal loss was seen in three or more volumes.
Infants were fed and wrapped and allowed to sleep naturally in the scanner without sedation. Pulse oximetry, electrocardiography and temperature were monitored. Flexible earplugs and neonatal earmuffs (MiniMuffs, Natus) were used for acoustic protection. All scans were supervised by a doctor or nurse trained in neonatal resuscitation. Structural images were reported by an experienced paediatric radiologist (A.J.Q.)

dMRI preprocessing
Diffusion images that were acquired on the MAGNETOM Verio scanner were denoised using a Marchenko-Pastur-principal component analysis (PCA)-based algorithm [45][46][47] ; eddy current distortion and head movement were corrected using outlier replacement; [48][49][50] bias field inhomogeneity correction was performed by calculating the bias field of the mean b0 volume and applying the correction to all the volumes. 51 FA and MD were calculated from the dMRI data.
The two dMRI acquisitions from the MAGNETOM Prisma scanner were first concatenated and then denoised using a Marchenko-Pastur-PCA-based algorithm [45][46][47] ; eddy current, head movement and EPI geometric distortions were corrected using outlier replacement and slice-to-volume registration; [48][49][50]52 bias field inhomogeneity correction was performed by calculating the bias field of the mean b0 volume and applying the correction to all the other volumes. 51 From the dMRI data, we calculated the three eigenvalues and eigenvectors of the water diffusion tensor, and NODDI (Bingham distribution) parametric maps using cuDIMOT [intracellular volume fraction (NDI) and the overall orientation dispersion index (ODI TOT )]. 12,13,53 FA and MD were calculated using single-shell data to match the Verio scanner.

The peak width of skeletonized water diffusion parameters
All the subjects were registered to the Edinburgh Neonatal Atlas (ENA50) using DTI-TK. 20 The diffusion tensor derived maps of each subject (FA and MD) were calculated after registration; NDI was then propagated to the template space using the previously calculated transformations. The data were skeletonized using the ENA50 skeleton and then multiplied by a custom mask. Finally, the peak width of the histogram of values computed within the skeletonized maps was calculated as the difference between the 95th and 5th percentiles. 17,20 Statistical analysis

Epigenome-wide association analyses
Unless otherwise stated, analysis was completed in R version 3.4.4. 54 An overview of the analysis pipeline is shown in Supplementary Fig. 1. Surrogate variable analysis (SVA) of the data matrix was carried out, to adjust for potential confounders such as batch, using the sva function in the sva package in R. 55,56 A fully adjusted model was specified before SVA to retain signal explained by biological variables of interest: CpG gestational age at birth + age at scan + birthweight Z-score + maternal smoking + sex + epithelial cells. SVA identified 17 significant surrogate variables (SVs) which were subsequently included in the analysis.
An epigenome-wide association study (EWAS) was performed using the limma package in R. 57 Beta values of each of 776 025 CpG sites were regressed (as dependent variables) on gestational age (GA) at birth using linear regression. Covariates were added to adjust for sex, birthweight Z-score, age at sample collection, maternal smoking, estimated epithelial cell proportions and 17 surrogate variables. A significance threshold of 3.6 × 10 −8 was selected, which represents genome-wide significance. 58

Differentially methylated region analysis
Differentially methylated regions (DMRs) were assessed using the dmrff function in the dmrff package in R. 59 Here, a differentially methylated region is a region containing two or more sites separated by ≤500 bp with EWAS analysis P ≤ 0.05 and methylation changes in a consistent direction. Following dmrff's subregion selection step, DMRs with Bonferroni-adjusted P ≤ 0.05 were significant.

Gene set testing
Gene set enrichment analysis was carried out using the gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) databases, and using the gometh function in missMethyl package, which controls for multiple probe bias. 60 We performed an analysis that included those sites that reached genome-wide significance in EWAS and a second analysis that also incorporated those sites that contributed to differentially methylated regions.

Principal component analysis
Principal component analysis (PCA) was conducted on CpG probes that reached genome-wide significance, using the prcomp function in R. CpGs were precorrected for the effects of biological covariates and surrogate variables via linear regression. The scree plot was visually inspected to select a principal component (eigenvalue . 1) to be carried forward for subsequent analysis.

Linear regression between DNAm and peak width skeletonized metrics
Pearson's correlation coefficient was used to assess the relationship between the first PC identified from PCA and gestational age at birth. This PC was used in linear regression models, as an independent variable, to test the associations between peak width of skeletonised mean diffusivity (PSMD), peak width of skeletonised neurite density index (PSNDI) and peak width of skeletonised fractional anisotropy (PSFA) with DNAm, conducted in R version 4.0.1. 54 In models testing PSMD and PSFA, MRI scanner was included as a binary covariate as MRI data from two scanners were included. PSNDI was only available from data acquired on one scanner, so it was not necessary to include the scanner as a covariate. We report standardized regression coefficients and P-values.

Data availability
The atlas with templates can be found at https://git.ecdf.ed. ac.uk/jbrl/ena and the code necessary to calculate histogrambased metrics is at https://git.ecdf.ed.ac.uk/jbrl/psmd. Requests for original image data will be considered through the BRAINS governance process: www.brainsimagebank.ac. uk. 61 DNAm data available upon request from Theirworld Edinburgh Birth Cohort, University of Edinburgh (james. boardman@ed.ac.uk or https://www.tebc.ed.ac.uk/2019/ 12/data-access-and-collaboration/). The DNAm and metadata are not publicly available due to them containing information that could compromise participant consent.

Cohort
DNAm data were collected from 311 neonates. Twenty-nine did not meet DNAm preprocessing QC criteria and were excluded. One participant with a congenital abnormality was removed, as were three participants whose sex predicted from DNAm data did not match their recorded sex. This group of 278 neonates included 20 sets of twins. After random removal of one twin from each set, there was no evidence of imbalance for birthweight between the twins that were removed and the twins that remained in the sample (t = −0.157, P = 0.88) or for sex (χ 2 = 0.417, P = 0.52).
The study group consisted of 258 neonates: 155 participants were preterm (gestational age range 23.28-34.84) and 103 were infants born at full term (gestational age range 36.42-42.14), see Table 1 for participant characteristics and Supplementary Fig. 2 for participant flow. Among the preterm infants, 38 (25%) had bronchopulmonary dysplasia (defined as need for supplementary oxygen ≥36 weeks gestational age), 9(6%) developed necrotizing enterocolitis requiring medical or surgical treatment and 31 (20%) had an episode of postnatal sepsis defined as either blood culture positivity with a pathogenic organism, or physician decision to treat for ≥5 days in the context of the growth of coagulase-negative staphylococcus from blood or a negative culture. Of the 258 participants with DNAm data, 214 also had MRI data.

Widespread differential saliva DNAm in association with gestational age at birth
We conducted an epigenome-wide association study whereby CpG methylation at 776 025 sites was regressed on gestational age at birth, controlling for birthweight Z-score, infant sex, age at sample collection, maternal smoking, estimated epithelial cell proportion and surrogate variables. The genomic inflation factor was 1.72 ( Supplementary Fig. 3). Differential methylation in relation to gestational age at birth was identified at 8870 CpG sites at genome-wide significance (P , 3.6 × 10 −8 ), Fig. 1. Of these, 4250 (47.9%) sites demonstrated a positive association with gestational age, while 4620 (52.1%) were negatively associated. Following Bonferroni adjustment, 1767 DMRs corresponding to 4664 CpG sites were significant at P , 0.05. Of these, 11 had 10 or more CpG sites contributing to the DMR. The largest DMR mapped to a genomic region that encodes two genes: NNAT and BLCAP. The 29 probes mapped to this region were all located within islands and positively associated with gestational age at birth, indicating that longer gestation corresponds to hypermethylation. Of the 10 most significant differentially methylated probes (DMPs), three probes were localized to the IRX4 gene, one probe to the GAL3ST4 gene, and one to LOXL4 (Table 2; Supplementary Fig. 4). The probes with the largest absolute magnitude effect size (top five hypermethylated and hypomethylated) were mapped to the following genes: IRX2, SMIM2, INTS1, HEATR2, ZBP1 and UBXN11 (Table 3).

Pathways implicated in functional testing
Based on the 8870 sites that met the genome-wide significance threshold (P , 3.6 × 10 −8 ), no KEGG terms remained significant following false discovery rate (FDR) correction for multiple comparisons. Two GO terms were enriched following FDR correction: one for anchoring (GO:0070161; q = 0.0062) and one for adherens junction (GO:0005912; q = 0.0062). In an analysis that incorporated all 11 752 distinct CpGs from both EWAS and DMR analysis, 14 GO terms were enriched (Table 4).

Discussion
By studying a unique database of DNA linked to brain imaging in a population of preterm and term infants, we have identified extensive differential methylation in association with gestational age at birth and revealed an association between the principal axis of methylation and brain dysmaturation within the same sample. Differentially methylated regions and probes were distributed widely across the genome, indicating that gestation duration has a global effect on DNAm. Gene enrichment analysis of changes associated with gestational age identified gene sets pertaining to cell contacts and cytoskeleton. A single principal component that explained 23.5% of the variance in differential DNAm linked to preterm birth was closely associated with markers of generalized dysconnectivity across the white matter skeleton. The data are consistent with studies that have reported associations between length of gestation and genome-wide variation in DNAm within foetal brain 22 and umbilical cord blood; 34,62 and widely distributed variation is reported in association with postmenstrual age at the sampling of preterm infants (a proxy for gestational age at birth). 31 The signature we identified in saliva sampled at term equivalent age included 233 probes that were previously shown to be differentially methylated in association with gestational age in a meta-analysis of umbilical cord blood samples that reported 8899 gestation-dependent CpGs. 34 The limited overlap could be explained by differences in the cellular composition of assessed tissues, different array types used to measure DNAm, or due to the time of sampling. Cord blood is collected at birth and so methylation changes at this time reflect foetal maturity and/or prenatal experience, whereas the methylation signature at term equivalent age reflects the allostatic load of early postnatal experiences, as well as the prenatal environment. We chose to sample at term equivalent age because postnatal co-morbidities of preterm birth and NICU care practices such as painful stress exposures alter DNAm profiles, and because cumulative DNAm variations over this time period may link exposure to behavioural outcome in preterm infants. 32,63,64 Functional analyses of DMPs identified two enriched GO terms, for adherens and anchoring junctions. When distinct probes that contributed to both DMPs and DMRs were combined, GO analysis identified an additional 14 terms related to cell-cell adhesion, cell adhesions with the extracellular matrix and signalling from the extracellular membrane; 12 of these were also identified in the meta-analysis of gestational age effects on DNAm obtained at birth from umbilical cord blood. 34 The most significant DMR mapped to a site encoding two genes: NNAT and BLCAP. NNAT encodes the neural fate initiator neuronatin, the expression of which decreases throughout development; 65 there was a positive association with increasing gestational age at birth. Hypomethylation of NNAT is associated with a corresponding increase in the expression of neuronatin. 66 BLCAP   encodes the bladder cancer-associated protein. This is a tumour suppressor that induces apoptosis, with high expression in brain and B lymphocyte. The candidacy of this locus, encoding NNAT and BLCAP, as a region of interest whose expression may be modified by perinatal exposures is supported by previous EWAS of gestational age, maternal body mass index, maternal smoking and schizophrenia. 22,29,62,[67][68][69][70] In addition, BLCAP was found to have reduced methylation in placental samples from women who had preeclampsia. 71 Probes that demonstrated the largest magnitude of effect in association with gestational age mapped to genes previously associated with gestational age or maternal risk factors in EWAS. 70 Hypermethylated probes were found in genes including ZBP1, which was identified in EWAS investigating gestational age and hypertensive disorders of pregnancy; 62,69,72 INTS1, which has been identified in EWAS of gestational age, hypertensive disorders of pregnancy, maternal body mass index, birthweight and breastfeeding duration. 22,30,62,67,72,73 Hypomethylated probes were found in genes including UBXN11, which was identified in studies of gestational age; 22,62 and IRX4. Three of the 10 most significant DMPs mapped to the IRX4 gene, all of which displayed a negative association with gestational age at birth. IRX4 is associated with cardiac development in vertebrates, including humans. 74 Its homologues have been implicated in retinal axon guidance in zebrafish, and neural patterning in Xenopus, 75,76 and it has been identified in previous EWAS of hypertensive disorders of pregnancy 72 and prenatal maternal stress. 77 The novel pathways and genes implicated by EWAS studies of gestational age could provide a strong empirical basis for the selection of genes in targeted analyses in association with neuroimaging. 27 For example, one of the genes identified in our EWAS has been implicated in neurodevelopmental disorders; biallelic mutations in INTS1 have been associated with a rare neurodevelopmental syndrome characterized by intellectual disability. [78][79][80] We used metrics of generalized white matter connectivity to assess relationships between DNAm and brain development because generalized dysconnectivity of structural networks is a hallmark of preterm brain dysmaturation. 10,14,81 PSNDI and PSMD were strongly associated with the first principal component of gestational age-dependent variation in DNAm but PSFA was not. This suggests that variations in DNAm could contribute to the higher variability in water Figure 2 Scatter plots with regression lines and 95% confidence intervals showing the relationships between gestational age at birth (weeks) and DNAm with PSMD and PSNDI, where peak width is the difference between the 95th and 5th centile of histogram values across the white matter skeleton. The associations between gestational age (weeks) and PSMD and PSNDI are shown in (A) and (B), respectively. The relationships between DNAm PC1 and PSMD and PSNDI are shown in (C) and (D), respectively. PS metrics are residualized for gestational age at scan; PSMD is additionally residualized for scanner variable. content and intra-axonal volume that characterize preterm brain dysmaturation. 20 We have previously reported that differential DNAm is associated with FA of the genu of the corpus callosum and tract shape of the right corticospinal tract; 26 it is most likely that we did not observe an association between gestational age or DNAm with PSFA because this metric is subject to histogram shift, 20 meaning that although there are groupwise differences in FA values across the skeleton, the spread of values is the same. Associations between gestational age at birth and both DNAm and image markers of dysconnectivity, and between DNAm and image features, suggest that differential DNAm contributes, in part, to the relationship between gestational age and brain network dysconnectivity in preterm infants. However, some of the DNAm signatures may be related to underlying causes of prematurity that operate in foetal life, such as infection or preeclampsia. We cannot rule out that such changes are involved in the aetiology of preterm birth, which might preclude them from mediating specific associations with brain development.
The strengths of this study are that we studied a population of preterm and term infants across the gestational age range of 23-42 weeks, who were uniquely phenotyped with DNAm and dMRI. We sampled after the period of NICU care to capture the allostatic load of preterm birth. We measured DNAm from neonatal saliva samples, which has consistency with brain DNAm patterns and is noninvasive. The Illumina EPIC platform provided extensive coverage of the methylome (850 000 sites) and we controlled for cell composition. Finally, we used an image phenotype that is robust to scanner variation. 17 There are some limitations. First, control for cell composition was based on estimation of cell proportions rather than measurement, so we cannot rule out the possibility that some of the signal identified was related to variation in cell composition. Second, mediation analysis to assess causation was not possible because the association between the DNAm PC1 and peak width skeletonized metrics might result from the DNAm PC being derived from CpG sites that are associated with gestational age, so variance attributable to the mediating variable cannot be assumed. This could be addressed by out-of-sample validation, which will require other neonatal cohorts with both saliva methylome and dMRI data. Cohorts, with such data, that recruit neonates from across the gestational age continuum would also provide valuable replication of analyses described here. There was also some evidence of inflation of test statistics based on the genomic inflation factor. However, the genomic inflation factor is thought to provide an overestimate of inflation and corrections based on it may   be overly conservative. 82 The value of the genomic inflation factor was also similar to those previously reported in neonatal cohorts. 34 In addition, visual inspection of our results via Manhattan plots suggests that our finding, of widespread differences in DNAm in relation to gestational age at birth, is in line with previous studies that have investigated this in cord blood and foetal brain tissue. 22,34 Finally, preterm birth itself can reflect maternal health or circumstance; in this cohort, mothers of preterm infants were more likely to have pregnancy-induced hypertension and the control mothers were more likely to have college or university qualifications. It is possible that these, or other maternal factors associated with preterm birth such as infection, chronic disease, poor nutrition and alcohol, tobacco or drug use, could influence the neonatal methylome at term equivalent age. A large prospective study with longitudinal sampling from birth would be required to test the hypothesis that maternal factors contribute to changes in the methylome apparent after neonatal intensive care.
In conclusion, gestational age at birth has a profound and widely distributed impact on the neonatal saliva methylome at term equivalent age, which reflects the allostatic load of preterm birth itself and postnatal exposures during neonatal intensive care. GO terms related to cell-cell contacts were enriched, indicating that cell contacts and organization are implicated in the phenotype. Associations between DNAm and PSMD and PSNDI suggest that variations in DNAm could contribute to white matter dysconnectivity commonly seen in preterm infants, and this analysis identified several genes and gene regions that could provide further insight into the molecular mechanisms by which early exposure to extrauterine life influences neurodevelopment.