Prenatal Lead (Pb) Exposure and Peripheral Blood DNA Methylation (5mC) and Hydroxymethylation (5hmC) in Mexican Adolescents from the ELEMENT Birth Cohort.

BACKGROUND: Gestational lead (Pb) exposure can adversely affect offspring health through multiple mechanisms, including epigenomic alterations via DNA methylation (5mC) and hydroxymethylation (5hmC), an intermediate in oxidative demethylation. Most current methods do not distinguish between 5mC and 5hmC, limiting insights into their individual roles.
OBJECTIVE: Our study sought to identify the association of trimester-specific (T1, T2, T3) prenatal Pb exposure with 5mC and 5hmC levels at multiple cytosine-phosphate-guanine sites within gene regions previously associated with prenatal Pb (HCN2, NINJ2, RAB5A, TPPP) in whole blood leukocytes of children ages 11-18 years of age.
METHODS: Participants from the Early Life Exposure in Mexico to Environmental Toxicants (ELEMENT) birth cohorts were selected (n=144) for pyrosequencing analysis following oxidative or standard sodium bisulfite treatment. This workflow directly quantifies total methylation (5mC+5hmC) and 5mC only; 5hmC is estimated by subtraction.
RESULTS: Participants were 51% male, and mean maternal blood lead levels (BLL) were 6.43±5.16μg/dL in Trimester 1 (T1), 5.66±5.21μg/dL in Trimester 2 (T2), and 5.86±4.34μg/dL in Trimester 3 (T3). In addition, 5hmC levels were calculated for HCN2 (mean±standard deviation(SD), 2.08±4.18%), NINJ2 (G/C: 2.01±5.95; GG: 0.90±3.97), RAB5A (0.66±0.80%), and TPPP (1.11±6.67%). Furthermore, 5mC levels were measured in HCN2 (81.3±9.63%), NINJ2 (heterozygotes: 38.6±7.39%; GG homozygotes: 67.3±9.83%), RAB5A (1.41±1.21%), and TPPP (92.5±8.03%). Several significant associations between BLLs and 5mC/5hmC were identified: T1 BLLs with 5mC in HCN2 (β=-0.37, p=0.03) and 5hmC in NINJ2 (β=0.49, p=0.003); T2 BLLs with 5mC in HCN2 (β=0.37, p=0.03) and 5hmC in NINJ2 (β=0.27, p=0.008); and T3 BLLs with 5mC in HCN2 (β=0.50, p=0.01) and NINJ2 (β=-0.35, p=0.004) and 5hmC in NINJ2 (β=0.45, p<0.001). NINJ2 5mC was negatively correlated with gene expression (Pearson r=-0.5, p-value=0.005), whereas 5hmC was positively correlated (r=0.4, p-value=0.04). DISCUSSION: These findings suggest there is variable 5hmC in human whole blood and that prenatal Pb exposure is associated with gene-specific 5mC and 5hmC levels at adolescence, providing evidence to consider 5hmC as a regulatory mechanism that is responsive to environmental exposures. https://doi.org/10.1289/EHP8507.

Introduction

The Developmental Origins of Health and Disease (DOHaD) hypothesis (Barker 2007) postulates that in utero exposures, including exposure to environmental contaminants such as lead (Pb), can permanently modify an organism’s molecular biology, physiology, and metabolism, potentially leading to myriad effects on cognition, growth, maturation, and metabolic risk (Bellinger et al. 2016). Pb is a widely abundant environmental pollutant known to be a potent neurotoxicant, even at low levels. Pregnancy is a vulnerable time period for Pb exposure because both current and past maternal Pb exposure may affect the developing fetus. Animal and human studies have provided evidence for the impact of early-life Pb exposure on the epigenome (Faulk et al. 2013; Montrose et al. 2020; Goodrich et al. 2015; Rygiel et al. 2020; Wu et al. 2017; Dou et al. 2019).

Epigenetic modifications are mitotically heritable molecular changes that regulate gene expression without altering the underlying DNA sequence. DNA 5-methylcytosine (5mC) is the addition of a methyl group covalently bound to the 5′-carbon of Cytosine; in mammals this typically occurs on a Cytosine adjacent to a Guanine, referred to as a CpG site (Illingworth et al. 2010; Ohlsson and Kanduri 2002). During demethylation, 5mC is oxidized into 5-hydroxymethylcytosine (5hmC), an intermediate of oxidative demethylation that can remain as a stable modification (Sadakierska-Chudy et al. 2015; Bachman et al. 2014). 5hmC perturbations in the brain have been reported in several early-life neurological disorders and later-life neurodegenerative disorders in both human and animal studies (Zhao et al. 2017; Wang et al. 2013; Zhubi et al. 2014; Chouliaras et al. 2013). Levels of 5hmC vary with tissue, where brain contains the highest levels and blood contains the lowest levels, but 5hmC is still detectable in blood at modest levels (Kochmanski et al. 2018a, 2018b; Globisch et al. 2010; Li and Liu 2011; Nestor et al. 2012). Recent studies suggest that the functional roles of 5hmC are distinct from 5mC (López et al. 2017). 5mC and 5hmC undergo dynamic changes during early gestation that disseminate through mitosis to new cells and developing organs, potentially persisting throughout the life span. These modifications are important for gene regulation, including in the nervous system, and have implications for learning and memory from early development, adolescence, and into adulthood (Jobe and Zhao 2017; Spiers et al. 2017; Zhao et al. 2017; Vogel Ciernia and LaSalle 2016). Most current methods for quantifying DNA methylation, including the gold standard sodium bisulfite sequencing, collectively measure 5mC and 5hmC without distinguishing between the two. This approach has limited our ability to identify whether environmental exposures alter 5mC, 5hmC, or both and what the implications are for early life neurodevelopment and long-term health.

Pregnancy is a key starting point to investigate the association between environmental contaminants and epigenetic modifications leading to health outcomes in offspring. An important mechanism by which developmental exposures can affect long-term disease risk is through the disruption of normal epigenetic processes, thereby affecting gene regulation and subsequent chronic outcomes (Dolinoy and Jirtle 2008). A study prenatally exposing dams, which are genetically invariant mice 93% identical to C57BL/6J strain, to physiologically relevant doses of Pb (, , and in water) 2 wk prior to mating through lactation until weaning at postnatal day 21 assessed DNA total methylation via pyrosequencing at four Intracisternal A particle (IAP) elements in the brain (Montrose et al. 2017). IAPs are a class of murine retrotransposons that are environmentally responsive (Waterland and Jirtle 2004). Prenatal Pb exposure reduced DNA methylation at three of the IAPs in the brain with dose–dependent and sex-specific effects in comparison with control mice. A study exposing rodents to 0.2% Pb acetate in water postnatally reported hypomethylation and initial transient 8-hydroxy-2′-deoxyguanosine () accumulation in the cerebral cortex, a neurodegeneration biomarker (Bolin et al. 2006). A U.S. prospective human pregnancy cohort, Project Viva, with low Pb exposure (averaging in erythrocytes) conducted an epigenome-wide analysis on 268 umbilical cord blood samples to evaluate the association between maternal Pb exposure and DNA methylation and identified sex-specific differentially methylated CpG sites (Wu et al. 2017). Using a similar epigenome-wide approach but with a moderately to highly exposed population, we identified differentially methylated genes in 89 umbilical cord blood DNA samples from the Early Life Exposure in Mexico to Environmental Toxicants (ELEMENT) cohort that varied by trimester of exposure (Rygiel et al. 2020). Although these and several other studies provide evidence for associations between Pb exposure and DNA methylation profiles in humans and animals, only total methylation was assessed.

