- Open Access
Interaction between prenatal pesticide exposure and a common polymorphism in the PON1 gene on DNA methylation in genes associated with cardio-metabolic disease risk—an exploratory study
Clinical Epigeneticsvolume 9, Article number: 35 (2017)
Prenatal environmental conditions may influence disease risk in later life. We previously found a gene-environment interaction between the paraoxonase 1 (PON1) Q192R genotype and prenatal pesticide exposure leading to an adverse cardio-metabolic risk profile at school age. However, the molecular mechanisms involved have not yet been resolved. It was hypothesized that epigenetics might be involved. The aim of the present study was therefore to investigate whether DNA methylation patterns in blood cells were related to prenatal pesticide exposure level, PON1 Q192R genotype, and associated metabolic effects observed in the children.
Whole blood DNA methylation patterns in 48 children (6–11 years of age), whose mothers were occupationally unexposed or exposed to pesticides early in pregnancy, were determined by Illumina 450 K methylation arrays.
A specific methylation profile was observed in prenatally pesticide exposed children carrying the PON1 192R-allele. Differentially methylated genes were enriched in several neuroendocrine signaling pathways including dopamine-DARPP32 feedback (appetite, reward pathways), corticotrophin releasing hormone signaling, nNOS, neuregulin signaling, mTOR signaling, and type II diabetes mellitus signaling. Furthermore, we were able to identify possible candidate genes which mediated the associations between pesticide exposure and increased leptin level, body fat percentage, and difference in BMI Z score between birth and school age.
DNA methylation may be an underlying mechanism explaining an adverse cardio-metabolic health profile in children carrying the PON1 192R-allele and prenatally exposed to pesticides.
A considerable part of modern pesticides has neurotoxic and/or endocrine disrupting properties [1–3] and therefore the potential to disturb development of neurobehavioral, neuroendocrine, and reproductive functions [4–8] especially if exposure occurs during vulnerable time periods in fetal life or early childhood. To investigate potential health effects of prenatal pesticide exposure, we have followed a cohort of children, whose mothers were employed in greenhouse horticulture in pregnancy. Some of the mothers were occupationally exposed to mixtures of pesticides in the first trimester before the pregnancy was recognized, and preventive measures were taken. Findings from this cohort include associations between maternal pesticide exposure and lower birth weight followed by increased body fat accumulation during childhood , impaired reproductive development in boys [10, 11], and earlier breast development  and impaired neurobehavioral function in girls .
The HDL-associated enzyme paraoxonase 1 (PON1) catalyzes the hydrolysis of a wide range of substrates including some organophosphate insecticides [14, 15]. It also protects lipoproteins from oxidative modifications and hence against development of atherosclerosis [16, 17]. A common polymorphism in the coding sequence of the PON1 gene substitutes glutamine (Q) to arginine (R) at position 192. This substitution seems to affect both properties of the enzyme, and several studies have indicated an increased risk of cardiovascular disease in R-allele carriers [17, 18]. To investigate if this polymorphism affected the sensitivity to prenatal pesticide exposure, the PON1 Q192R genotype was determined in the children. We found a marked interaction between prenatal pesticide exposure and the PON1 Q192R genotype. At school age, exposed children with the R-allele had significantly higher BMI, body fat percentage, abdominal circumference, and blood pressure compared to unexposed children with the same genotype. In the group of children with the QQ genotype, there was no effect of prenatal pesticide exposure on these parameters . In addition, serum concentrations of leptin, glucagon, and plasminogen activator inhibitor type-1 (PAI-1) were enhanced in prenatally pesticide exposed children with the R-allele, also after adjusting for BMI  which also indicates disturbance of metabolic pathways related to development of metabolic syndrome [21–23]. In addition, leptin seemed to be a mediator of the increased fat accumulation during childhood related to prenatal pesticide exposure in children with the PON1 192R-allele . Thus, the obtained results indicate a gene-environment interaction between pesticide exposure and PON1 gene heterogeneities already in early prenatal life that might enhance the risk of cardio-metabolic diseases later in life.
The mechanism behind this interaction is not yet understood but might be mediated by epigenetic alterations depending on both genotype and prenatal exposure. Epigenetic marks, including DNA methylation and covalent histone modifications, are dynamic and can adapt to a variety of external stimuli . Furthermore, during fetal development extensive de- and re-methylation events are taking place making this period highly vulnerable for epigenetic changes caused by environmental conditions . Indeed, emerging evidence in experimental animals and in humans associate altered DNA methylation patterns with a variety of prenatal exposures including dietary factors, parental care, infections, smoking, and environmental pollutants [26–31]. In experimental animals, early life changes in DNA methylation have been associated with diet-induced obesity and insulin resistance . Recently, also human studies have suggested that DNA methylation patterns at birth are related to birth weight and fat mass later in childhood [33, 34]. The aim of this exploratory study was to investigate whether methylation patterns in blood samples of school children were related to prenatal pesticide exposure, PON1 Q192R genotype, and adverse health outcomes already observed in the children. We hypothesized that the health effects associated with early prenatal pesticide exposure were related to differential epigenetic modifications in children with the QQ-genotype and children carrying the R-allele.
This study is a part of an ongoing prospective study including 203 children born between 1996 and 2001 by female greenhouse workers. The children were examined for the first time at 3 months of age  and followed-up at school age when 44 new age-matched controls were included , and the PON1 genotype was determined for 141 children . For this exploratory study, 48 pre-pubertal (Tanner Stage 1) children, whose mothers reported not to have smoked during pregnancy, were selected equally distributed between the PON1 192QQ and QR/RR genotype. The QR/RR genotype group consisted of 3 children with the RR genotype and 21 with the QR genotype. After excluding children of mothers who smoked in pregnancy, the number of unexposed controls within each genotype was low, 20 with the QQ genotype and 16 with the QR/RR genotype. DNA qualified for methylation analysis was only available for 11 and 12 of these children, respectively. For each genotype, we then used individual matching to select one exposed child of same sex and age for each of the controls. For the QQ-genotype, two exposed children were selected for each of two controls to obtain 24 children. Thus, in total we used data from 13 exposed and 11 unexposed children with the QQ genotype, and 12 exposed and 12 unexposed children with the QR/RR genotype (Table 1).
Recruitment, characteristics, exposure categorization, and clinical examinations of the children have previously been described in detail [9, 11, 19]. Briefly, we recruited pregnant women working in greenhouses and referred to the local Department of Occupational Health for risk assessment of their working conditions and guidance for safe work practices during pregnancy. Detailed information about working conditions inclusive pesticide use for the previous 3 months was obtained from maternal interview at enrollment (gestational weeks 4–10) and supplemented by telephone contact to the employers. For all women, re-entry activities (such as moving or packing potted plants or nipping cuttings) constituted their main work functions. Approximately 20% of the women reported having been directly involved in applying pesticides, mainly by irrigating fungicides or growth retardants. Only few (6%) of the women had applied insecticides. The women were categorized as occupationally exposed if pesticides were applied in the working area more than once a month, and the women handled treated plants within 1 week after treatment and/or the women were directly involved in applying pesticides. The women were categorized as occupationally unexposed if none of the above criteria was fulfilled. All exposure assessments and categorization of the mothers as pesticide exposed or unexposed were performed independently by two toxicologists before the first examination of the children. Women categorized as pesticide exposed went on paid leave or were moved to work functions with less or no pesticide exposure shortly after enrollment. Hence, the exposure classification relates to the early weeks of the first trimester before study enrollment.
The exposure situation was complex since the use of specific pesticides varied with time and location, both within the same company and between companies, depending on the plant production and the type of pest to be controlled. Out of 124 different active pesticide ingredients used in the greenhouses were 59 insecticides (17 organophosphates, 12 pyrethroids, 9 carbamates, and 21 others), 40 fungicides, 11 growth regulators, and 14 herbicides. Some were used only in few greenhouses or in short periods, whereas others were used more often. Organophosphate insecticides were used to some extent in the working areas for 91% of the exposed mothers in the entire cohort, and for 24 out of the 25 exposed mothers whose children were included in this study. The most used organophosphates were dichlorvos, dimethoate, and chlorpyrifos. Other frequently used pesticides were the pyrethroid insecticides deltamethrin and fenpropathrin; the carbamate insecticides methiocarb, pirimicarb, and methomyl, and the fungicides fenarimol, prochloraz, tolclofos-methyl, vinclozolin, iprodion, and chlorothalonil. In general, the time interval between applying insecticides and working in the treated areas was longer (1–3 days) than for fungicides and growth regulators (often a few hours). Because of the complexity of the exposure situation and because most of the women at enrollment had been off work for some days while the risk assessment of their working conditions was performed, biomonitoring of the exposure was not feasible. A complete list of the pesticides used in the greenhouses can be obtained from the corresponding author.
At follow-up at age 6 to 11 years, 177 children underwent a standardized clinical examination in which systolic and diastolic blood pressure, pubertal staging, height, weight, thickness of skin folds, and other anthropometric parameters were measured . The same pediatrician performed all clinical examinations blinded to information about maternal pesticide exposure during pregnancy.
Venous non-fasting blood samples were collected (between midmorning and late afternoon) in EDTA-coated and uncoated vials (Venoject). After centrifugation at 2000 g for 10 min at 20 °C, buffy coat for genotyping and epigenetic analysis was separated from the EDTA-treated samples. Buffy coat and serum from the uncoated vials were stored at −80 °C until analysis.
As previously described , C-108T (rs705379) and Q192R (rs662) polymorphisms of the PON1 gene was determined by the Taqman-based allele discrimination using the ABI Prism 7700 Sequence Detection System, serum activity of PON1 was determined by spectrophotometry with paraoxon as substrate, and insulin (proinsulin and insulin) and leptin concentrations in serum were determined by commercial ELISA hormone kits from RayBio.
Genotyping and all serum analyses were performed blinded to both exposure information and examination outcomes.
DNA from buffy coat samples was extracted using QIAamp DNA Blood Mini Kit (Qiagen, Hilden, Germany). The blood spin protocol was applied according to manufacturer’s instructions. Samples were eluted in 100 μl elution buffer. DNA samples were bisulfite-converted using the EZ DNA methylation kit from Zymo according to manufacturer’s instructions. Successful bisulfite conversion was checked using a bisulfite-specific PCR of an amplicon in the SALL3 gene (see Additional file 1 for primer sequences). Only samples showing an intense band on agarose gel were further analyzed by the 450 K methylation array. As a negative control non-converted gDNA was used.
DNA methylation and data preprocessing
The Infinium HumanMethylation450 BeadChip array (Illumina, San Diego, CA, USA) was used to measure DNA methylation genome-wide. 4 μL of bisulfite-converted DNA from each sample was amplified, fragmented, precipitated, resuspended, and subsequently hybridized onto the BeadChips. After overnight incubation of the BeadChips, unhybridized fragments were washed away, while hybridized fragments were extended using fluorescent nucleotide bases. Finally, the BeadChips were scanned using the Illumina iScan system to obtain raw methylation intensities for each probe.
We used the R package RnBeads to preprocess the Illumina 450 K methylation data . Cg-probes were filtered before normalization based on following criteria: probes containing a SNP within 3 bp of the analyzed CpG site, bad quality probes based on an iterative Greedycut algorithm where a detection p value of 0.01 was set as a threshold for an unreliable measurement, and probes with missing values in at least one sample. After filtering these cg-probes, beta values (ratio of methylated probe intensity versus total probe intensity) were within-array normalized using the beta mixture quantile dilation (BMIQ) method . Another filtering step was performed after normalization based on the following criteria: probes measuring methylation not at CpG sites and probes on sex chromosomes. The two filtering steps removed a total of 20,338 cg-probes and ended up with a data set of normalized methylation values for 465,239 cg-probes. Beta values were transformed to M values (M = log2(β/(1−β))) prior to further analyses. Principal component analysis (PCA) was conducted to detect possible batch effects. Associations between the first eight principal components and possible batch effect covariates were measured. The Kruskall-Wallis test was used to find associations with sentrix_ID (BeadChip), while the two-sided Wilcoxon sum rank test was used for associations with the processing date, exposure and PON1 Q192R genotype. Significant associations between principal component 2 and sentrix_ID (BeadChip) and processing date were suggestive for batch effects and were therefore corrected using the ComBat function in the SVA R package  (Additional files 2 and 3). Raw and normalized array data were uploaded to the Gene Expression Omnibus (GEO) database and have accession number: GSE90177.
For each sample, the relative cell type contribution was measured using the approach described by Houseman et al. . Reference methylomes of each blood cell type (granulocyte, CD4+ T-cell, CD8+ T-cell, B-cell, monocyte, NK-cell) were obtained from the study of Reinius et al. using the FlowSorted.Blood.450 K R package . The analysis was limited to the 100,000 most variable sites. The top 500 cg-probes associated with the cell types were used to estimate the relative cell type composition in each sample. One-way ANOVA was used to determine differences in relative cell type composition between the exposed and the unexposed children and between children with the QQ and QR/RR genotype. Associations between relative cell type composition and health outcomes (percentage body fat, delta BMI z-scores from birth to school age, and BMI Z scores), leptin levels and age were analyzed using simple linear regression.
Differential methylation was analyzed both at the single CpG site level and at the region level (Fig. 1). At the single CpG site level, multiple linear regression (Matlab version 2014b, The Mathworks®, Natick, MA, USA) was performed in which methylation was the dependent variable and PON1 Q192R genotype and prenatal pesticide exposure (yes/no) were the independent variables. Our statistical approach was designed to explain—at the level of methylation—the previously reported gene-environment interaction between the paraoxonase 1 (PON1) Q192R genotype and prenatal pesticide exposure leading to an adverse cardio-metabolic risk profile at school age among children carrying the R-allele . Thus, our primary interest was to identify methylation marks associated with exposure that were more altered in R-allele carriers than in QQ-homozygotes. Two statistical models were included in our statistical approach. In the first model, effect modification (interaction) of exposure by PON1 Q192R genotype was allowed by including main effects (exposure and genotype) and cross-product terms (exposure*genotype) in the models. Statistical significant effects of exposure in the PON1 192QR/RR group were defined as follows: P value interaction term ≤0.1 and P value of exposure in the QR/RR group ≤ 0.001. This model allows studying synergistic effects where the combined effect of prenatal exposure and in the QR/RR group is greater than the sum of the effects of each factor alone. In the second model, effect modification of exposure by PON1 Q192R genotype was not assumed (no cross product term included). Statistical significant effects of exposure were defined as follows: P value of exposure ≤0.001, P value of PON1 genotype ≤ 0.1. In this model the combined effect of exposure and being R-allele carrier is equal to the sum of the effect of each factor separately. For both models, the associations were adjusted for child sex. To identify probes that were most aberrant in the exposed QR/RR group, we set an additional filter for both models in which we defined that the prenatally exposed QR/RR group should either be highest or lowest methylated (based on mean methylation level) as compared to the other three groups (exposed QQ, unexposed QR/RR and unexposed QQ). These sites are defined as significantly differentially methylated positions (sig-DMPs) in the remainder of this text. Sig-DMPs were annotated using the HumanMethylation450 v1.2 manifest file. The freely available EpiExplorer tool was used to add further annotation including chromatin state segmentation and histone modifications based on the UCSC hg19 browser . Genomic locations of transcription factor binding sites (TFBS) were directly downloaded from the UCSC h19 genome browser. Enrichment or depletion of sig-DMPs in a particular genomic region was determined using the Fisher’s exact test.
Differentially methylated regions (DMRs) were detected using the limma-based DMRcate R package . We only looked for regions differentially methylated between the exposed QR/RR group and one of the other groups (exposed QQ, unexposed QR/RR and unexposed QQ). In line with identification of sig-DMPs, significant regions (P adj value < 0.05) were selected in which the exposed R-allele carriers showed either the highest or the lowest methylation state which are called sig-DMRs in the remainder of this text. P values were corrected for multiple testing using the Benjamini-Hochberg method (P adj).
We used bisulfite pyrosequencing to further verify the methylation differences observed in the methylation array. We selected regions in four genes that are known to be involved in metabolism: LEP, GPR39, PPARG, and OPCML (Additional file 4). LEP DNA methylation has been associated with BMI, birth weight, and cholesterol levels [42–44]. Also, maternal conditions have an effect on the methylation status of the LEP promoter [45–48]. GPR39 belongs to the ghrelin receptor family and was shown to be associated with obesity . PPARG is a nuclear receptor involved in regulation of lipid and metabolism as well as a target for some obesogenic endocrine disruptors [20, 50–53]. Furthermore, PPARγ is directly involved in the regulation of PON1 gene expression [54–56]. OPCML (Opioid Binding Protein/Cell Adhesion Molecule-like) is a member of the IgLON family. A SNP in the OPCML gene was associated with coronary artery calcified plaque in African Americans with type 2 diabetes . A mouse and human GWAS analysis identified an OPCML SNP associated with obesity traits and visceral adipose/subcutaneous adipose ratio, respectively [58, 59]. 1 μg DNA from each sample was bisulfite-converted using the EpiTect Fast bisulfite Conversion Kit (Qiagen, Hilden, Germany) according to manufacturer’s instructions. 15 ng of bisulfite-treated DNA was subsequently used in PCR amplification using the PyroMark PCR Kit (Qiagen, Hilden, Germany). Reverse primers were biotinylated to get biotin-labeled PCR products. Finally, DNA sequences were pyrosequenced using the PyroMark Q24 Advanced instrument (Qiagen, Hilden, Germany). First, streptavidin-coated Sepharose beads (High Performance, GE Healthcare, Uppsala, Sweden) were used to immobilize the biotin-labeled PCR products. Subsequently, PCR products were captured by the PyroMark vacuum Q24 workstation, washed and denaturated. The single stranded PCR products were mixed and were annealed with their corresponding sequencing primer. After the pyrosequencing run was finished, the results were analyzed using the PyroMark Q24 Advanced software (Qiagen, Hilden, Germany). Biotinylated-reverse, forward, and sequencing primers were designed using the PyroMark Assay Design 2.0 software (Qiagen, Hilden, Germany) (Additional file 1).
For a subset of sig-DMPs and sig-DMRs we analyzed (1) whether methylation is a mediator between exposure in PON1 192R-allele carriers and leptin levels; and (2) whether methylation is a mediator between exposure in PON1 192R-allele carriers and body fat accumulation (using delta BMI-score (from birth to school age), and percentage body fat as endpoints). Mediation analysis was restricted to the subset of the methylation data that overlap between the list of sig-DMPs (interaction model) and sig-DMRs. The analysis was performed by the procedure described by Baron and Kenny (1986) . Leptin concentrations were logarithmically (ln) transformed prior to analysis. In mediation analysis considering body fat percentage and leptin, the models were adjusted for sex. As sex was already considered when calculating BMI Z score, associations considering mediation between pesticide exposure and BMI Z score were not adjusted for sex.
To demonstrate mediation, four requirements must be met: (model 1) the dependent outcome variable (leptin or a body fat measure) should be significantly associated with pesticide exposure (independent variable); (model 2) the DNA methylation mark (mediator) should be significantly associated with pesticide exposure; (model 3) the dependent variable should be significantly associated with the DNA methylation mark; and (model 4) the DNA methylation mark should be a significant predictor of the outcome variable, while controlling for pesticide exposure. The estimated exposure-related change in the outcome variables in model 4 should be less than in model 1 to demonstrate partial mediation, and drop to zero to demonstrate full mediation. A P value below 0.05 was used as a cut-off for statistical significance in each of the models.
Functionally relevant mediators, i.e., mediators that have been reported to be involved in development of weight gain/obesity, insulin resistance/diabetes, cardiovascular disease, and/or fetal growth retardation were subjected to further statistical analysis. R-package “mediation” was used to calculate the significance of the causal mediation effect using a bootstrapping approach . It should be noted that the age of the children varied between 6 and 11 years at the follow-up examination where blood was collected. As child age might affect methylation levels, the exposed and unexposed children selected for this study were age-matched within each genotype.
Ingenuity Pathway Analysis (IPA, Ingenuity Systems®) was used for biological interpretation. The overlap between sig-DMPs and sig-DMRs was determined and used as input for canonical pathway analysis. A Fisher’s exact test was used to determine whether the gene lists include more genes associated with a given pathway as compared to random chance (P value ≤ 0.05).
The DisGeNet platform (http://www.disgenet.org/) was used to screen for gene disease associations . The database (currently) contains 429111 gene disease associations for which the platform provides a reliability score (DisGeNET Score). This score ranges from 0 to 1 and takes into account the number and type of sources (level of curation, organisms), and the number of publications supporting the association (for further details we refer to the DisGeNet website). For this manuscript, we extracted the associations with a score above 0.1. By this criterion, 34180 gene disease associations remain in the database. Associated diseases were mapped to the overlapping list of genes between sig-DMPs and sig-DMRs.
Descriptive statistics of the study population
Characteristics, inclusive anthropometric data, for the 48 children (6–11 years of age) are presented in Table 1. In accordance with the findings for the whole cohort , birth weights were significantly lower and measures of body composition (abdominal circumference, skin fold thickness), increase in BMI Z score from birth to school age (delta BMI Z score), diastolic blood pressure, and leptin and insulin concentrations at school age were significantly higher in the exposed PON1 192QR/RR group compared with the unexposed QR/RR group. For children with the QQ genotype, none of the variables was significantly affected by prenatal pesticide exposure (P > 0.05).
Prenatal pesticide exposure-induced methylation changes at CpG sites enriched in promoter regions in PON1 192R-allele carriers
Genome-wide DNA methylation in whole blood samples from the children was determined by Illumina 450 K methylation arrays and differential methylation patterns related to prenatal pesticide exposure and PON1 Q192R genotype were analyzed. First differential methylation was detected at the single CpG level using two multiple linear regression models (Fig. 1). Because relative cell type composition was not associated with pesticide exposure and PON1 Q192R genotype (Additional file 5), differences in cellular composition were not further considered in the workflow of statistical analysis. Allowing effect modification by PON1 Q192R genotype, 767 sig-DMPs were identified of which 128 were hypermethylated and 639 hypomethylated in prenatally exposed PON1 192R allele carriers. When effect modification was not assumed, and the interaction term between exposure and PON1 genotype was removed from the models, 70 sig-DMPs of which 44 were hypermethylated and 26 hypomethylated in prenatally exposed PON1 192R-allele carriers were identified. Hierarchical clustering of the samples using all the sig-DMPs demonstrated a clear cluster of exposed PON1 192R-allele carriers (Fig. 2). Confidence in detection of differentially methylated genes was increased by further analysis showing that the changes in methylation were not restricted to single CpGs, but were often located in regions or so called differentially methylated regions (DMRs). 5002 sig-DMRs were identified, of which 2264 were hypermethylated and 2738 hypomethylated in the exposed PON1 192R carrier group compared to the other groups. Allowing interaction between exposure and PON1 Q192R genotype to determine sig-DMPs, 547 out of 767 sites (71.3%) were overlapping with the list of sig-DMRs (Additional file 6). When effect modification was not considered, 57 out of 70 sites (81.4%) were overlapping (Additional file 7).
The pyrosequencing methylation percentages confirmed the robustness of Illumina results. They showed significant positive correlations with the Illumina 450 K beta values for all measured CpG probes (Fig. 3), except for two probes in the LEP gene (cg00840332 and cg26814075) which were borderline significant (P value: 0.07 and 0.16, respectively). The reason for this less strong correlation between the Illumina and the pyrosequencing LEP methylation is probably the lower inter-individual methylation variability in this region compared to GPR39 and PPARG.
In accordance with the Illumina results, the pyrosequencing LEP methylation values were not associated with pesticide exposure and/or PON1 Q192R genotype. Furthermore, the serum leptin concentrations were not correlated with LEP methylation status (data not shown). For GPR39, the region analyzed with pyrosequencing contained three Illumina cg-probes (cg17172683, cg11552903, and cg18444763), which showed a high correlation (r > 0.78) between the Illumina beta values and the pyrosequencing methylation percentages. For most CpGs in the pyrosequencing region, we could verify a significant exposure effect, and in each CpG site, prenatally exposed children with the QR/RR genotype had the lowest mean methylation value (Additional files 8 and 9). In the PPARG promoter, a region was selected containing one Illumina cg-probe (cg01412654). Also here, the correlation between the 450 K Illumina beta values and the pyrosequencing methylation percentages was strong. However, DNA methylation in this region was not associated with pesticide exposure and/or PON1 Q192R genotype and did not correlate with PON1 activity (data not shown). A region in the OPCML gene was found to be higher methylated in prenatal pesticide-exposed children carrying the PON1 192R-allele. The significant interaction effect between pesticide exposure and PON1 Q192R genotype could be successfully verified by pyrosequencing. The pyrosequencing methylation values were significantly higher methylated in exposed children compared to unexposed children carrying the PON1 192R-allele for most of the CpG sites in the region (Additional file 9).
Next, we questioned whether the sig-DMPs were enriched or depleted in a specific genomic location (Fig. 4). Sig-DMPs for which interaction between exposure and PON1 Q192R genotype was seen, were enriched in promoter regions (200 and 1500 bp upstream of transcription start sites) and depleted in gene bodies, 3′UTRs and intergenic regions. This was also evident when we overlapped the sig-DMPs with different chromatin states, where we observed enrichment in active and poised promoters, while DMPs were depleted in regions like transcriptional elongation, weak transcribed, and heterochromatin regions. Furthermore, DMPs were significantly more located in CpG islands and less observed in CpG poor regions. Sig-DMPs found in the models without an interaction term were not enriched or depleted in a particular genomic region.
We also looked for enrichment in TFBS using available ChIP ENCODE data from the UCSC genome browser. Thirty-nine of the 161 TFBS were significantly enriched for the model with interaction (Bonferroni adjusted P value < 0.05) while no enrichment was found for the sig-DMPs found in the model without interaction (Additional file 10).
DNA methylation differences were enriched for genes involved in neuro-endocrine signaling pathways
Overlapping the list of sig-DMPs with the list of sig-DMRs we obtained a robust and a high confidence list of differentially methylated genes (N = 446). This list was used as an input for ingenuity pathway analysis. The top enriched canonical pathways (based on P value) were dopamine-DARPP32 feedback cAMP signaling, corticotrophin releasing hormone signaling, nNOS signaling in neurons, CDK5 signaling, and neuregulin signaling (Table 2). In the context of this manuscript, other significantly enriched pathways such as mTOR signaling (rank 9, −log(P value) = 1.85) and type II diabetes mellitus signaling (rank 16, −log(P value) = 1.51) are also highly relevant.
DNA methylation (partially) mediates associations between pesticide exposure and higher leptin concentrations, body fat content, and delta BMI Z scores
The list of genes that overlaps between sig-DMPs (as identified by the interaction model) and sig-DMRs was also used as input for mediation analysis. We identified, respectively, 20, 31, and 45 candidate methylation marks that (partly) mediate the effect between pesticide exposure and serum leptin concentrations; delta BMI Z score; and body fat content. Based on applied cut-off criteria, we were not able to identify methylation marks that mediate the effect on BMI Z score. Currently known gene disease associations allowed to extract mediators that were reported to be involved in development of weight gain/obesity, insulin resistance/diabetes, cardiovascular disease, and/or fetal growth retardation. This subset of mediators is given in Table 3. Based on Baron and Kenny’s steps to analyze mediation, the association between pesticide exposure and delta BMI Z score was partially mediated by hypomethylation of UQCRC2, MTNR1B and GRIN2A, and by hypermethylation of FABP4 and LRP8. Methylation of UQCRC2 and LRP8 was also a partial mediator in the association between pesticide exposure and body fat percentage. LRP8 was also found to mediate the association between pesticide exposure and serum leptin concentration. The P value for significance of the causal mediation effect is included in Table 3 and was below 0.1 for all mediators except for UQCRC2 and GRIN2A. Irrespective of disease association of interest, the full list of potential mediators is provided in Additional file 11 which also includes the outcome of the statistical analysis.
DNA methylation at the PON1 promoter is affected by the PON1-108CT SNP (rs705379) and negatively correlated with paraoxonase 1 activity
Beside the genome-wide DNA methylation effects of the PON1 Q192R genotype, we also observed a wide variation in DNA methylation in the PON1 promoter itself for nine Illumina cg-probes. Prenatal pesticide exposure and/or PON1 Q192R genotype did not affect PON1 promoter methylation status. However, another polymorphism (rs705379, PON1 -108CT) in the promoter region of PON1 could explain a large extent of this variation (Fig. 5). Individuals homozygous for the T-allele showed higher methylation values compared with the homozygous C-allele carriers. As expected, heterozygous individuals had an intermediate methylation value. Furthermore, the paraoxonase 1 activity was significantly associated with DNA methylation in the PON1 promoter region, with higher methylation values resulting in lower paraoxonase 1 activity (Fig. 6). PON1 Q192R genotype had the strongest effect on PON1 activity, while variation in PON1 promoter methylation led to a smaller but significant effect on PON1 activity.
We found that prenatal pesticide exposure was associated with a differential DNA methylation profile in children carrying the PON1 192R-allele compared to children with the PON1 192QQ genotype and unexposed children. 767 sig-DMPs were identified of which 128 were hypermethylated and 639 hypomethylated in prenatally exposed PON1 192R-allele carriers. The profiles of PON1 192R-allele carriers are clustered together. As far as we know, our study is the first one to demonstrate a link between epigenetics and genetic susceptibility towards pesticide exposure in fetal life. Our study supports a linkage of a differential methylation pattern and higher body fat content and serum leptin concentrations in school age children dependent on both PON1 Q192R genotype and prenatal pesticide exposure.
The majority of the detected sig-DMPs were hypomethylated in exposed children with the PON1 192QR/RR genotype. Interestingly, these DMPs were mainly located in gene promoters, CpG islands and transcription factor-binding sites, suggesting a possible direct link with gene expression. To increase the confidence of our findings, we also screened for differentially methylated regions. Most of the single CpG sites were part of a DMR suggesting that these were independent of technical variation and could be considered as reliable.
Technical reliability of the outcomes from the 450 K Illumina methylation array was successfully confirmed by bisulfite pyrosequencing of corresponding CpG probe regions of four selected genes, i.e., LEP, PPARG, GPR39, and OPCML for which corresponding probes were available.
LEP was chosen because we previously found leptin to be a potential mediator of the association between prenatal pesticide exposure and body fat accumulation in children with the PON1 192R-allele . In addition, multiple studies demonstrated associations between LEP DNA methylation and BMI, birth weight, and cholesterol concentrations [42–44]. LEP was also found to be differentially methylated in the offspring of mothers suffering from the Dutch winter famine . However, our pyrosequencing results did not demonstrate a correlation between leptin DNA methylation and leptin serum concentrations, and prenatal pesticide exposure was not associated with changes in leptin DNA methylation. This suggests that the higher leptin concentration observed in exposed children with the R-allele is not due to a direct effect on DNA methylation of the leptin gene itself. Another gene whose methylation was confirmed by pyrosequencing was PPARG, a nuclear receptor controlling the expression of genes involved in lipid storage and glucose metabolism and target for obesogenic compounds [50–53]. Furthermore, PPARγ is involved in the regulation of PON1 expression [54–56]. However, we did not find a correlation between PPARG DNA methylation and PON1 activity (data not shown). In our dataset, prenatal pesticide exposure did not seem to change PPARG methylation levels irrespective of PON1 Q192R genotype.
Reduced GPR39 DNA methylation observed in prenatally pesticide exposed R-allele carriers was confirmed with pyrosequencing. GPR39 is receptor for obestatin (belonging to the ghrelin receptor family), involved in regulation of appetite and glucose homeostasis [63, 64] and associated with obesity . Furthermore, GPR39 knock-out mice showed an increased fat accumulation due to changes in lipolysis and energy expenditure . So, mis-regulation of this gene due to methylation changes might lead to an obese phenotype. To our knowledge, no other study has yet reported methylation differences in this region associated with obesity or metabolic disorders or showed links with pesticide exposure.
The higher methylation values of the OPCML DMR in exposed children carrying the PON1 192R-allele could be confirmed by pyrosequencing. OPCML encodes for a protein belonging to the IgLON family. OPCML was shown to be a tumor suppressor and inactivated by DNA methylation in a variety of cancer types [65–68]. There is also a link with metabolic diseases, as SNPs in this gene were found to be associated with obesity traits, coronary artery calcified plaque, and visceral adipose/subcutaneous adipose ratio [57–59].
Further analysis revealed that the differences in DNA methylation were most pronounced in genes involved in neuro-endocrine signaling pathways, including “dopamine-DARPP32 feedback in cAMP signaling”, “corticotropin releasing hormone signaling”, “nNOS signaling in neurons”, and “CDK5 signaling”. These pathways are important in the control of food intake and energy balance. Dopamine signaling, for example, is one of the key players in the reward pathway, also controlling food intake and preferences. Reduced dopamine signaling is assumed to induce overeating [69, 70]. In mice, a high-fat diet during pregnancy resulted in altered gene expression and DNA methylation of the dopamine transporter gene in the offspring, leading to an increased preference for sucrose and fat . Another study found similar results in prenatally stressed rats given a high fat-sucrose diet . These studies suggest that prenatal and early life conditions may influence food intake and food preferences later in life through modulation of the dopamine pathway [73–77]. Organophosphate insecticides have been shown to modulate dopamine signaling . Furthermore, low-dose exposure of neonatal rats caused metabolic dysfunction resembling prediabetes, and in adulthood, exposed animals gained excess weight when fed a high fat diet compared to unexposed rats on the same diet .
Corticotropin-releasing hormone (CRH) is a neuropeptide secreted in response to stress. However, a role for CRH in regulating energy balance and food intake has also been described [80–82] including a relation to the action of leptin .
Also NOS1 neurons are involved in energy balance and food intake [84–86]. Knock-out of NOS1 in leptin receptor- and NOS1-expressing hypothalamic neurons results in hyperphagic obesity, decreased energy expenditure, and hyperglycemia in mice . Interestingly, organophosphates have been shown to alter NOS1-expressing neurons during development in mice [87, 88].
Neureguline 1 treatment in rodents has been shown to increase serum leptin concentrations, prevent weight gain, and lower food intake. Hence, affecting this pathway may also change food intake and energy metabolism [89, 90].
A limitation of this study is that the methylation profile is measured at the same time as health outcomes and causality as such cannot be proven. Some of the genes that relate to the sig-DMPs are involved in neuro-endocrine pathways that regulate appetite and energy balance, but this study cannot rule out if these sig-DMPs are a consequence of alterations of food habits and physical activity among the exposed children with the PON1 192R-allele or an underlying mechanism. However, the mediation analysis suggested that some of the differentially methylated marks are on the mechanistic pathway between prenatal pesticide exposure and the measured outcomes. This result suggests that, at least in some, CpG sites a change in methylation might contribute to metabolic disturbances later in life. Furthermore, the association was not significant between pesticide exposure and BMI Z score as such, but between pesticide exposure and delta BMI Z score which integrates fat accumulation from birth and onwards to school age.
Interestingly, some of the mediator marks could be linked to specific genes that were reported earlier to play a role in the development of weight gain/obesity, insulin resistance/diabetes, cardiovascular disease, and/or fetal growth retardation: UQCRC2, MTNR1B, GRIN2A, FABP4, and LRP8. FABP4 encodes for a member of the fatty acid-binding protein family regulating lipid trafficking, signaling, and metabolism. Different studies have demonstrated the role of this protein in obesity, type 2 diabetes and atherosclerosis development [91–93]. In ApoE deficient mice with hyperhomocysteine FABP4 DNA methylation is reduced in the aorta compared to wild type mice, leading to a higher gene expression [94, 95]. UQCRC2 encodes a protein which is part of the ubiquinol-cytochrome c reductase complex in the mitochondria. UQCRC2 was shown to be downregulated in individuals who were susceptible to weight gain and obesity development . The melatonin receptor 1B (MTNR1B) has a main function in regulating circadian rhythm. Interestingly, several polymorphisms in the MTNR1B gene are associated with type 2 diabetes, fasting glucose concentration, and insulin secretion [97–99]. GRIN2A encodes for a NMDA glutamate receptor subunit. Polymorphisms in the GRIN2A gene are associated with epilepsy and different neurological and mental disorders [100–104]. A decreased gene expression of GRIN2A in rats after intrauterine growth retardation suggests a possible role for this gene in fetal growth and development . LRP8 encodes for a member of the LDL receptor family. Common polymorphisms in the LRP8 gene are associated with coronary artery disease, myocardial infarction, and high birth weight [106–110]. Thus, the mediation analysis suggests a mechanistic role of epigenetics in the development of an adverse metabolic risk profile among the prenatally exposed children with the PON1 R-allele as previously reported for these children  and confirmed in the selected subset of children.
A few studies have investigated associations between PON1 genotype and metabolic disturbances in children. A recent study showed a higher risk of insulin resistance (HOMA-IR) in Mexican children with the RR-genotype as compared to children with the QQ or QR genotypes although BMI did not differ between the groups . Among Mexican-American children from an agricultural community in California, a trend of increased BMI Z scores with increased number of PON1 192Q alleles was seen . However, potential interactions between PON1 genotype and prenatal exposure to pesticides, or other environmental contaminants, were not investigated in these studies. In our cohort, unexposed QQ-homozygote children also tended to have higher body fat content than unexposed R-carriers, but prenatally pesticide exposed children with the R-allele accumulated more fat during childhood and had a more unhealthy metabolic risk profile at school age than unexposed children and exposed children with the QQ genotype .
We also demonstrated that methylation in the PON1 promoter itself is affected by a SNP (PON1 -108CT, rs705379). In addition, PON1 methylation values were negatively associated with paraoxonase 1 activity. These results are in agreement with the outcome of a recent study from Huen and colleagues . They found methylation in the same nine CpG sites to be associated with the PON1 -108CT polymorphism and also reported an inverse association with AREase activity as a measure of PON1 expression, both in newborns and 9-year-old children. Furthermore, they demonstrated that PON1 methylation mediates the relationship between PON1 expression and the promoter -108 genotype. However, the effect of prenatal pesticide exposure on the health outcomes shown in Table 1 was not modulated by PON1 -108CT genotype (data not shown).
Our findings indicate that the higher vulnerability among children with the R-allele towards prenatal pesticide exposure might be mediated by genotype-specific epigenetic alterations. However, a limitation of this study is that we cannot identify individual pesticides related to these findings, since the study design did not allow bio-monitoring of pesticide exposure in the mothers, and the exposure classification of the mothers encompassed more than 100 pesticides used in different mixtures .
However, the existence of mixed exposure is a real-world situation, and the longitudinal design, the blinded exposure classification, and the blinded clinical examinations, and genotyping minimized the possible impact of exposure misclassification and bias.
Since PON1 is known to detoxify some organophosphate insecticides (e.g., chlorpyrifos), and these substances were frequently applied in the mothers’ working areas, organophosphate insecticides could be assumed to be responsible for the observed effects. However, the mechanism is unclear and does not seem to be related to the hydrolysis efficiency, since R-carriers have higher paraoxonase activity than QQ homozygotes. Besides, at relatively low exposure levels, as in this study, the capacity to detoxify organophosphates is considered to be independent of the PON1 Q192R genotype , and furthermore, serum PON1 activity was reported to be low in newborns and may be even lower before birth, as indicated by lower activity in premature compared to term babies [115, 116]. Thus, differences in fetal detoxification of pesticides related to PON1 genotype might not be a likely explanation of the exposure-related difference in methylation pattern between children with the QR/RR and QQ genotype.
Another limitation of the study is that DNA methylation analyses were performed in white blood cells as surrogates for the target tissues. We do not know whether the differences in DNA methylation patterns found in blood mirror a similar change in adipose tissue, for example. A recent study from Huang et al. demonstrated several potential limitations in using methylation profiles in blood to mirror the corresponding profile in target tissues by comparing paired blood and adipose tissue methylation profiles . Furthermore, the composition of blood cell types may be variable and might affect the DNA methylation analyses. In our dataset, prenatal pesticide exposure and/or PON1 Q192R genotype did not affect the relative blood cell counts determined by the reference-based method of Houseman. Cell counts were not included in the models due to the small sample size of the study. Since we found that some of the health effects (mainly leptin) were associated with cell type count (Additional file 12), we cannot exclude that the results of the mediation analysis were biased by differences in cell type composition. Based on the data of Reinius et al. , methylation of only two CpG probes (cg18202502 and cg15134033) in Table 3 were slightly associated with cell types (data not shown). Methylation in the other CpG probes in Table 3 was not significantly different between the blood cell types.
Finally, the small number of subjects included in this exploratory study is a clear weakness because of the limited statistical power. Despite these constraints, our findings suggest that DNA methylation might be a link between prenatal pesticide exposure and cardio-metabolic risk profile in children carrying the PON1 192R-allele. The findings deserve further investigation in a larger study with quantitative data on pesticide exposure. Whether this DNA methylation pattern is unique to pesticide exposure or is shared by other adverse prenatal environmental factors also needs further investigation.
In summary, our data indicate that DNA methylation may be an underlying mechanism explaining an adverse cardio-metabolic risk profile in prenatally pesticide-exposed children carrying the PON1 192R-allele.
Differentially methylated posistion
Differentially methylated region
Ingenuity Pathway Analysis
Plasminogen activator inhibitor type-1
Transcription factor binding site
Bjorling-Poulsen M, Andersen HR, Grandjean P. Potential developmental neurotoxicity of pesticides used in Europe. Environ Health. 2008;7:50.
Andersen HR, Vinggaard AM, Rasmussen TH, Gjermandsen IM, Bonefeld-Jorgensen EC. Effects of currently used pesticides in assays for estrogenicity, androgenicity, and aromatase activity in vitro. Toxicol Appl Pharmacol. 2002;179:1–12.
Orton F, Rosivatz E, Scholze M, Kortenkamp A. Widely used pesticides with previously unknown endocrine activity revealed as in vitro anti-androgens. Environ Health Perspect. 2011;119:794–800.
Grandjean P, Landrigan PJ. Neurobehavioural effects of developmental toxicity. Lancet Neurol. 2014;13:330–8.
London L, Beseler C, Bouchard MF, Bellinger DC, Colosio C, et al. Neurobehavioral and neurodevelopmental effects of pesticide exposures. Neurotoxicology. 2012;33:887–96.
Li AA, Baum MJ, McIntosh LJ, Day M, Liu F, et al. Building a scientific framework for studying hormonal effects on behavior and on the development of the sexually dimorphic nervous system. Neurotoxicology. 2008;29:504–19.
Gore AC. Neuroendocrine targets of endocrine disruptors. Hormones (Athens). 2010;9:16–27.
Jacobsen PR, Christiansen S, Boberg J, Nellemann C, Hass U. Combined exposure to endocrine disrupting pesticides impairs parturition, causes pup mortality and affects sexual differentiation in rats. Int J Androl. 2010;33:434–42.
Wohlfahrt-Veje C, Main KM, Schmidt IM, Boas M, Jensen TK, et al. Lower birth weight and increased body fat at school age in children prenatally exposed to modern pesticides: a prospective study. Environ Health. 2011;10:79.
Wohlfahrt-Veje C, Andersen HR, Jensen TK, Grandjean P, Skakkebaek NE, et al. Smaller genitals at school age in boys whose mothers were exposed to non-persistent pesticides in early pregnancy. Int J Androl. 2012;35:265–72.
Andersen HR, Schmidt IM, Grandjean P, Jensen TK, Budtz-Jorgensen E, et al. Impaired reproductive development in sons of women occupationally exposed to pesticides during pregnancy. Environ Health Perspect. 2008;116:566–72.
Wohlfahrt-Veje C, Andersen HR, Schmidt IM, Aksglaede L, Sorensen K, et al. Early breast development in girls after prenatal exposure to non-persistent pesticides. Int J Androl. 2012;35:273–82.
Andersen HR, Debes F, Wohlfahrt-Veje C, Murata K, Grandjean P. Occupational pesticide exposure in early pregnancy associated with sex-specific neurobehavioral deficits in the children at school age. Neurotoxicol Teratol. 2015;47:1–9.
Costa LG, Cole TB, Furlong CE. Polymorphisms of paraoxonase (PON1) and their significance in clinical toxicology of organophosphates. J Toxicol Clin Toxicol. 2003;41:37–45.
Mackness B, Durrington P, Povey A, Thomson S, Dippnall M, et al. Paraoxonase and susceptibility to organophosphorus poisoning in farmers dipping sheep. Pharmacogenetics. 2003;13:81–8.
Aviram M, Rosenblat M, Bisgaier CL, Newton RS, Primo-Parmo SL, et al. Paraoxonase inhibits high-density lipoprotein oxidation and preserves its functions. A possible peroxidative role for paraoxonase. J Clin Invest. 1998;101:1581–90.
Durrington PN, Mackness B, Mackness MI. Paraoxonase and atherosclerosis. Arterioscler Thromb Vasc Biol. 2001;21:473–80.
Seo D, Goldschmidt-Clermont P. The paraoxonase gene family and atherosclerosis. Curr Atheroscler Rep. 2009;11:182–7.
Andersen HR, Wohlfahrt-Veje C, Dalgard C, Christiansen L, Main KM, et al. Paraoxonase 1 polymorphism and prenatal pesticide exposure associated with adverse cardiovascular risk profiles at school age. PLoS One. 2012;7:e36830.
Jorgensen A, Nellemann C, Wohlfahrt-Veje C, Jensen TK, Main KM, et al. Interaction between paraoxonase 1 polymorphism and prenatal pesticide exposure on metabolic markers in children using a multiplex approach. Reprod Toxicol. 2015;51:22–30.
Patel SB, Reams GP, Spear RM, Freeman RH, Villarreal D. Leptin: linking obesity, the metabolic syndrome, and cardiovascular disease. Curr Hypertens Rep. 2008;10:131–7.
Meas T, Deghmoun S, Chevenne D, Gaborit B, Alessi MC, et al. Plasminogen activator inhibitor type-1 is an independent marker of metabolic disorders in young adults born small for gestational age. J Thromb Haemost. 2010;8:2608–13.
Huang KC, Lin RC, Kormas N, Lee LT, Chen CY, et al. Plasma leptin is associated with insulin resistance independent of age, body mass index, fat mass, lipids, and pubertal development in nondiabetic adolescents. Int J Obes Relat Metab Disord. 2004;28:470–5.
Ho SM, Johnson A, Tarapore P, Janakiram V, Zhang X, et al. Environmental epigenetics and its implication on disease risk and health outcomes. ILAR J. 2012;53:289–305.
Faulk C, Dolinoy DC. Timing is everything: the when and how of environmentally induced changes in the epigenome of animals. Epigenetics. 2011;6:791–7.
Saffery R, Novakovic B. Epigenetics as the mediator of fetal programming of adult onset disease: what is the evidence? Acta Obstet Gynecol Scand. 2014;93:1090–8.
Reynolds RM, Jacobsen GH, Drake AJ. What is the evidence in humans that DNA methylation changes link events in utero and later life disease? Clin Endocrinol (Oxf). 2013;78:814–22.
Perera F, Herbstman J. Prenatal environmental exposures, epigenetics, and disease. Reprod Toxicol. 2011;31:363–73.
Chmurzynska A. Fetal programming: link between early nutrition, DNA methylation, and complex diseases. Nutr Rev. 2010;68:87–98.
Knopik VS, Maccani MA, Francazio S, McGeary JE. The epigenetics of maternal cigarette smoking during pregnancy and effects on child development. Dev Psychopathol. 2012;24:1377–90.
Casati L, Sendra R, Sibilia V, Celotti F. Endocrine disrupters: the new players able to affect the epigenome. Front Cell Dev Biol. 2015;3:37.
Skinner MK. Environmental epigenomics and disease susceptibility. EMBO Rep. 2011;12:620–2.
Relton CL, Groom A, St Pourcain B, Sayers AE, Swan DC, et al. DNA methylation patterns in cord blood DNA and body size in childhood. PLoS One. 2012;7:e31821.
Godfrey KM, Sheppard A, Gluckman PD, Lillycrop KA, Burdge GC, et al. Epigenetic gene promoter methylation at birth is associated with child’s later adiposity. Diabetes. 2011;60:1528–34.
Assenov Y, Muller F, Lutsik P, Walter J, Lengauer T, et al. Comprehensive analysis of DNA methylation data with RnBeads. Nat Methods. 2014;11:1138–40.
Teschendorff AE, Marabita F, Lechner M, Bartlett T, Tegner J, et al. A beta-mixture quantile normalization method for correcting probe design bias in Illumina Infinium 450 k DNA methylation data. Bioinformatics. 2013;29:189–96.
Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics. 2012;28:882–3.
Houseman EA, Accomando WP, Koestler DC, Christensen BC, Marsit CJ, et al. DNA methylation arrays as surrogate measures of cell mixture distribution. BMC Bioinformatics. 2012;13:86.
Reinius LE, Acevedo N, Joerink M, Pershagen G, Dahlen SE, et al. Differential DNA methylation in purified human blood cells: implications for cell lineage and studies on disease susceptibility. PLoS One. 2012;7:e41361.
Halachev K, Bast H, Albrecht F, Lengauer T, Bock C. EpiExplorer: live exploration and global analysis of large epigenomic datasets. Genome Biol. 2012;13:R96.
Peters TJ, Buckley MJ, Statham AL, Pidsley R, Samaras K, et al. De novo identification of differentially methylated regions in the human genome. Epigenetics Chromatin. 2015;8:6.
Obermann-Borst SA, Eilers PH, Tobi EW, de Jong FH, Slagboom PE, et al. Duration of breastfeeding and gender are associated with methylation of the LEPTIN gene in very young children. Pediatr Res. 2013;74:344–9.
Garcia-Cardona MC, Huang F, Garcia-Vivas JM, Lopez-Camarillo C, Del Rio Navarro BE, et al. DNA methylation of leptin and adiponectin promoters in children is reduced by the combined presence of obesity and insulin resistance. Int J Obes (Lond). 2014;38:1457–65.
Houde AA, Legare C, Biron S, Lescelleur O, Biertho L, et al. Leptin and adiponectin DNA methylation levels in adipose tissues and blood cells are associated with BMI, waist girth and LDL-cholesterol levels in severely obese men and women. BMC Med Genet. 2015;16:29.
Tobi EW, Lumey LH, Talens RP, Kremer D, Putter H, et al. DNA methylation differences after exposure to prenatal famine are common and timing- and sex-specific. Hum Mol Genet. 2009;18:4046–53.
Jousse C, Parry L, Lambert-Langlais S, Maurin AC, Averous J, et al. Perinatal undernutrition affects the methylation and expression of the leptin gene in adults: implication for the understanding of metabolic syndrome. FASEB J. 2011;25:3271–8.
Lesseur C, Armstrong DA, Paquette AG, Koestler DC, Padbury JF, et al. Tissue-specific Leptin promoter DNA methylation is associated with maternal and infant perinatal factors. Mol Cell Endocrinol. 2013;381:160–7.
Lesseur C, Armstrong DA, Paquette AG, Li Z, Padbury JF, et al. Maternal obesity and gestational diabetes are associated with placental leptin DNA methylation. Am J Obstet Gynecol. 2014;211(654):e651–9.
Petersen PS, Jin C, Madsen AN, Rasmussen M, Kuhre R, et al. Deficiency of the GPR39 receptor is associated with obesity and altered adipocyte metabolism. FASEB J. 2011;25:3803–14.
Janani C, Ranjitha Kumari BD. PPAR gamma gene--a review. Diabetes Metab Syndr. 2015;9:46–50.
Androutsopoulos VP, Hernandez AF, Liesivuori J, Tsatsakis AM. A mechanistic overview of health associated effects of low levels of organochlorine and organophosphorous pesticides. Toxicology. 2013;307:89–94.
Pillai HK, Fang M, Beglov D, Kozakov D, Vajda S, et al. Ligand binding and activation of PPARgamma by Firemaster(R) 550: effects on adipogenesis and osteogenesis in vitro. Environ Health Perspect. 2014;122:1225–32.
Grimaldi M, Boulahtouf A, Delfosse V, Thouennon E, Bourguet W, et al. Reporter cell lines for the characterization of the interactions between human nuclear receptors and endocrine disruptors. Front Endocrinol (Lausanne). 2015;6:62.
Khateeb J, Gantman A, Kreitenberg AJ, Aviram M, Fuhrman B. Paraoxonase 1 (PON1) expression in hepatocytes is upregulated by pomegranate polyphenols: a role for PPAR-gamma pathway. Atherosclerosis. 2010;208:119–25.
Camps J, Garcia-Heredia A, Rull A, Alonso-Villaverde C, Aragones G, et al. PPARs in regulation of paraoxonases: control of oxidative stress and inflammation pathways. PPAR Res. 2012;2012:616371.
Khateeb J, Kiyan Y, Aviram M, Tkachuk S, Dumler I, et al. Urokinase-type plasminogen activator downregulates paraoxonase 1 expression in hepatocytes by stimulating peroxisome proliferator-activated receptor-gamma nuclear export. Arterioscler Thromb Vasc Biol. 2012;32:449–58.
Divers J, Palmer ND, Lu L, Register TC, Carr JJ, et al. Admixture mapping of coronary artery calcified plaque in African Americans with type 2 diabetes mellitus. Circ Cardiovasc Genet. 2013;6:97–105.
Parks BW, Nam E, Org E, Kostem E, Norheim F, et al. Genetic control of obesity and gut microbiota composition in response to high-fat, high-sucrose diet in mice. Cell Metab. 2013;17:141–52.
Fox CS, Liu Y, White CC, Feitosa M, Smith AV, et al. Genome-wide association for abdominal subcutaneous and visceral adipose reveals a novel locus for visceral fat in women. PLoS Genet. 2012;8:e1002695.
Baron RM, Kenny DA. The moderator-mediator variable distinction in social psychological research: conceptual, strategic, and statistical considerations. J Pers Soc Psychol. 1986;51:1173–82.
Tingley D, Yamamoto T, Hirose K, Keele L, Imai K. mediation: R package for causal mediation analysis. J Stat Soft. 2014;59:1-38.
Pinero J, Queralt-Rosinach N, Bravo A, Deu-Pons J, Bauer-Mehren A, et al. DisGeNET: a discovery platform for the dynamical exploration of human diseases and their genes. Database (Oxford). 2015;2015:bav028.
Verhulst PJ, Lintermans A, Janssen S, Loeckx D, Himmelreich U, et al. GPR39, a receptor of the ghrelin receptor family, plays a role in the regulation of glucose homeostasis in a mouse model of early onset diet-induced obesity. J Neuroendocrinol. 2011;23:490–500.
Zhang JV, Ren PG, Avsian-Kretchmer O, Luo CW, Rauch R, et al. Obestatin, a peptide encoded by the ghrelin gene, opposes ghrelin’s effects on food intake. Science. 2005;310:996–9.
Cui Y, Ying Y, van Hasselt A, Ng KM, Yu J, et al. OPCML is a broad tumor suppressor for multiple carcinomas and lymphomas with frequently epigenetic inactivation. PLoS One. 2008;3:e2990.
Sellar GC, Watt KP, Rabiasz GJ, Stronach EA, Li L, et al. OPCML at 11q25 is epigenetically inactivated and has tumor-suppressor function in epithelial ovarian cancer. Nat Genet. 2003;34:337–43.
Wu Y, Davison J, Qu X, Morrissey C, Storer B, et al. Methylation profiling identified novel differentially methylated markers including OPCML and FLRT2 in prostate cancer. Epigenetics. 2016;11:247–58.
Li C, Tang L, Zhao L, Li L, Xiao Q, et al. OPCML is frequently methylated in human colorectal cancer and its restored expression reverses EMT via downregulation of smad signaling. Am J Cancer Res. 2015;5:1635–48.
Volkow ND, Wang GJ, Baler RD. Reward, dopamine and the control of food intake: implications for obesity. Trends Cogn Sci. 2011;15:37–46.
Murray S, Tulloch A, Gold MS, Avena NM. Hormonal and neural mechanisms of food reward, eating behaviour and obesity. Nat Rev Endocrinol. 2014;10:540–52.
Vucetic Z, Kimmel J, Totoki K, Hollenbeck E, Reyes TM. Maternal high-fat diet alters methylation and gene expression of dopamine and opioid-related genes. Endocrinology. 2010;151:4756–64.
Paternain L, Batlle MA, De la Garza AL, Milagro FI, Martinez JA, et al. Transcriptomic and epigenetic changes in the hypothalamus are involved in an increased susceptibility to a high-fat-sucrose diet in prenatally stressed female rats. Neuroendocrinology. 2012;96:249–60.
Ong ZY, Gugusheff JR, Muhlhausler BS. Perinatal overnutrition and the programming of food preferences: pathways and mechanisms. J Dev Orig Health Dis. 2012;3:299–308.
Palmer AA, Brown AS, Keegan D, Siska LD, Susser E, et al. Prenatal protein deprivation alters dopamine-mediated behaviors and dopaminergic and glutamatergic receptor binding. Brain Res. 2008;1237:62–74.
Vucetic Z, Totoki K, Schoch H, Whitaker KW, Hill-Smith T, et al. Early life protein restriction alters dopamine circuitry. Neuroscience. 2010;168:359–70.
Wright TM, Fone KC, Langley-Evans SC, Voigt JP. Exposure to maternal consumption of cafeteria diet during the lactation period programmes feeding behaviour in the rat. Int J Dev Neurosci. 2011;29:785–93.
Teegarden SL, Scott AN, Bale TL. Early life exposure to a high fat diet promotes long-term changes in dietary preferences and central reward signaling. Neuroscience. 2009;162:924–32.
Torres-Altoro MI, Mathur BN, Drerup JM, Thomas R, Lovinger DM, et al. Organophosphates dysregulate dopamine signaling, glutamatergic neurotransmission, and induce neuronal injury markers in striatum. J Neurochem. 2011;119:303–13.
Slotkin TA. Does early-life exposure to organophosphate insecticides lead to prediabetes and obesity? Reprod Toxicol. 2011;31:297–301.
Richard D, Huang Q, Timofeeva E. The corticotropin-releasing hormone system in the regulation of energy balance in obesity. Int J Obes Relat Metab Disord. 2000;24 Suppl 2:S36–9.
Richard D, Lin Q, Timofeeva E. The corticotropin-releasing factor family of peptides and CRF receptors: their roles in the regulation of energy balance. Eur J Pharmacol. 2002;440:189–97.
Sharma R, Banerji MA. Corticotropin releasing factor (CRF) and obesity. Maturitas. 2012;72:1–3.
Uehara Y, Shimizu H, Ohtani K, Sato N, Mori M. Hypothalamic corticotropin-releasing hormone is a mediator of the anorexigenic effect of leptin. Diabetes. 1998;47:890–3.
Sutton AK, Pei H, Burnett KH, Myers Jr MG, Rhodes CJ, et al. Control of food intake and energy expenditure by Nos1 neurons of the paraventricular hypothalamus. J Neurosci. 2014;34:15306–18.
Leshan RL, Greenwald-Yarnell M, Patterson CM, Gonzalez IE, Myers Jr MG. Leptin action through hypothalamic nitric oxide synthase-1-expressing neurons controls energy balance. Nat Med. 2012;18:820–3.
Sansbury BE, Hill BG. Regulation of obesity and insulin resistance by nitric oxide. Free Radic Biol Med. 2014;73:383–99.
Naseh M, Vatanparast J, Baniasadi M, Hamidi GA. Alterations in nitric oxide synthase-expressing neurons in the forebrain regions of rats after developmental exposure to organophosphates. Neurotoxicol Teratol. 2013;37:23–32.
Naseh M, Vatanparast J. Enhanced expression of hypothalamic nitric oxide synthase in rats developmentally exposed to organophosphates. Brain Res. 2014;1579:10–9.
Snodgrass-Belt P, Gilbert JL, Davis FC. Central administration of transforming growth factor-alpha and neuregulin-1 suppress active behaviors and cause weight loss in hamsters. Brain Res. 2005;1038:171–82.
Ennequin G, Boisseau N, Caillaud K, Chavanelle V, Etienne M, et al. Neuregulin 1 affects leptin levels, food intake and weight gain in normal-weight, but not obese, db/db mice. Diabetes Metab. 2015;41:168–72.
Wu G, Li H, Zhou M, Fang Q, Bao Y, et al. Mechanism and clinical evidence of lipocalin-2 and adipocyte fatty acid-binding protein linking obesity and atherosclerosis. Diabetes Metab Res Rev. 2014;30:447–56.
Hotamisligil GS, Bernlohr DA. Metabolic functions of FABPs—mechanisms and therapeutic implications. Nat Rev Endocrinol. 2015;11:592–605.
Furuhashi M, Saitoh S, Shimamoto K, Miura T. Fatty acid-binding protein 4 (FABP4): pathophysiological insights and potent clinical biomarker of metabolic and cardiovascular diseases. Clin Med Insights Cardiol. 2014;8:23–33.
Jiang Y, Ma S, Zhang H, Yang X, Lu GJ, et al. FABP4-mediated homocysteine-induced cholesterol accumulation in THP-1 monocyte-derived macrophages and the potential epigenetic mechanism. Mol Med Rep. 2016;14:969-76.
Yang AN, Zhang HP, Sun Y, Yang XL, Wang N, et al. High-methionine diets accelerate atherosclerosis by HHcy-mediated FABP4 gene demethylation pathway via DNMT1 in ApoE(-/-) mice. FEBS Lett. 2015;589:3998–4009.
Marrades MP, Gonzalez-Muniesa P, Arteta D, Martinez JA, Moreno-Aliaga MJ. Orchestrated downregulation of genes involved in oxidative metabolic pathways in obese vs. lean high-fat young male consumers. J Physiol Biochem. 2011;67:15–26.
Mussig K, Staiger H, Machicao F, Haring HU, Fritsche A. Genetic variants in MTNR1B affecting insulin secretion. Ann Med. 2010;42:387–93.
Nagorny C, Lyssenko V. Tired of diabetes genetics? Circadian rhythms and diabetes: the MTNR1B story? Curr Diab Rep. 2012;12:667–72.
Karamitri A, Renault N, Clement N, Guillaume JL, Jockers R. Minireview: Toward the establishment of a link between melatonin and glucose homeostasis: association of melatonin MT2 receptor variants with type 2 diabetes. Mol Endocrinol. 2013;27:1217–33.
Liu R, Dang W, Du Y, Zhou Q, Liu Z, et al. Correlation of functional GRIN2A gene promoter polymorphisms with schizophrenia and serum D-serine levels. Gene. 2015;568:25–30.
Endele S, Rosenberger G, Geider K, Popp B, Tamer C, et al. Mutations in GRIN2A and GRIN2B encoding regulatory subunits of NMDA receptors cause variable neurodevelopmental phenotypes. Nat Genet. 2010;42:1021–6.
Yoo HJ, Cho IH, Park M, Yang SY, Kim SA. Family based association of GRIN2A and GRIN2B with Korean autism spectrum disorders. Neurosci Lett. 2012;512:89–93.
Carvill GL, Regan BM, Yendle SC, O’Roak BJ, Lozovaya N, et al. GRIN2A mutations cause epilepsy-aphasia spectrum disorders. Nat Genet. 2013;45:1073–6.
Lal D, Steinbrucker S, Schubert J, Sander T, Becker F, et al. Investigation of GRIN2A in common epilepsy phenotypes. Epilepsy Res. 2015;115:95–9.
Ninomiya M, Numakawa T, Adachi N, Furuta M, Chiba S, et al. Cortical neurons from intrauterine growth retardation rats exhibit lower response to neurotrophin BDNF. Neurosci Lett. 2010;476:104–9.
Wang L, Wang X, Laird N, Zuckerman B, Stubblefield P, et al. Polymorphism in maternal LRP8 gene is associated with fetal growth. Am J Hum Genet. 2006;78:770–7.
Shen GQ, Li L, Wang QK. Genetic variant R952Q in LRP8 is associated with increased plasma triglyceride levels in patients with early-onset CAD and MI. Ann Hum Genet. 2012;76:193–9.
Shen GQ, Girelli D, Li L, Olivieri O, Martinelli N, et al. Multi-allelic haplotype association identifies novel information different from single-SNP analysis: a new protective haplotype in the LRP8 gene is against familial and early-onset CAD and MI. Gene. 2013;521:78–81.
Shen GQ, Girelli D, Li L, Rao S, Archacki S, et al. A novel molecular diagnostic marker for familial and early-onset coronary artery disease and myocardial infarction in the LRP8 gene. Circ Cardiovasc Genet. 2014;7:514–20.
Guo T, Yin RX, Yao LM, Huang F, Pan L, et al. Integrative mutation, haplotype and G x G interaction evidence connects ABGL4, LRP8 and PCSK9 genes to cardiometabolic risk. Sci Rep. 2016;6:37375.
Alegria-Torres JA, Garcia-Dominguez ML, Cruz M, Aradillas-Garcia C. Q192R polymorphism of paraoxonase 1 gene associated with insulin resistance in Mexican children. Arch Med Res. 2015;46:78–83.
Huen K, Harley K, Beckman K, Eskenazi B, Holland N. Associations of PON1 and genetic ancestry with obesity in early childhood. PLoS One. 2013;8:e62565.
Huen K, Yousefi P, Street K, Eskenazi B, Holland N. PON1 as a model for integration of genetic, epigenetic, and expression data on candidate susceptibility genes. Environ Epigenet. 2015;1:dvv003.
Coombes RH, Meek EC, Dail MB, Chambers HW, Chambers JE. Human paraoxonase 1 hydrolysis of nanomolar chlorpyrifos-oxon concentrations is unaffected by phenotype or Q192R genotype. Toxicol Lett. 2014;230:57–61.
Cole TB, Jampsa RL, Walter BJ, Arndt TL, Richter RJ, et al. Expression of human paraoxonase (PON1) during development. Pharmacogenetics. 2003;13:357–64.
Costa LG, Giordano G, Cole TB, Marsillach J, Furlong CE. Paraoxonase 1 (PON1) as a genetic determinant of susceptibility to organophosphate toxicity. Toxicology. 2013;307:115–22.
Huang YT, Chu S, Loucks EB, Lin CL, Eaton CB, et al. Epigenome-wide profiling of DNA methylation in paired samples of adipose tissue and blood. Epigenetics. 2016;11:227–36.
We are grateful to the families for their participation in the greenhouse cohort study. We thank Mariann Bøllund and the greenhouse cohort study team for the skilled help with child examinations and the database. We thank Karen Hollanders (VITO, Belgium) for DNA extraction from blood samples.
The study was funded by The Danish Environmental Protection Agency (project number 667-00164). The funding organization had no role in study design, data collection and analysis, interpretation of the results, or preparation of the manuscript.
Availability of data and materials
The datasets generated and/or analyzed during the current study are available in the Gene Expression Omnibus (GEO) repository with accession number GSExxxx.
KD, SR, GS, WVB, and HRA conceived and designed the experiments. KD, SR, and HRA performed the experiments and analyzed the data. CWV, KMM, GVC, GS, WVB, and HRA contributed reagents/materials/analysis tools. KD, SR, GS, WVB, and HRA wrote the paper. KD, SR, CWV, KMM, GVC, GS, WVB, and HRA evaluated the manuscript text. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
Ethics approval and consent to participate
The study was conducted according to the Helsinki II Declaration with written informed consent by all parents and oral consent by the children as approved by The Regional Scientific Ethical Committees for Southern Denmark (S-20070068) and the Danish Data Protection Agency.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Primer sequences. (XLSX 10 kb)
PCA before and after batch effect correction for Sentrix_ID and processing date using ComBat. (TIFF 240 kb)
Associations between the first eight principal components and covariates before and after ComBat batch correction. Associations between principal components and Sentrix_ID were measured using the Kruskal-Wallis test. Associations between principal components and processing date, exposure and PON1 Q192R genotype were measured using the two-sided Wilcoxon sum rank test. (TIFF 102 kb)
Genomic location of the pyrosequencing assays represented as a UCSC genome browser track. The first track indicates the sequence analyzed by pyrosequencing (Seq_to_analyse). Other custom tracks include: CpG islands, Dnase I hypersensitivity clusters, H3K27ac histone marks, transcription factor-binding sites, and the Illumina 450 K methylation probes. A) LEP assay B) GPR39 assay C) PPARG assay and D) OPCML assay. (TIFF 3014 kb)
Relative cell type contribution estimated by the Houseman approach. Differences in cell type composition between the exposure groups were measured using one-way ANOVA. (TIFF 259 kb)
Sig-DMPs overlapping with DMRs (interaction model). For each DMP P values are given for the interaction between pesticide exposure and PON1 genotype (P.Value.int_EXP:PON1), and for the exposure effect PON1 R-allele carrier group (P.Value.EXP_when PON1 QR/RR). The mean beta values in each exposure group are listed. (XLSX 82 kb)
Sig-DMPs overlapping with DMRs (model without interaction). For each DMP P values are given for the PON1 effect (P.Value.PON1) and exposure effect (P.Value.EXP). The mean beta values in each exposure group are listed. (XLSX 16 kb)
Outcome of GPR39 DMR pyrosequencing. Boxplots showing methylation differences between the exposure groups in the GPR39 pyrosequencing region. P values shown are those of the exposure effect. (TIFF 193 kb)
Outcome of GPR39 DMR pyrosequencing. (XLSX 11 kb)
Enrichment of TFBS for DMPs significant in the interaction model. P value from the Fisher’s exact test were adjusted using the Bonferroni correction. (XLSX 10 kb)
Mediation analysis. Outcome statistics and gene disease associations of (partial) mediators between pesticide exposure and body fat measures in PON1 R-allele carriers. (XLSX 47 kb)
Association between estimated blood cell counts and health outcomes. Simple linear regression was used to determine associations between the relative blood cell type composition and the health outcomes (body fat, BMI Z score, delta BMI Z score, leptin levels, and age). (XLSX 11 kb)
About this article
- DNA methylation
- Prenatal pesticide exposure
- Paraoxonase 1
- PON1 Q192R genotype
- Illumina 450 K methylation array
- Cardio-metabolic health