Only one study has considered intermediates (5hmC) in adverse environmental associations with prenatal Pb, which may have implications for understanding the association between prenatal exposures and adverse health outcomes. That study provided some evidence for associations between prenatal Pb and 5hmC in umbilical cord blood at the regional level through a modified epigenome-wide method called hMeDIP-450K chip (Sen et al. 2015). The investigators identified both sex-independent and sex-specific differentially methylated and hydroxymethylated regions, where sex-dependent associations were more common in 5mC in comparison with 5hmC. There is a mechanistic basis by which Pb may affect 5hmC levels. Pb-induced oxidative stress results in the accumulation of () (Tretter and Adam-Vizi 2005), a cofactor for ten-eleven translocation (TET) enzymes, which are involved in the oxidation of 5mC to 5hmC (Chia et al. 2011; Coulter et al. 2013). Thus, Pb may increase activity of TET enzymes and increase 5hmC across the genome (Chia et al. 2011).

We hypothesized that prenatal Pb exposure would alter epigenetic programming of 5mC and 5hmC within genes involved in neurological function and that this can be detected in samples from adolescents. We leveraged the ELEMENT longitudinal birth cohort with rich data on prenatal Pb exposure and offspring follow-up through adolescence, including whole blood leukocyte DNA and prenatal Pb exposure biomarkers, to investigate whether Pb is associated with programming of neurocognitive-related genes between the ages of 11–18 y in offspring whole blood leukocytes. We used oxidative-bisulfite (oxBS) pyrosequencing to profile both 5mC and 5hmC in adolescent blood leukocyte DNA at four neurocognitive-related genes—HCN2, NINJ2, RAB5A, and TPPP—that were associated with prenatal Pb exposure in ELEMENT in our previous epigenome-wide study of total methylation in umbilical cord blood samples (Rygiel et al. 2020).

Methods

Study Population

Participants are from the second and third cohorts of the ELEMENT project, a series of longitudinal birth cohorts originally designed to investigate the influence of Pb exposure—in utero and in childhood —on child growth and neurodevelopment. Mothers were recruited between 1997 and 2000 (second cohort) and 2001–2005 (third cohort) from the Mexican Social Security Institute hospital in Mexico City. Eligibility and exclusion criteria are as previously described (Perng et al. 2019; Hu et al. 2006). Data collected included sex, gestational age, household socioeconomic status, anthropometry, and other environmental exposures at multiple follow-up visits from participants’ infancy through adolescence. For the current study, 144 of the ELEMENT participants who were followed-up through the 2015 study visit and who provided an adolescent blood sample for DNA isolation were selected (Figure 1). These participants with prenatal and adolescent Pb measures and/or previous cord blood DNA methylation analysis were prioritized for the adolescent epigenetic analysis (Rygiel et al. 2020).

Figure 1.

Data collection timeline from the ELEMENT birth cohort in Mexico City, Mexico, between 1997–2000 and 2001–2005 for pregnancies and through adolescence into 2018. Maternal whole blood samples were collected during the first trimester (T1), second trimester (T2), and third trimester (T3) and analyzed for blood Pb concentrations. Maternal bone Pb was measured 1 month postpartum as an indicator of cumulative exposure over the course of the gestational period. Adolescent whole blood samples were collected in offspring at a follow-up visit occurring once between the ages of 11 and 18 y for Pb measures and DNA 5mC and 5hmC profiling. Covariate data including demographics and anthropometry were obtained at each stage with sample collection. Note: BLL, blood lead levels.

Prior to participation, all mothers were informed about the study; those who agreed to participate read and signed a letter of informed consent about the original study. Children also provided informed assent, and their mothers provided informed consent prior to participation in the adolescent follow-up visits. The research protocol and all amendments to the study were approved by the ethics committees of the National Institutes of Public Health of Mexico, participating hospitals, and the internal review boards at participating institutions, including the University of Michigan.

Pb Exposure Assessment and Genomic DNA Isolation

Cohort 2 maternal venous blood lead levels (BLLs) from each trimester were measured using inductively coupled plasma mass-spectrometry (ICP-MS; Thermo Finnigan) at the University of California, Santa Cruz, as described previously (Lamadrid-Figueroa et al. 2006). Cohort 3 trimester-specific maternal BLLs were measured using graphite furnace atomic absorption spectrometry (instrument model 3000; PerkinElmer) at the Trace Metal Laboratory of the American British Cowdry Hospital. Perinatal maternal bone Pb levels were also measured in the left patella (trabecular bone) and mid-shaft of the left tibia (cortical bone) 1 month postpartum as an indicator of cumulative Pb exposure during pregnancy using a spot-source K-shell X-ray fluorescence (K-XRF) instrument. The technical specifications and validation of this instrument are described in detail elsewhere (Aro et al. 1994). Tibia and patella bond Pb measures were dropped if their associated uncertainty levels were greater than () and (), respectively. Next, any negative tibia () and patella () bone Pb measurements were recoded as positive values with random numbers in a uniform distribution between 0 and the limit of detection. Whole blood samples were collected during the 2015 follow-up visit conducted when the children were between the ages of 11 and 18 y old, and blood was stored frozen at . DNA was isolated from blood leukocytes using Qiagen kits and standard protocols for blood DNA isolation. Nucleic acid yield and purity were assessed first using a NanoDrop spectrophotometer (ThermoFisher Scientific), and double-stranded DNA was also quantified via a Qubit fluorometer. All DNA samples were stored at .

Candidate Gene Selection

CpG sites that mapped back to four genes relevant to neurological function were selected as candidates from an epigenome-wide study of prenatal Pb exposure conducted with a subset of ELEMENT children in umbilical cord blood (Rygiel et al. 2020) and adolescent whole blood DNA (unpublished): HCN2 (cg06657917), NINJ2 (cg19692784, cg14911689, cg05578102), RAB5A (cg17138393), and TPPP (cg25353752; probe IDs from the Infinium MethylationEPIC BeadChip). RAB5A DNA methylation in umbilical cord blood (Rygiel et al. 2020) and TPPP DNA methylation (unpublished) in adolescent whole blood were inversely associated with T1 maternal BLLs and were selected for the current targeted analysis. HCN2 and NINJ2 were selected because T1 BLLs were associated with greater than 5% change in methylation (unpublished; hypermethylated in HCN2 and hypomethylated in NINJ2) in both umbilical cord blood and adolescent whole blood leukocytes, and these genes have been previously shown to be involved in early-life neurological and/or cognitive development, for which mutations or disruptions in function have been associated with neuronal activity (Zhong et al. 2018), and neurite growth and regeneration (Araki and Milbrandt 2000). EXT1 and LRFN1 were additionally reported as significantly associated with prenatal Pb exposure biomarkers in the original umbilical cord blood epigenome-wide study but were not included in this study because a) we were unable to design primers to amplify the desired region in EXT1 and b) LRFN1 was not associated with first trimester BLLs, which was the exposure time period of most interest.

DNA 5mC and 5hmC Quantification

Pyrosequencing primers for HCN2, NINJ2, RAB5A, and TPPP were designed using the PyroMark Assay Design Software 2.0 Methylation Analysis (CpG) Assay (Table S1). Primers were designed to target a specific CpG site within a region reported to be differentially methylated by first trimester Pb exposure in umbilical cord blood and/or adolescent whole blood DNA (see “Candidate Gene Selection” section). ELEMENT genomic DNA samples were oxidative and/or bisulfite treated according to the NuGen TrueMethyl oxBS Module protocol (Cat. No. 0414-32) and Zymo EZ DNA Methylation kit (Cat. No. D5003). Briefly, of input genomic DNA was dissolved in nuclease-free water to , and each genomic DNA sample was divided into two aliquots (Booth et al. 2013). Each aliquot underwent independent, parallel treatments and were either oxidative bisulfite–converted with the NuGen TrueMethyl oxBS Module or Zymo EZ DNA Methylation kit. The yield and purity of treated samples were quantified using a NanoDrop spectrophotometer.

The target loci within HCN2, NINJ2, RAB5A, and TPPP were amplified in both bisulfite and oxidative bisulfite–converted samples. Polymerase chain reaction (PCR) products were verified using the QIAxcel automated DNA electrophoresis. DNA methylation levels were quantified using the PyroMark Q96 ID instrument (Qiagen). Targeted pyrosequencing captured 6 CpG sites within a CpG island of HCN2, 8 within a DNase hypersensitivity region of NINJ2, 11 within a DNase hypersensitivity region in the first exon of RAB5A, and 4 within the third exon of TPPP (Table S2). The sixth CpG site covered by the NINJ2 assay was dropped because 90% of samples failed at this location, leaving seven CpG sites for NINJ2 from downstream analysis. Sequencing the bisulfite converted samples quantifies the total level of , whereas sequencing the oxidative bisulfite–treated samples quantifies total levels of 5mC. Thus, 5hmC levels were quantified by subtracting the results from the oxidative bisulfite–converted samples (5mC) from the bisulfite-converted sample (). It should be noted that because 5hmC is based on a calculation, the difference can sometimes be negative at sites with zero or low levels of methylation as a consequence of random noise (Field et al. 2015; Stewart et al. 2015). A threshold for calling 5hmC was set to , which corresponds to the 95th percentile of all negative 5hmC values observed across all samples and sites. It should be noted that fewer than 2.0% of all CpG sites were less than , for which those values were dropped from subsequent analyses. For quality control, plates were run with Qiagen EpiTect bisulfite-converted unmethylated (0%) and methylated (100%) human methylation standards (Cat. No. 59665 and No. 59655). OxBS triplicates of eight samples were included for quality control. Standard deviations of triplicate measures following OxBS treatment of NINJ2 (variably methylated) and TPPP (highly methylated) for 5mC averaged 2.11 and 2.28, respectively, at each CpG site captured. Measures of 5mC were precise, with coefficient of variation, whereas the measures of methylation using the standard bisulfite method were more precise, with coefficient of variation. DNA samples were randomized across experimental batches that consisted of four plates for each gene. Paired oxBS and BS conversion samples were always in the same batch. All data that failed internal quality control checks within the PyroMark software were excluded from analysis (HCN2: 5mC , 5hmC ; NINJ2: 5mC , 5hmC ; RAB5A: 5mC , 5hmC ; and TPPP: 5mC , 5hmC ).

Estimates of cell-type composition [ and T lymphocytes, B cells, natural killer (NK) cells, monocytes, granulocytes] for each sample were performed using an established method based on adult cell type–specific differentially methylated regions using data from the Infinium EPIC array for each sample (Houseman et al. 2012).

A methylation quantitative trait locus (MeQTL) was identified within NINJ2 and was found to be correlated with a single-nucleotide polymorphism (SNP) (C/G; rs34038797) within the pyrosequenced region (Gaunt et al. 2016). Genotyping was performed on genomic DNA from all 144 samples using the PyroMark Q96 ID instrument (Qiagen). Pyrosequencing SNP primers for NINJ2 were designed using the PyroMark Assay Design Software 2.0 Genotyping (SNP) Assay (dbSNP, https://www.ncbi.nlm.nih.gov/snp/) (Table S1). Primers were designed to target the rs34038797 SNP.

Gene Expression Analysis via RNA Sequencing

For a subset of participants (), next-generation sequencing of RNA (RNA-seq) was conducted to obtain relative expression data for all genes, and we used data from the four genes that are the focus of this paper. Following collection of whole blood into EDTA-containing tubes, white blood cells were extracted by centrifugation, preserved in RNALater, and stored frozen () until further processing. RNA was isolated via the All-Prep kit (Qiagen). Quality and quantity were assessed via a Bioanalyzer Tapestation (Agilent). Libraries were prepared with Universal Plus mRNA-Seq with Human globin AnyDeplete (NuGEN Technologies, Inc.), which removes globin transcripts that are highly abundant in blood samples. Library preparation and sequencing were performed at the University of Michigan Advanced Genomics Core. Paired-end 50 cycle sequencing on an Illumina HiSeq 4000 was performed. Quality of the raw reads data for each sample was checked using FastQC (version 0.11.3, Babraham Bioinformatics). The Tuxedo Suite software package was used for alignment (Trapnell et al. 2009, 2013; Langmead et al. 2009). Briefly, reads were aligned to the reference mRNA transcriptome (hg19) using TopHat (version 2.0.13) and Bowtie2 (version 2.2.1., both from John Hopkins University, Center for Computational Biology) followed by a second round of post-alignment quality control in FastQC, which allows at most three mismatched values. One sample was dropped due to low alignment rates. All samples used in downstream analysis () had at least good quality aligned reads with alignment rates averaging 60%. Prior to analysis, read counts were normalized by the trimmed mean of M-values method (Robinson and Oshlack 2010).

Statistical Analysis

All statistical analyses were performed in R (version 3.6.1; R Foundation for Statistical Computing). Summary statistics were first calculated. Pb variables analyzed include maternal BLLs at each trimester (T1, T2, T3), bone Pb levels in maternal patella and tibia, and BLLs measured at the follow-up visit, which were all treated as continuous variables. We performed a Wilcox signed-rank test to compare the current subset of ELEMENT with the entire ELEMENT population from the same cohorts (2 and 3). Univariate analyses between Pb exposure variables and covariates [e.g., Pearson correlation coefficient analysis: cell-type proportions, maternal age at birth, adolescent age at time of sample collection, height-for-age z-score, body mass index (BMI)-for-age z-score, weight, education; -test: adolescent sex; analysis of variance (ANOVA): socioeconomic status] were performed (Table S3A–C). We also estimated Pearson correlations between DNA methylation and cell-type proportions [T lymphocytes (CD4T, CD8T), B cells, NK cells, monocytes, granulocytes], where cell-types were not correlated with either DNA 5mC or 5hmC (Table S3D). Potential confounding variables that were significantly associated with both trimester-specific Pb–adolescent age at follow-up and adolescent BLL were identified for inclusion in final statistical models of site-specific DNA methylation data. Although not associated with Pb, sex was also included in statistical models because of its biological effect on DNA methylation. Experimental plate (i.e., batch) was also added as a covariate because of its effect on the technical measurement of DNA methylation. A random intercept for each participant was included in all models to account for autocorrelation from matched 5mC and 5hmC percentages for each individual. Last, in models of NINJ2 methylation, C/C individuals had 0% regional methylation as expected, because the SNP converts the CpG site to a CpC site, resulting in the inability to methylate the locus and subsequent regional methylation changes. Thus, C/C individuals () were dropped from analysis and genotype (G/G or G/C) was included as a covariate.

A mixed-effects model was run, treating 5mC and 5hmC values as repeated measures of a single outcome variable—DNA methylation, given that 5mC and 5hmC are both measured at each CpG site, and the values for these two marks are dependent on each other both biologically and statistically (Kochmanski et al. 2019). To determine whether Pb exposure modifies the balance between the DNA modification categories (5mC and 5hmC), an interaction term between prenatal Pb measures and a dichotomous variable signifying whether the outcome measure was 5hmC or 5mC were included in the statistical model. Because there was evidence for an interaction between Pb and type of DNA methylation, we next used separate mixed-effects regression models treating CpG sites as repeated measures to estimate associations between prenatal Pb and 5mC and associations between prenatal Pb and 5hmC, separately. To assess potential sex-specific effects on the association between prenatal Pb and DNA 5mC and 5hmC, identical models were run as described above but stratified by sex. The lme4 and lmerTest packages within the statistical program R were used for these analyses, and a was considered significant. Finally, RNAseq normalized read counts of HCN2, NINJ2, RAB5A, and TPPP were log-transformed and compared with 5mC% and 5hmC% from pyrosequencing. We used Pearson’s product-moment correlation to quantify the strength of the relationship between these expression data 5mC or 5hmC. were considered significant. The ggplot2 R package was used to plot 5mC and 5hmC results.

Due to outliers in many of the Pb biomarkers, a sensitivity analysis was performed. In this analysis, we reran the models excluding outliers for BLLs or bone Pb measures. Outliers were defined as from the mean. We excluded 4, 5, 1, 3, and 5 outliers that were all 3 SD greater than the means for T1 BLL, T2 BLL, T3 BLL, tibia bone and patella bone Pb, respectively. We compared results with and without the outliers.

Results

Population Parameters and Phenotypic Data

Among the 144 children in the analytic sample, 73 (51%) were male (Table 1). The mean gestational age at birth was 39.1 () wk, with the minimum being 36.0 wk and the maximum being 42.0 wk. The mean age of the offspring was 14.0 () y, ranging from 11.1 y to 17.7 years of age. At the adolescent follow-up visit, mean weight was 53.2 () kg. The height for age z-score and the BMI for age z-score were () and 0.51 (), respectively. Among the 144 mothers, the mean age at offspring birth was 26.7 () y. About 77% of the mothers were from households of low-middle, middle, or middle-high socioeconomic status. Mean maternal education was 11.0 () y, with the minimum being 3.0 y and the maximum being 21.0 y. We compared mother–offspring pairs included in this analysis with all ELEMENT mother–offspring pairs from cohorts 2 and 3 and found that demographic characteristics and Pb biomarker concentrations were not statistically different between the subset and entire ELEMENT population, with the exception of gestational age and adolescent age at follow-up (Table 1).

Table 1

Characteristics of ELEMENT Mexican mother–offspring pairs with adolescent blood leukocyte DNA methylation data compared to all ELEMENT mother–infant pairs.

Characteristics	No.	ELEMENT Subsetb				All ELEMENTc
		Mean±SD or n (%)	Median	Min	Max	No.	Mean±SD or n (%)	Min	Max	p-Valuea
Mothers
Age at birth (y)	144	26.7±5.42	26.5	14.0	39.0	1458	26.1±5.35	14.0	44.0	0.20
Blood lead (μg/dL)
First trimester	127	6.43±5.16	5.30	0.90	35.8	594	5.77±4.03	0.00	36.0	0.34
Second trimester	130	5.66±5.21	4.40	0.00	38.2	616	5.20±3.93	0.00	38.0	0.43
Third trimester	131	5.86±4.34	4.90	0.00	34.0	575	5.54±4.13	0.00	38.0	0.29
Average all trimesters	133	6.02±4.41	5.01	0.43	33.1	643	5.53±3.51	0.00	33.0	0.41
Bone lead (μg/g)
Patella	140	11.4±9.45	9.90	0.11	48.0	1252	11.0±8.64	0.11	68.0	0.79
Tibia	122	9.19±6.88	8.28	0.12	19.0	1104	9.20±7.23	0.00	44.0	0.79
Education (y)	144	11.0±2.90	12.0	3.00	21.0	1212	10.8±3.01	0.00	24.0	0.40
Socioeconomic statusd	144					554				0.13
High		13 (9%)					35 (6%)
Medium/high		26 (18%)					95 (17%)
Medium		50 (35%)					165 (30%)
Low/medium		35 (24%)					178 (32%)
Low		6 (4%)					23 (4%)
Very low		14 (10%)					58 (10%)
Adolescent children
Gestational age at birth (wk)	143	39.1±1.15	39.0	36.0	42.0	1196	38.8±1.57	27.0	42.0	0.02a
Age at follow-up (y)	144	14.0±1.96	13.7	11.0	18.0	549	14.53±2.10	10.7	18.1	0.002a
Male sex	144	73 (51%)				546	268 (49%)			0.80
Blood Pb (μg/dL)	144	3.29±4.44	2.15	0.00	41.0	404	3.09±3.25	0.00	41.0	0.54
Weight (kg)	144	53.2±12.2	53.0	25.9	95.0	554	54.9±13.3	25.9	109.0	0.27
Height for age z-score	144	−0.23±0.88	−0.13	−2.30	−2.06	546	−0.29±0.93	−2.84	3.32	0.43
BMI for age z-score	144	0.51±1.25	0.46	−2.57	3.40	546	0.50±1.25	−3.81	3.40	0.95

Note: BLL, blood lead levels; BMI, body mass index; Max, maximum; Min, minimum; SD, standard deviation; SES, socioeconomic status.

using Wilcox signed-rank test comparing entire ELEMENT cohort to the subset used in this current research. Here we are comparing their characteristics to those of all ELEMENT mother–child pairs from the same cohorts (2 and 3) with available data for each variable.

Study population (): full-term singleton live births in Mexico City, Mexico, from Cohorts 2 (1997–2000) and three (2001–2005) to mothers who have at least one measure of prenatal Pb, archived whole blood samples for DNA methylation analysis, and with nonmissing covariates (i.e., sex, adolescent BLL, and adolescent age).

ELEMENT longitudinal birth cohort (): mother–infant pairs sequentially recruited at pregnancy or delivery from maternity hospitals in Mexico City, Mexico, between 1997–2000 (Cohort 2) or 2001–2005 (Cohort 3).

The AMAI scale (Asociación Mexicana de Agencias de Investigación de Mercado) was used as an indicator of SES.

The mean maternal BLLs averaged over all three trimesters was 5.98 , with a first trimester (T1) mean of 6.43 , second trimester (T2) mean of 5.66 , and a third trimester mean of 5.86 . Maternal BLLs between the trimesters were highly correlated (, ) according to a Spearman’s rank order correlation test. Bone Pb measures represent long-term Pb accumulation because the half-life of Pb in bone can be decades. Average Pb concentration measured in 1-month postpartum maternal patella was 11.4 , and tibia was 9.19 . Current adolescent BLL was 3.29 , which was moderately correlated (, ) with maternal trimester-specific Pb exposures. Genotyping of study participants at rs34038797 SNP revealed 27 individuals with G/G genotype, 83 G/C, and 34 C/C. Participants included in this analysis were compared with all ELEMENT mother–infant pairs from the same cohorts (2 and 3). Characteristics of offspring and Pb biomarker concentrations were not statistically different between the subset included in this analysis and the entire sample (Table 1), with the exception of adolescent age at follow-up, which was slightly lower in the subset.

5mC and 5hmC within Candidate Regions

To measure the levels of DNA 5mC and 5hmC at HCN2, NINJ2, RAB5A, and TPPP loci, we performed parallel BS and oxBS conversion of DNA samples from the adolescent ELEMENT samples. Both 5mC and 5hmC were measured at each gene loci in all adolescent whole blood DNA samples. The mean percentages of 5hmC across genetic regions were measured at 2.00% () in HCN2, 2.01% () for G/C samples in NINJ2 and 0.90% () in G/G NINJ2, 0.65% () in RAB5A, and 1.11% () in TPPP (Table S2). Point measurements of 5hmC in blood at each locus were detectable but much lower than for 5mC, with 13 participants displaying average 5hmC measures at or above 10% in at least one gene. 5mC was measured at 81.3% () in HCN2, 38.6% () for G/C samples NINJ2 and 67.3% () in G/G NINJ2, 1.41% () in RAB5A and 92.5% () in TPPP.

Association of Prenatal Pb with 5mC and 5hmC

Next, an interaction term was used to determine whether prenatal Pb exposure is associated with 5mC and 5hmC in the same or a different manner. In a repeat measures model with the paired 5mC and 5hmC measures at each CpG site in each gene as outcomes, several interaction terms were statistically significant. Because there was evidence for differences in the association between Pb and DNA methylation by type of methylation at several genes, the association between each Pb exposure biomarker and 5mC or 5hmC at each gene was tested separately. In one example, the 5mC and 5hmC proportions in HCN2 shift toward increased 5mC at higher Pb exposure measures when analyzing T1 BLLs, T2 BLLs, T3 BLLs, and patella bone Pb, which was reflected in the significant positive interaction terms (p 0.001). In the models stratified by type of methylation, for each unit increase () in maternal BLL during T1, there was a 0.43 [95% confidence interval (CI): 0.07, 0.78] increase in percent 5mC in HCN2 (), and T2 BLLs demonstrated a positive association with 5mC levels within HCN2 [ (95% CI: 0.04, 0.69); ] (Table 2). Similarly, T3 maternal BLLs were associated with an increase in HCN2 5mC [ (95% CI: 0.10, 0.90); )] However, patella Pb and tibia Pb were not significantly associated with 5mC or 5hmC in HCN2.

Table 2

Results from three different models assessing associations between prenatal Pb biomarkers and DNA methylation. For each Pb biomarker and gene combination, the first two rows display the effect estimates for Pb and for an interaction term between the prenatal Pb biomarker and a dichotomous variable signifying whether the outcome measure is for 5mC or 5hmC; these estimates are from a mixed-effects model with repeated measures of 5hmC and 5mC at multiple CpG sites for each gene as the outcome variables. The third row is a mixed-effects model of 5mC data only; in this model outcomes are repeat measures of 5mC at each CpG site within the gene, and the effect estimate for the Pb biomarker is shown. The fourth row is a model of 5hmC measures only. All three models adjusted for the following covariates: sex, adolescent age, adolescent BLL, and batch. In a sensitivity analysis, the same models for 5mC only and 5hmC only were run with outliers of Pb exposure excluded.

	All observations				Outliers excluded
	n	Betab (95% CI)	p-Value	p-Intc	n	Betab (95% CI)	p-Value
HCN2
T1 Pb
Main effect	120	0.05 (−0.12, 0.21)	0.58			NA	NA
Interaction	120	0.30 (0.16, 0.45)		<0.001		NA	NA
5mC	124	0.43 (0.07, 0.78)	0.02		122	0.58 (0.01, 1.06)	0.02
5hmC	120	−0.04 (−0.18, 0.12)	0.65		118	−0.03 (−0.23, 0.18)	0.79
T2 Pb
Main effect	123	0.06 (−0.13, 0.24)	0.54			NA	NA
Interaction	123	0.28 (0.12, 0.45)		0.001		NA	NA
5mC	127	0.37 (0.04, 0.69)	0.03		123	0.41 (−0.16, 0.97)	0.16
5hmC	123	−0.10 (−0.29, 0.09)	0.30		119	0.01 (−0.28, 0.25)	0.93
T3 Pb
Main effect	124	0.06 (−0.12, 0.25)	0.50			NA	NA
Interaction	124	0.32 (0.16, 0.49)		<0.001		NA	NA
5mC	128	0.50 (0.10, 0.90	0.01		127	0.39 (−0.06, 0.83)	0.09
5hmC	124	−0.06 (−0.16, 0.05)	0.55		123	0.01 (−0.20, 0.22)	0.93
Tibia
Main effect	116	−0.03 (−0.15, 0.08)	0.55			NA	NA
Interaction	116	0.04 (−0.078, 0.15)		0.54		NA	NA
5mC	120	0.03 (−0.21, 0.26)	0.83		117	−0.02 (−0.29, 0.25)	0.88
5hmC	116	−0.06 (−0.16, 0.05)	0.32		113	−0.07 (−0.20, 0.06)	0.29
Patella
Main effect	133	−0.01 (−0.09, 0.06)	0.72			NA	NA
Interaction	133	0.13 (0.05, 0.21)		0.001		NA	NA
5mC	137	0.10 (−0.06, 0.26)	0.21		133	0.12 (−0.08, 0.32)	0.25
5hmC	133	−0.01 (−0.09, 0.06)	0.76		129	−0.01 (−0.10, 0.00)	0.91
NINJ2a
T1 Pb
Main effect	92	0.41 (−0.02, 0.84)	0.06			NA	NA
Interaction	92	−0.82 (−1.06, −0.57)		<0.001		NA	NA
5mC	92	−0.26 (−0.59, 0.07)	0.13		90	−0.29 (−1.22, 0.64)	0.54
5hmC	92	0.49 (0.17, 0.80)	0.003		90	0.36 (0.01, 0.71)	0.04
T2 Pb
Main effect	93	0.38 (0.24, 0.52)	<0.001			NA	NA
Interaction	93	−0.68 (−0.87, −0.49)		<0.001		NA	NA
5mC	93	−0.19 (−0.40, 0.02)	0.07		92	0.01 (−1.00, 1.01)	0.99
5hmC	93	0.27 (0.07, 0.47)	0.01		92	0.38 (0.03, 0.73)	0.03
T3 Pb
Main effect	93	0.46 (0.29, 0.63)	<0.001			NA	NA
Interaction	93	−0.82 (−1.04, −0.60)		<0.001		NA	NA
5mC	93	−0.35 (−0.60, −0.11)	0.004		92	−0.01 (−0.79, 0.77)	0.98
5hmC	93	0.45 (0.21, 0.68)	<0.001		92	0.32 (0.03, 0.61)	0.03
Tibia
Main effect	87	0.13 (0.01, 0.26)	0.04			NA	NA
Interaction	87	−0.28 (−0.45, −0.10)		0.002		NA	NA
5mC	87	−0.08 (−0.55, 0.39)	0.74		85	−0.17 (−0.75, 0.40)	0.56
5hmC	87	0.13 (−0.03, 0.28)	0.11		85	0.11 (−0.08, 0.30)	0.27
Patella
Main effect	101	0.20 (0.12, 0.28)	<0.001			NA	NA
Interaction	101	−0.40 (−0.51, −0.29)		<0.001		NA	NA
5mC	101	−0.27 (−0.56, 0.03)	0.08		98	−0.36 (−0.73, 0.01)	0.06
5hmC	101	0.17 (0.07, 0.27)	<0.001		98	0.19 (0.07, 0.32)	0.003
RAB5A
T1 Pb
Main effect	126	−0.01 (−0.02, 0.01)	0.44			NA	NA
Interaction	126	0.01 (−0.01, 0.028)		0.42		NA	NA
5mC	126	0.01 (−0.01, 0.02)	0.45		122	−0.01 (−0.03, 0.02)	0.52
5hmC	126	−0.01 (−0.04, 0.02)	0.61		122	0.00 (−0.04, 0.05)	0.91
T2 Pb
Main effect	129	−0.01 (−0.02, 0.01)	0.42			NA	NA
Interaction	129	0.01 (−0.01, 0.03)		0.25		NA	NA
5mC	129	0.01 (−0.01, 0.02)	0.34		124	0.00 (−0.04, 0.05)	0.82
5hmC	129	−0.01 (−0.03, 0.02)	0.63		124	0.00 (−0.04, 0.04)	0.97
T3 Pb
Main effect	130	0.00 (−0.02, 0.02)	0.90			NA	NA
Interaction	130	0.01 (−0.02, 0.03)		0.64		NA	NA
5mC	130	0.01 (−0.01, 0.03)	0.42		129	0.01 (−0.01, 0.03)	0.41
5hmC	130	−0.01 (−0.04, 0.03)	0.85		129	−0.01 (−0.05, 0.03)	0.69
Tibia
Main effect	122	0.01 (−0.00, 0.02)	0.19			NA	NA
Interaction	122	−0.01 (−0.03, −0.004)		0.16		NA	NA
5mC	122	0.00 (−0.01, 0.01)	0.83		119	0.00 (−0.01, 0.01)	0.95
5hmC	122	0.01 (−0.01, 0.02)	0.59		119	0.01 (−0.02, 0.03)	0.51
Patella
Main effect	139	0.00 (−0.01, 0.01)	0.72			NA	NA
Interaction	139	0.00 (−0.01, −0.01)		0.98		NA	NA
5mC	139	0.00 (−0.01, 0.01)	0.92		134	0.00 (−0.01, 0.01)	0.61
5hmC	139	0.00 (−0.02, 0.01)	0.70		134	0.01 (−0.01, 0.02)	0.50
TPPP
T1 Pb
Main effect	124	−0.07 (−0.22, 0.08)	0.39			NA	NA
Interaction	124	0.18 (−0.03, 0.36)		0.09		NA	NA
5mC	125	0.11 (−0.09, 0.30)	0.28		121	0.15 (−0.13, 0.42)	0.31
5hmC	124	−0.01 (−0.20, 0.26)	0.88		120	−0.03 (−0.30, 0.25)	0.85
T2 Pb
Main effect	127	−0.04 (−0.20, 0.11)	0.58			NA	NA
Interaction	127	0.11 (−0.10, 0.32)		0.31		NA	NA
5mC	128	0.04 (−0.16, 0.23)	0.72		123	0.03 (−0.27, 0.33)	0.85
5hmC	127	0.01 (−0.18, 0.19)	0.94		122	−0.03 (−0.33, 0.27)	0.84
T3 Pb
Main effect	128	−0.13 (−0.30, 0.03)	0.12			NA	NA
Interaction	128	0.28 (0.06, 0.50)		0.01		NA	NA
5mC	129	0.17 (−0.05, 0.39)	0.13		128	0.22 (−0.02, 0.45)	0.07
5hmC	128	−0.01 (−0.27, 0.13)	0.50		127	−0.16 (−0.38, 0.06)	0.15
Tibia
Main effect	120	0.04 (−0.05, 0.14)	0.38			NA	NA
Interaction	120	−0.11 (−0.24, −0.03)		0.13		NA	NA
5mC	121	−0.06 (−0.19, 0.06)	0.34		118	0.38 (0.04, 0.71)	0.03
5hmC	120	0.05 (−0.07, 0.17)	0.40		117	0.12 (−0.04, 0.36)	0.14
Patella
Main effect	137	−0.06 (−0.12, 0.02)	0.06			NA	NA
Interaction	137	0.11 (0.03, 0.20)		0.01		NA	NA
5mC	138	0.05 (−0.04, 0.13)	0.26		133	0.21 (−0.28, 0.70)	0.40
5hmC	137	−0.05 (−0.13, 0.02)	0.18		132	−0.06 (−0.7, 0.05)	0.26

Note: BLL, blood lead levels; CI, confidence interval; NA, -int, -value for intercept effect estimates; T1, first trimester; T2, second trimester; T3, third trimester. In each analysis, the outcome is percentages of 5mC or 5hmC at all quantified CpG sites treated as repeated measures, and all models adjust for sex, adolescent age, adolescent BLL, and batch. Outliers are defined as those with Pb measures greater than three standard deviations above the mean Pb concentration. NA represents the Model not run in the sensitivity analysis.

NINJ2 models were controlled for genotype at rs34038797 (G/C) and C/C samples () were dropped from the analysis.

Per change in maternal blood Pb at T1, T2, or T3, and change in maternal bone Pb measured in tibia and patella 1 month postpartum.

Estimates were derived from a mixed-effects model with an interaction term between prenatal Pb exposures and a dichotomous variable signifying whether the outcome was DNA 5mC or 5hmC. A positive value represents an increase in the slope toward greater 5mC with increasing Pb exposure, whereas a negative value represents an increase toward greater 5hmC with increasing Pb exposure.

Within NINJ2, as prenatal Pb exposure increased, there was a decrease in the slope between 5mC and 5hmC in whole blood, signifying a shift toward a greater proportion of 5hmC compared with 5mC. This shift toward 5hmC was reflected in the significant negative interaction terms between Pb exposure and DNA modification category ( for all Pb exposure biomarkers) (Table 2). For each unit increase () in maternal BLL during T1, there was a 0.49 (95% CI: 0.17, 0.80) increase in percent 5hmC within the NINJ2 locus () (Figure 2; Table 2). T2 BLLs demonstrated a positive association with 5hmC levels within NINJ2 [ (95% CI: 0.07, 0.47); ]. T3 BLLs were also associated with an increase in NINJ2 5hmC levels [ (95% CI: 0.21, 0.68); ], but a decrease in NINJ2 5mC [ (95% CI: , ); ]. A increase in patella bone Pb was associated with a significant increase in 5hmC [ (95% CI: 0.07, 0.27); ] and a nonsignificant decrease in 5mC [ (95% CI: , 0.03); ] in NINJ2 (Table 2).

Figure 2.

Association of DNA average percent 5mC and 5hmC in adolescent human whole blood within the regions of HCN2 (A, C, and E) and NINJ2 (B, D, and F) by T1 (A and B), T2 (C and D), and T3 (E and F) maternal Pb exposure measured in whole blood. The lines represent the simple (unadjusted) relationship between trimester-specific Pb measures with adolescent 5mC or 5hmC in males and females separately. 5mC and 5hmC values for each sample are depicted as the average among all CpG sites captured within the gene region. This figure includes selected results only; numeric estimates for associations between all exposures and genes are reported using mixed-effects models treating CpG sites as repeated measures adjusting for sex, adolescent age, adolescent BLL, and batch (and genotype rs34038797 (G/C) within NINJ2 models only) and are represented in Table 2 for all children and in Table S4 for males and females. Note: T1 Pb by HCN2 5mC () and 5hmC (); T2 Pb by HCN2 5mC () and 5hmC (); T3 Pb by HCN2 5mC () and 5hmC (); T1 Pb by NINJ2 5mC () and 5hmC () without rs34038797 SNP C/C (); T2 Pb by NINJ2 5mC () and 5hmC () without rs34038797 SNP C/C (); T3 Pb by NINJ2 5mC () and 5hmC () without rs34038797 SNP C/C (). BLL, blood lead levels; Pb, lead; T1, first trimester; T2, second trimester; T3, third trimester.

An interesting finding was that TPPP 5mC and 5hmC proportions shifted toward 5mC at higher patella bone Pb levels (), but there was no significant shift with tibia bone Pb measures (). There was also a shift toward 5mC at higher T3 BLLs (). However, there were no observed significant associations between trimester-specific Pb measures or bone Pb measures and TPPP 5mC or 5hmC. There was no significant interaction between Pb and type of DNA methylation and also no significant main effect of Pb in the model of RAB5A. Following these findings, RAB5A DNA 5mC or 5hmC stratified models showed no significant associations between trimester-specific Pb measures or tibia or patella bone Pb and percent 5mC or 5hmC.

The results of the sensitivity analysis excluding outliers for trimester-specific BLLs and bone Pb measures demonstrated associations with DNA 5mC and 5hmC in the same direction for each exposure-outcome relationship that were statistically significant in models including all subjects (Table 2). After outliers were removed, associations remained significant for T1, T2, and T3 BLLs and patella bone Pb with 5hmC in NINJ2; however, the estimate for T3 BLLs and 5mC in NINJ2 changed from a significant inverse association to a null association (; 95% CI: , 0.77) after one outlier for this exposure was excluded from the analysis. The association between T1 BLLs and 5mC in HCN2 remained significant, whereas associations of T2 and T3 BLLs with 5mC in HCN2 remained positive but were no longer statistically significant within the sensitivity analysis.

Sex-Specific Estimates

We performed sex-stratified models to estimate effect sizes for the associations between Pb exposure biomarkers and 5mC and/or 5hmC at HCN2, NINJ2, RAB5A, and TPPP loci (Table S4) separately in males and females. For HCN2, associations between increases in T2 and T3 BLLs and 5mC were positive for females and null for males [e.g., T2: (95% CI: 0.25, 1.38) vs. (95% CI: , 0.45), respectively], whereas associations with 5hmC were negative for females and positive for males [e.g., T3 BLLs: (95% CI: , 0.03) vs. (95% CI: , 0.36), respectively]. Associations with 5mC and 5hmC in HCN2 were not notably different by sex. For NINJ2, T2 BLLs were positively associated with 5hmC in females, whereas the association was close to the null for males [ (95% CI: 0.10, 0.86) vs. (95% CI: , 0.34), respectively]. Associations of 5mC in NINJ2 with T1 BLLs and patella Pb were negative for females and close to the null for males, but confidence intervals overlapped; otherwise, effect estimates were similar by sex. Consistent with the results for RAB5A in males and females combined, associations with 5mC and 5hmC were consistently null for both groups. For TPPP, patella Pb was positively associated with 5mC and negatively associated with 5hmC in males [ (95% CI: 0.05, 0.27) and (95% CI: , ), respectively], whereas corresponding estimates were close to the null in females [ (95% CI: , 0.08) and (95% CI: , 0.10), respectively].

The sex-stratified sensitivity analysis results excluding outliers for trimester-specific BLLs and bone Pb measures demonstrated some effect estimate changes (Table S5). In HCN2, the magnitude of the association between increases in T1 BLLs and 5mC became more positive in females by 0.42% [ (95% CI: 0.08, 1.59)], whereas males remained unchanged. All other associations with 5mC and 5hmC in HCN2 were not notably different in the sensitivity analysis. In NINJ2, T1 BLLs were more negatively associated with 5mC by 0.64%, whereas the association among females was diminished [ (95% CI: , 0.68) and (95% CI: (, 1.43), respectively]. The association between T2 BLLs and 5mC became positive in females [ (95% CI: , 1.48)], whereas the magnitude of males with 5hmC increased by 0.31% [ (95% CI: , 0.99)], but neither were statistically significant. The associations between patella Pb with 5mC and 5hmC were statistically significant after excluding outliers [ (95% CI: , ), (95% CI: 0.08, 0.44), respectively]; otherwise, effect estimates were similar when outliers were included.

5mC and 5hmC Correlations with Gene Expression

To compare quantified methylation levels with gene expression (normalized read counts), we compared DNA 5mC and 5hmC measures from the pyrosequencing method with gene expression from a subset of individuals who had available RNA-seq data (used as normalized read counts; Figure 3 and Figure S1). DNA 5mC was negatively correlated with NINJ2 expression ( ), whereas increased DNA 5hmC was positively correlated with NINJ2 expression (; ) (Figure 3). RAB5A had similar directionality but was not statistically significant; 5mC was negatively correlated with expression (; ), and 5hmC was positively correlated (; ). HCN2 (5mC: , ; 5hmC: , ) and TPPP (5mC: , ; 5hmC: , ) did not have statistically significant correlations.

Figure 3.

Correlations between gene expression (normalized read counts from RNA-seq) and 5mC% (left column) and 5hmC% (right column) from pyrosequencing in NINJ2 within a subset of individuals (). Regression line indicates the association between log-transformed NINJ2 expression by average 5mC or 5hmC among all CpG sites captured within the gene region of NINJ2 for each individual assessed. Note: Figure includes selected results only. Results for the other genes are reported in Figure S1.

Discussion

Prenatal and early-life Pb exposures have been associated with altered DNA methylation patterns, which serve as potential mechanistic links for Pb-induced health effects (Dou et al. 2019; Montrose et al. 2017; Wu et al. 2017; Goodrich et al. 2015; Faulk et al. 2014; Faulk et al. 2013; Bakulski et al. 2012). The majority of studies, however, quantify total DNA methylation and fail to distinguish between 5mC and 5hmC. These include methods based on sodium bisulfite conversion of DNA, which is unable to discriminate between 5mC and 5hmC, with important implications for the interpretation of published epigenetic studies focused on Pb exposure to date. To our knowledge, this is the first epigenetic study to investigate the independent associations of prenatal Pb exposure with DNA 5mC and 5hmC in human blood samples. Our study was strengthened by the inclusion of multiple biomarkers of early-life Pb exposure.

First, by using oxBS pyrosequencing, we provide evidence for the presence of variable 5hmC in human whole blood. Second, we estimated the effect of prenatal Pb exposure on both 5mC and 5hmC of four genes with neurological functions selected from an epigenome-wide study of prenatal Pb exposure (HCN2, NINJ2, RAB5A, and TPPP). We observed associations between T1, T2, and T3 BLLs with HCN2 5mC, in addition to T3 BLL and patella bone Pb with NINJ2 5mC, and showed that T1, T2, T3 BLLs and cumulative gestational Pb exposure measured in patella bone are associated with 5hmC levels in NINJ2. This finding supports the idea that prenatal Pb exposure stably alters gene-specific 5mC and 5hmC levels that can be detected into adolescence. Further, male- and female-specific associations of prenatal Pb exposure with both 5mC and 5hmC were not notably different by sex. Last, we identified correlations between 5mC and 5hmC with gene expression in one of the genes, NINJ2.

Our study provided estimates of 5mC and 5hmC in whole blood at four genes with varying methylation profiles. Levels of 5hmC in blood were relatively low, with wide variability between individuals, consistent with prior research that estimated the highest levels of 5hmC are in brain tissue, with lower levels estimated in tissues such as blood (Nestor et al. 2012; Li and Liu 2011; Globisch et al. 2010). A cross-sectional study of female neonates () and French-Canadian women ages 25–30 y, 70–75 y, and ( in each age group) reported that mean 5hmC levels were highest in cord blood DNA and lowest in peripheral blood DNA from women in the two oldest groups, suggesting declining levels with age (Buscarlet et al. 2016). Of note, we quantified average 5hmC levels across CpG sites as high as 10% within HCN2, NINJ2, and TPPP genes for 13 individuals. However, we also detected negative 5hmC values, which is an inherent limitation of the oxBS subtraction method when assessing gene-specific 5hmC via pyrosequencing. Although brain is the tissue of interest in terms of Pb toxicity, surrogate epigenetic measures such as in blood are typically necessary in longitudinal epidemiological studies. Research has shown that there are widespread changes in 5hmC occurring during human brain development, with sex-differences in 5hmC levels in fetal brain (Spiers et al. 2017). Whether 5hmC levels in blood correlate with crucial 5hmC patterning in brain and if so, at which genes, is currently unknown but needs to be characterized to understand the implications of exposure–5hmC associations in surrogate tissues (Bakulski et al. 2016). The Toxicant Exposures and Responses by Genomic and Epigenomic Regulators of Transcription (TaRGET II) Consortium is one research program that is conducting analyses of 5mC, 5hmC, and more in target and surrogate tissues (including brain and blood, respectively) from animals prenatally exposed to common environmental toxicants to understand the role that environment plays in disease susceptibility as a function of tissue-specific epigenomic perturbations (Wang et al. 2018). Multiple publicly available databases compared matched blood and brain tissue DNA methylation profiles that allow researchers to assess the level of similarity between human whole blood and various brain tissues at CpG sites available on the 450K (Edgar et al. 2017; Hannon et al. 2015) and EPIC arrays (Braun et al. 2019). They each concluded that correlations are highly variable among blood and brain DNA methylation profiles, with some CpG sites highly correlated between tissues and others exhibiting weak or no correlation; therefore, research must consider the genomic region of interest when planning to use whole blood as a surrogate for brain.

By using adolescent-age DNA methylation profiles, we observed evidence that Pb exposure may modify the balance between 5mC and 5hmC within total methylation. We detected a positive association between trimester-specific BLLs, but not bone, with 5mC in HCN2. Further, increasing prenatal Pb exposure was associated with increases in 5hmC but decreases in 5mC in whole blood in NINJ2. For example, a increase in T3 BLLs was associated with lower 5mC but higher 5hmC in NINJ2. However, the association between T3 BLLs and 5mC was null after an outlier was excluded. Metal exposures, such as exposures to Pb, cause increased oxidative DNA damage that inhibits the ability of methyltransferases to interact with DNA, and these can result in hypomethylation at some loci (Sanchez et al. 2017; Pizzino et al. 2014; Castellani et al. 2006; Valinluck et al. 2004). TET proteins are involved in the oxidation of 5mC to 5hmC, and free radicals are known to interact with these proteins (Chia et al. 2011; Coulter et al. 2013). Pb-induced oxidative stress inhibits () dehydrogenase in the mitochondria, resulting in the accumulation of , a cofactor for TET enzymes (Tretter and Adam-Vizi 2005). Thus, we hypothesize that Pb exposure would increase activity of TET enzymes, resulting in increases in 5hmC levels. Out of the four genes assessed, only the findings with NINJ2 reflect this hypothesis where increased 5hmC was associated with all Pb exposure measures, with or without Pb biomarker outliers included. In contrast, we identified a shift toward higher 5mC as Pb exposure increased for HCN2. More studies are needed to understand the mechanisms by which Pb and other metals lead to locus-specific changes in DNA methylation in either direction.

Although both males and females are vulnerable to the adverse effects of Pb exposure, there is significant evidence that sex can influence the severity of Pb neurotoxicity (Varma et al. 2017; Schneider et al. 2016). We observed some evidence for sex-specific associations, though based on magnitude, direction, and 95% CI overlap between male- and female-specific effect estimates, it is not entirely conclusive. Our study would benefit from a larger sample size to draw such conclusions around sex-specific effects. These data are inconsistent with previous studies in humans and mice in which more Pb-associated methylated CpG sites were identified in females compared with males (Singh et al. 2018). It is known that global DNA methylation profiles are dimorphic by sex when assessed in fetal, child, and adult human brain tissue (Numata et al. 2012) and influenced by hormones when assessed in brain tissue in rats (Schwarz et al. 2010). Thus, it is important to continue evaluating whether environmental impacts on the epigenome differ by sex.

Though we were limited to analyzing whole blood samples, results from this research provide preliminary evidence that Pb exposure may affect regulation of both 5mC and 5hmC. Challenges in current developmental toxicology studies include determining whether persistent epigenetic alterations are due to global changes across all cells, altered cell-type proportions, or changes in specific cell types (Bakulski et al. 2020; Jaffe and Irizarry 2014). At this time, the stability of 5hmC within different blood cell-types is unknown (Houseman et al. 2012). Therefore, it is possible that associations between prenatal Pb and gene-specific 5hmC alterations documented here are simply reflecting blood cell–type proportion shifts. The associations reported here between Pb and blood DNA methylation could be under or overrepresenting potential effects on the target tissue of interest for Pb toxicity: brain. However, additional studies have shown DNA methylation conservation among different tissue and cell types (Skinner 2016), suggesting that DNA methylation from whole blood can be a valid biomarker for tissues of interest at many some genes. Future studies should consider cell and tissue type-specific impacts of perinatal Pb exposure to explain mechanisms and reveal biological pathways underlying the adverse outcomes of Pb exposure. A limitation of this study revolves around the subtraction method of the 5hmC percentages (i.e., negative values). This method produced a wide variability among the 5hmC measures, for which the negative values might be capturing random noise within the technical measurements. This approach might explain why roughly a quarter of the participants had negative values. To address this issue, future research studies should employ methods that directly measure 5hmC [i.e., HMeDIP-seq (Tan et al. 2013)] to avoid negative measurements. Our study had a relatively wide range of Pb exposures, giving us the opportunity to investigate associations relevant to populations with a range of exposures. Larger studies should incorporate oxBS-based methodologies to study the sex-specific impact of Pb and other environmental exposures on 5hmC and 5mC. Although this project measured 5mC and 5hmC levels in only a few selected loci, it demonstrates the utility of oxBS pyrosequencing as a quantitative and high-throughput method to estimate 5hmC levels at specific regions that does not require antibodies.

In terms of gene regulation, we were able to detect statistically significant correlations within NINJ2, in which increases in gene expression were associated with hypomethylation and hyper-hydroxymethylation. This finding follows previous research that has shown that increased levels of 5mC are associated with decreased transcription factor binding at promoter/enhancer sites and suppression of transcription (Medvedeva et al. 2014), whereas 5hmC is associated with increases in gene expression (Stroud et al. 2011; Wu et al. 2011). We were unable to detect these same correlations within HCN2, RAB5A, and TPPP. Despite nearly half the samples having RNA available for RNA-seq, the use of these samples allowed for direct measurement of expression of each gene analyzed in this study. To further explore the idea that altered 5hmC in human blood affects gene regulation, future studies should include a larger sample size and should investigate expression in individual cell types using technologies such as single-cell RNA-seq.

In conclusion, these data suggest that a) prenatal Pb exposure differentially influences 5mC and 5hmC, and b) 5hmC is present and detectable in whole blood samples at variable levels. We provided little evidence for differences in these associations by sex, though power was limited in the sex-stratified analyses. These modification-specific associations would not have been detected using methods based on standard sodium bisulfite treatment, demonstrating the utility of incorporating oxBS treatment in cohort studies. Prenatal Pb exposure results in a shift from 5mC to 5hmC in NINJ2, with the opposite shift seen in HCN2 (5hmC to 5mC), and these relationships were detected using all blood Pb biomarker measures. Prenatal Pb exposure has been previously associated with a decrease in IQ and adverse effects on neurobehavioral outcomes, but whether epigenetic mechanisms contribute to these long-term effects is not well characterized. Given the potential importance of Pb-induced DNA methylation perturbations on processes related to early neurodevelopment and potential late-onset delays (i.e., adulthood) in cognition (Bakulski et al. 2012), careful examination of Pb-induced 5mC and 5hmC alterations may eventually improve our knowledge of the epigenetic pathways involved in neurodevelopment. Our results suggest the prenatal Pb exposure may modify both 5mC and 5hmC, and the relationship between these changes and neurodevelopment outcomes merits further study. This study and others are providing additional evidence around DNA 5hmC as being an environmentally responsive modification (Kochmanski et al. 2018a). Altogether, it is vital for environmental studies to consider 5hmC as a regulatory response mechanism to environmental exposures.

Click here for additional data file.

Click here for additional data file.