Placental Expression of Imprinted Genes, Overall and in Sex-Specific Patterns, Associated with Placental Cadmium Concentrations and Birth Size.

BACKGROUND: Prenatal cadmium (Cd) exposure has been recognized to restrict growth, and male and female fetuses may have differential susceptibility to the developmental toxicity of Cd. Imprinted genes, which exhibit monoallelic expression based on parent of origin, are highly expressed in placental tissues. The function of these genes is particularly critical to fetal growth and development, and some are expressed in sex-specific patterns.
OBJECTIVES: We aimed to examine whether prenatal Cd associates with the expression of imprinted placental genes, overall or in fetal sex-specific patterns, across two independent epidemiologic studies.
METHODS: We tested for Cd–sex interactions in association with gene expression, then regressed the placental expression levels of 74 putative imprinted genes on placental log-Cd concentrations while adjusting for maternal age, sex, smoking history, and educational attainment. These models were performed within study- and sex-specific strata in the New Hampshire Birth Cohort Study (NHBCS; [Formula: see text]) and the Rhode Island Child Health Study (RICHS; [Formula: see text]). We then used fixed-effects models to estimate the sex-specific and overall associations across strata and then examine heterogeneity in the associations by fetal sex.
RESULTS: We observed that higher Cd concentrations were associated with higher expression of distal-less homeobox 5 (DLX5) ([Formula: see text]), and lower expression of h19 imprinted maternally expressed transcript (H19) ([Formula: see text]) and necdin, MAGE family member (NDN) ([Formula: see text]) across study and sex-specific strata, while three other genes [carboxypeptidase A4 (CPA4), growth factor receptor bound protein 10 (GRB10), and integrin-linked kinase (ILK)] were significantly associated with Cd concentrations, but only among female placenta ([Formula: see text]). Additionally, the expression of DLX5, H19, and NDN, the most statistically significant Cd-associated genes, were also associated with standardized birth weight z-scores. DISCUSSION: The differential regulation of a set of imprinted genes, particularly DLX5, H19, and NDN, in association with prenatal Cd exposure may be involved in overall developmental toxicity, and some imprinted genes may respond to Cd exposure in a manner that is specific to infant gender. https://doi.org/10.1289/EHP4264.

Introduction

Cadmium (Cd) is an ubiquitous environmental contaminant that bioaccumulates in plants grown in contaminated soils, resulting in human exposure when those plants are consumed (Khan et al. 2017). Though Cd is most commonly recognized for its roles in kidney damage, bone and joint problems, and various cancers (ATSDR 2012), it is also a developmental and reproductive toxicant (Thompson and Bannigan 2008). Maternal Cd not transferred to the fetus tends to accumulate in the placenta, making the placenta a useful biomarker of maternal Cd exposures and burden during pregnancy (Piasek et al. 2014). Cd increases oxidative stress, interrupts cell cycles, induces apoptosis and cell death (Rani et al. 2014), and acts as an endocrine disruptor, which may be particularly important for reproductive and developmental toxicity, even at relatively low exposure levels (Knazicka et al. 2015). Thus, despite the protective role of the placenta in limiting direct fetal exposure to Cd, toxic effects on this critical developmental organ may in turn elicit adverse impacts on fetal development.

Numerous epidemiologic studies of prenatal Cd exposure have observed associations with anthropometric measures at birth (Al-Saleh et al. 2014, 2015; Wang et al. 2016) and/or fetal growth restriction (Llanos and Ronco 2009; Wang et al. 2018), impaired cognition in childhood (Kippler et al. 2012b), preeclampsia (Laine et al. 2015; Wang et al. 2018), and impaired zinc (Zn) transport (Kippler et al. 2010). Additionally, multiple studies have implicated epigenetic dysregulation (Everson et al. 2018; Kippler et al. 2013; Mohanty et al. 2015; Sanders et al. 2014) and/or altered expression of genes related to oxidative and inflammatory response (Everson et al. 2016, 2017) as components underlying how Cd might impact the functions of developmental tissues and lead to these adverse pregnancy and birth outcomes.

Recent publications have also suggested that imprinted genes, which have monoallelic expression patterns dependent on the parent of origin, may be particularly susceptible to Cd-associated differential methylation and/or expression. Xu et al. (2017) used a mouse model to identify that the placental expression of paternally expressed gene 10 (PEG10) can be significantly down-regulated and cyclin dependent kinase inhibitor 1C (CDKN1C) significantly up-regulated in response to higher prenatal Cd exposures (Xu et al. 2017). Cowley et al. (2018) performed a genome-wide study in human maternal and cord blood samples to identify differentially methylated regions (DMRs) that associated with prenatal Cd exposure, and found that Cd-associated DMRs were significantly enriched for imprinting control regions (ICR) in both maternal and fetal blood (Cowley et al. 2018). These findings are particularly important because imprinted genes play critical roles in the fetal growth and development, and their dysregulation can result in numerous developmental disorders (Cassidy and Charalambous 2018).

Though studies have begun to implicate imprinted genes as potential targets of prenatal Cd toxicity, evidence of this in humans is still quite limited. In fact, most studies in human populations have focused on exposure associated variation in the DNA methylation of ICRs and not on the functional expression of those genes. The current study aims to address these gaps by examining how placental Cd concentrations are associated with the placental expression of both known and putative imprinted genes across two large independent U.S. samples, the New Hampshire Birth Cohort Study (NHBCS, https://www.dartmouth.edu/∼childrenshealth/) and the Rhode Island Child Health Study (RICHS, https://www.rhodeislandkidshealth.com/). There are also data suggesting that the influence of Cd exposure on molecular activities, as well as pregnancy and birth outcomes, differ by the sex of the fetus. Cd-associated reductions in anthropometric measures at birth appeared to be more common among female infants in two prior studies (Kippler et al. 2012a; Taylor et al. 2016), and fetal epigenetic response to prenatal Cd exposure also appeared to differ by fetal sex in one study (Kippler et al. 2013). Thus, to further elucidate whether prenatal Cd elicits sex-specific responses, we assessed whether Cd-associated variations in placental imprinted expression differed by fetal sex.

Methods

New Hampshire Birth Cohort Study

The NHBCS is an ongoing birth cohort initiated in 2009 (Gilbert-Diamond et al. 2016). This study includes pregnant women who used an unregulated private well as a primary source of drinking water, were between 18 and 45 y of age, and that attended one of the study clinics in New Hampshire for prenatal care. All participants provided written informed consent in accordance with the requirements of the Committee for the Protection of Human Subjects at Dartmouth College. For the current study, we included placental samples from participants who were enrolled between February 2012 and September 2013, for whom placenta were selected for genetic, epigenetic, and exposure biomarker assays. From this sample, we included all participants who had complete data on imprinted gene expression and placental concentrations of Cd (). Sociodemographic, lifestyle, and anthropometric data were obtained via self-administered questionnaires and medical record abstraction.

Rhode Island Child Health Study

The RICHS includes mother–infant pairs with healthy and successful pregnancies at the Women and Infants’ Hospital in Providence, Rhode Island, United States, that were between September 2010 and February 2013 (Appleton et al. 2016). Mothers younger than 18 y of age, pregnancies resulting in preterm birth, or infants born with congenital or chromosomal abnormalities were excluded from the RICHS. To examine contrasts between babies born large, adequate, or small for gestational age (LGA, AGA, or SGA), infants born SGA ( birth weight percentile) or LGA ( birth weight percentile) were oversampled, then infants born AGA (between the 10th and 90th birth weight percentiles) were enrolled who matched on gestational age and maternal age. All protocols were approved by the institutional review boards at the Women and Infants Hospital of Rhode Island, Dartmouth College, and Emory University, and all participants provided written informed consent. The current study included all mother–infant pairs for whom placental imprinted gene expression and Cd concentrations were performed (). Sociodemographic, lifestyle, and anthropometric data were obtained via interviewer-administered questionnaires and medical record abstraction.

Covariates

We adjusted for variables that have been associated with Cd exposure and that are recognized as influencing pregnancy and birth outcomes. Smoking status, highest educational attainment, and age were observed to be the strongest predictors of urinary Cd concentrations in a recent study of women and/or children across 16 European countries (Berglund et al. 2015). Highest educational attainment was defined as mothers who obtained an education beyond the high school level vs. those who did not go beyond high school. Maternal smoking during pregnancy (MSDP) was defined as any self-reported smoking during any point of pregnancy vs. those who reported no smoking during pregnancy. Maternal age was included as a continuous variable. All models were stratified on fetal sex and were adjusted for smoking status during pregnancy, highest educational attainment, and maternal age. In RICHS, we also adjusted for maternal race as a dichotomous variable (white vs. nonwhite) since this sample was racially diverse, whereas in NHBCS, the sample was homogenous for race, and additional adjustments were not necessary.

Placental RNA Sampling

Both cohorts followed a similar placental sampling protocol, obtaining biopsies within 2 h of delivery from the fetal side adjacent to the cord insertion site and free of maternal decidua. Samples were placed in RNAlater (Life Technologies), then rinsed, frozen, and stored at . RNA was extracted with the RNAeasy mini kit (#74,106; Qiagen) and subsequently stored at . RNA quality and integrity were assessed via the NanoDrop™ ND 1000 spectrophotometer (Thermo Fisher) and the Agilent Bioanalyzer 2100 (Agilent).

Imprinted Gene Expression

A set of 108 candidate known or putative imprinted genes were measured via a custom nCounter® (NanoString Technologies) panel in RICHS; the probe selection, design, and normalization methods are described in detail elsewhere (Kappil et al. 2015). Of these genes, 74 were also measured in the NHBCS using the same technology and normalization techniques; annotations for these genes are provided in Excel Table S1. Briefly, the raw nCounter® data were normalized against the geometric mean of positive controls, and negative controls were used to determine background noise levels. Samples with expression values that were indistinguishable from noise for more than half of the candidate genes, or those with expression levels standard deviations for more than 10 genes were excluded. Normalized counts were transformed, and batch effects were removed using ComBat from the sva package (Leek et al. 2012). We also compared the expression levels measured via NanoString nCounter to those measured via RNA-Seq, which we have previously published on (Everson et al. 2018); 69 of the 74 genes mapped to the same gene IDs in the RNAseq data. Overall, 56 of the 69 genes exhibited positive and statistically significant correlations between mRNA measured via NanoString and RNA-Seq (σ ranged from 0.14 to 0.75, with a mean rho of 0.40).

Cadmium Quantification

Placental concentrations of Cd were quantified in both cohorts at the Dartmouth Trace Elements Analysis Core via inductively coupled plasma mass spectrometry; details of the processing are described elsewhere (Punshon et al. 2016). Only three NHBCS samples and none of the RICHS samples yielded undetectable Cd concentrations. Samples with nondetectable concentrations were assigned a value equal to the minimum detectable values. Placental Cd concentrations were log transformed to better approximate a normal distribution.

Statistical Analyses

All statistical analyses were conducted in R (version 3.4.4; R Project, http://www.R-project.org/). We used the limma package with standard errors estimated via an empirical Bayes method (Ritchie et al. 2015) to assess the relationships between log-Cd and placental imprinted expression. First, we tested the for the presence of an interaction between log-Cd and fetal sex by regressing imprinted gene expression on log-Cd, sex, and an interaction term (). Overall estimates of interaction across RICHS and NHBCS were estimated with inverse variance–weighted fixed-effects models using the metafor package (Viechtbauer 2010), and models yielding interaction terms with were determined to be statistically significant. The linear relationships between log-Cd and -imprinted gene expression were then assessed for all 74 genes within four separate strata (by study and sex). Sex-specific associations, estimated across NHBCS and RICHS, and overall associations, estimated across all four strata, were estimated with inverse variance–weighted fixed-effects models using the metafor package (Viechtbauer 2010). To control for multiple comparisons, we implemented the Benjamini-Hochberg false discovery rate (FDR) adjustment and considered associations that yielded FDR-adjusted to be statistically significant. Venn diagrams were produced to show the concordance and the discordance in Cd-associated imprinted gene expression between male and female samples. Those genes that yielded at least nominally significant interaction terms (raw ) and FDR-significant associations within male or female strata were determined to have sex-specific associations with Cd. Those genes that yielded FDR-significant overall associations were determined to have overall relationships with Cd.

We performed a sensitivity analysis to further interrogate the potential confounding effects of maternal smoking. For this analysis, we reproduced the cohort- and sex-specific models between log-Cd and -imprinted gene expression, while excluding all mothers who reported any smoking during pregnancy. We examined the correlation between regression coefficients from models that included and excluded smokers to determine the impact of MSDP on our overall findings. We performed a second sensitivity analysis to test for potential interactions with maternal race, which was restricted to the RICHS sample. For this analysis, we generated two additional models, first regressing imprinted gene expression on log-Cd and maternal race (white vs. nonwhite) while including an interaction term () for these variables, and second additional model testing for the three-way interaction with log-Cd, sex, and race ().

We additionally explored whether any of the Cd-associated genes identified in the current study produced nominally significant () associations at epigenetic loci from our prior epigenome-wide association study (EWAS) (Everson et al. 2018). We used a Fisher’s exact test to examine whether the number of nominally significant loci within these genes was greater than expected, and we considered those genes with to be enriched for Cd-associated differential methylation. For this test, we included all CpGs that were of the start and end coordinates of the imprinted genes that were associated with Cd, and we were thus focused on the Cd-associated variations in cis-acting DNA methylation, as opposed to trans-acting CpGs or variations in the ICR for which consensus regions have not been defined for some of the genes that we are studying. We used all other loci as the background level of nominal significance.

Finally, we estimated overall and sex-specific associations between imprinted gene expression with z-scores for birth weight, birth length, and head circumference at the Cd-associated genes from the above analysis. z-Scores were calculated via Fenton growth curves and were standardized by sex and gestational age (Fenton and Kim 2013). Estimates of association were obtained for either sex-specific or nonsex-specific strata based on whether the Cd models yielded sex-specific or nonsex-specific associations. These estimates of association were obtained via inverse variance–weighted fixed-effects models using the metafor package.

Results

The RICHS cohort had slightly higher median Cd concentrations than the NHBCS, with females in RICHS having the highest reported concentrations () and females in the NHBCS having the lowest (), while the overall distributions of Cd concentrations between males and females were similar. The RICHS participants also tended to have lower educational attainment, higher proportions of MSDP, and greater racial/ethnic heterogeneity, while the distributions of these variables within study-specific strata were similar across the sexes, and all four strata were very similar in terms of average maternal age and gestational duration (Table 1).

Table 1

Frequencies (percentages) and of demographic characteristics and the distributions of placental Cd concentrations from the New Hampshire Birth Cohort Study (NHBCS) and the Rhode Island Child Health Study (RICHS), stratified by male and female newborns.

Maternal and fetal characteristics stratified by study and fetal sex		NHBCS (n=326)		RICHS (n=211)
		Female (n=153)	Male (n=173)	Female (n=100)	Male (n=111)
Maternal educational attainment	>HS	138 (90)	153 (88)	75 (75)	88 (79)
	≤HS	15 (10)	20 (12)	25 (25)	23 (21)
Maternal smoking during pregnancy	None	141 (92)	148 (86)	74 (74)	77 (69)
	Any	12 (8)	25 (15)	26 (26)	34 (31)
Maternal race	Nonwhite	0 (0)	0 (0)	28 (28)	25 (23)
	White	153 (100)	173 (100)	72 (72)	86 (78)
Gestational weeks	Mean±SD	39.41±1.41	39.35±1.69	39.26±0.97	39.3±0.97
Maternal age (years)	Mean±SD	31.39±4.83	31.68±4.99	30.37±5.74	29.56±5.52
Placental Cd (ng/g)	Minimum	0.19	0.19	1.22	1.06
	25th percentile	2.07	2.14	2.94	2.79
	Median	3.00	3.27	4.38	3.97
	Mean±SD	3.47±1.94	3.74±2.57	4.65±2.65	4.44±2.44
	75th percentile	4.74	4.62	5.27	5.43
	Maximum	10.43	22.31	16.42	17.99

Note: Cd, cadmium; HS, High school; NHBCS, New Hampshire Birth Cohort Study; RICHS, Rhode Island Child Health Study. Maternal variables including education, smoking behavior, race, and age were collected prior to delivery for NHBCS and postdelivery but prior to discharge for RICHS.

Tests for Interactions between Cd and Fetal Sex

First, we tested for potential interactions between sex and log-Cd while adjusting for maternal age, maternal educational attainment, and MSDP in both cohorts, and adjusted for maternal race/ethnicity in RICHS (Excel Table S2). We identified five imprinted genes that yielded nominally significant interaction terms (): integrin-linked kinase (ILK), carboxypeptidase A4 (CPA4), sarcoglycan epsilon (SGCE), growth factor receptor bound protein 10 (GRB10), and zinc finger and BTB domain containing 8B (ZBTB8B).

Study- and Sex-Stratified Results

We then assessed the linear relationships between placental log-Cd and imprinted expression levels using linear models, stratified by study and fetal sex (Excel Table S3). Among males, the strongest association in NHBCS was observed at coatomer protein complex subunit gamma 2 (COPG2) [; 95% confidence interval (CI): 0.069, 0.654; ], which did not yield an association in RICHS (; 95% CI: , 0.468; ), while the strongest association in RICHS was observed at necdin, MAGE family member (NDN) (; 95% CI: , ; ), which did yield a similar association in NHBCS (; 95% CI: , 0.008; ). Among females, the strongest association in NHBCS was observed at the imprinted maternally expressed transcript, H19 (; 95% CI: , ; ), which did not yield an association in RICHS but did yield somewhat attenuated associations in both male-specific strata (NHBCS: ; 95% CI: , 0.035; ; RICHS: ; 95% CI: , ; ). The strongest association among females in RICHS was observed at SGCE (; 95% CI: 0.172, 0.766; ), which did not yield similar associations within the other strata. Overall, we observed a greater proportion of nominally significant associations among females [NHBCS: 17/74 (23%), RICHS: 13/74 (18%)] when compared with males [NHBCS: 4/74 (5%), RICHS: 8/74 (11%)].

Meta-Analysis for Sex-Specific and Overall Associations

We then estimated sex-specific associations for males and females, meta-analyzed across studies, and identified nine nominally associated genes among males (Excel Table S4), two of which yielded FDR-significant associations: NDN and distal-less homeobox 5 (DLX5). We identified 17 nominally associated genes among females (Excel Table S5), four of which yielded FDR-significant associations: GRB10, ILK, CPA4, and thrombospondin type 1 domain containing 7A (THSD7A). Those genes that yielded nominally significant interaction terms and FDR-significant associations with log-Cd within the male or female strata were defined as having sex-specific associations with Cd. These included GRB10, ILK, and CPA4, all of which were specific to female placenta (Table 2). Three genes were differentially expressed with placental Cd concentrations at nominal significance levels among both males and females: DLX5, H19, and insulin-like growth factor 2 antisense 1 (IGF2-AS) (Figure 1A). We then estimated overall associations, meta-analyzed across all four strata for all 74 genes (Excel Table S6), and found that seven genes yielded FDR-significant associations with log-Cd (Table 3), with the most statistically significant Cd-associated genes being DLX5 (; 95% CI: 0.207, 0.566; ), H19 (; 95% CI: , ; ), and NDN (; 95% CI: , ; ). We examined whether the Cd-associated genes exhibited absolute differential expression between male and female placentae, and found that their overall expression did not differ by sex (Figure 1B; Excel Table S7). There were numerous moderate-to-strong correlations in the expression levels of these Cd-associated genes, with the strongest positive correlations between H19 and IGF2 (; ) and the strongest inverse correlation between delta-like noncanonical notch ligand 1 (DLK1) and SGCE (; ) (Figure 1C). We then examined the consistency in the associations across strata using forest plots of the FDR-significant sex-specific (Figure 2) and overall associations (Figure 3). Among the genes yielding overall associations with Cd, DLX5 demonstrated the greatest homogeneity in its association with log-Cd across all four strata, consistently yielding parameter estimates between 0.239 and 0.484, with only one CI overlapping the null, while H19 and IGF2-AS yielded more heterogenous associations with Cd across study- and sex-specific strata.

Table 2

Parameter estimates for sex-specific associations between placental log-cadmium (Cd) and gene expression levels that were adjusted for maternal age, maternal educational attainment, and maternal smoking during pregnancy in all models and additionally adjusted for maternal race/ethnicity in RICHS. This includes genes that yielded interaction and FDR-significant associations within meta-analyzed male- or female-specific strata (FDR ). Beta coefficients represent the sex-specific inverse variance–weighted fixed-effects estimates across NHBCS and RICHS for males and females, and estimates for each gene according to study and sex. Results for all other genes can be found in Excel Tables S3, S4, and S5.

Gene	Strata		Model results
	Study	Sex	Beta	p-Value	FDR q-value	95% CI
CPA4	Sex specific	Females	−0.56	0.00085	0.021	−0.89, −0.23
		Males	−0.015	0.92	0.975	−0.30, 0.27
	NHBCS	Females	−0.44	0.057	—	−0.89, 0.014
		Males	−0.083	0.67	—	−0.46, 0.30
	RICHS	Females	−0.70	0.0054	—	−1.19, −0.21
		Males	0.079	0.73	—	−0.37, 0.53
GRB10	Sex specific	Females	0.42	0.00028	0.018	0.19, 0.64
		Males	0.060	0.59	0.899	−0.16, 0.28
	NHBCS	Females	0.41	0.042	—	0.015, 0.80
		Males	−0.0032	0.99	—	−0.37, 0.36
	RICHS	Females	0.42	0.0035	—	0.14, 0.70
		Males	0.10	0.49	—	−0.18, 0.38
ILK	Sex specific	Females	0.23	0.00049	0.018	0.10, 0.36
		Males	0.015	0.84	0.942	−0.13, 0.16
	NHBCS	Females	0.22	0.023	—	0.032, 0.42
		Males	0.040	0.72	—	−0.18, 0.26
	RICHS	Females	0.24	0.010	—	0.058, 0.42
		Males	−0.0035	0.97	—	−0.20, 0.19

Note: —, no data; Cd, cadmium; CI, confidence interval; CPA4, carboxypeptidase A4; FDR, false discovery rate; GRB10, growth factor receptor bound protein 10; ILK, integrin-linked kinase; NHBCS, New Hampshire Birth Cohort Study; RICHS, Rhode Island Child Health Study.

Figure 1.

(A) Venn diagram showing the overlap in nominally significant associations across females (orange) and males (blue) [data for cadmium (Cd)-associated expression for males and females are available in Excel Tables S4 and S5, respectively], as well as (B) the average expression levels ( counts) among male and female placenta (data for expression levels among male and female placenta are available in Excel Table S7), and (C) correlations (, ) in the expression patterns across these cadmium (Cd)-associated genes.

Table 3

Parameter estimates for overall associations between placental log-Cd and gene expression levels that were adjusted for maternal age, maternal educational attainment, and maternal smoking during pregnancy in all models and additionally adjusted for maternal race/ethnicity in RICHS. This includes genes that yielded FDR-significant associations when meta-analyzed across all four sex- and study-specific specific strata (FDR ). Beta coefficients represent the inverse variance–weighted fixed-effects estimates across sex- and study-specific specific strata, and estimates for each gene according to study and sex. Results for all other genes can be found in Excel Tables S3, S4, and S6.

Gene	Strata		Model results
	Study	Sex	Beta	p-Value	FDR q-value	95% CI
DLX5	Overall		0.39	0.000025	0.0018	0.21, 0.57
	NHBCS	Females	0.48	0.012	—	0.11, 0.86
		Males	0.42	0.024	—	0.056, 0.78
	RICHS	Females	0.24	0.17	—	−0.10, 0.58
		Males	0.43	0.023	—	0.061, 0.81
H19	Overall		−0.34	0.00027	0.010	−0.52, −0.16
	NHBCS	Females	−0.91	0.000083	—	−1.35, −0.45
		Males	−0.37	0.073	—	−0.78, −0.035
	RICHS	Females	0.012	0.94	—	−0.30, 0.32
		Males	−0.40	0.17	—	−0.74, −0.056
NDN	Overall		−0.37	0.00064	0.016	−0.59, −0.16
	NHBCS	Females	−0.47	0.035	—	−0.90, −0.033
		Males	−0.46	0.054	—	−0.93, 0.0076
	RICHS	Females	0.10	0.66	—	−0.34, 0.54
		Males	−0.60	0.0027	—	−0.99, −0.21
IGF2-AS	Overall		−0.35	0.00097	0.018	−0.56, −0.14
	NHBCS	Females	−0.54	0.0081	—	−0.93, −0.14
		Males	−0.39	0.041	—	−0.77, −0.017
	RICHS	Females	−0.025	0.93	—	−0.58, 0.53
		Males	−0.28	0.18	—	−0.70, 0.13
IGF2	Overall		−0.36	0.0019	0.029	−0.59, −0.13
	NHBCS	Females	−0.79	0.0041	—	−1.33, −0.26
		Males	−0.54	0.036	—	−1.04, −0.035
	RICHS	Females	0.029	0.89	—	−0.40, 0.46
		Males	−0.36	0.093	—	−0.78, 0.061
GRB10	Overall		0.24	0.0034	0.042	0.078, 0.39
	NHBCS	Females	0.41	0.042	—	0.015, 0.79
		Males	−0.0031	0.98	—	−0.37, 0.36
	RICHS	Females	0.42	0.0035	—	0.14, 0.70
		Males	0.10	0.49	—	−0.18, 0.38
THSD7A	Overall		0.36	0.0041	0.043	0.12, 0.61
	NHBCS	Females	0.39	0.13	—	−0.12, 0.89
		Males	0.27	0.29	—	−0.23, 0.77
	RICHS	Females	0.69	0.0062	—	0.20, 1.17
		Males	0.084	0.74	—	−0.42, 0.58

Note: —, no data; Cd, cadmium; CI, confidence interval; DLX5, distal-less homeobox 5; FDR, false discovery rate; GRB10, growth factor receptor bound protein 10; H19, h19 imprinted maternally expressed transcript; IGF2, insulin-like growth factor 2; IGF2-AS, insulin-like growth factor 2 antisense 1; NDN, necdin, MAGE family member; NHBCS, New Hampshire Birth Cohort Study; RICHS, Rhode Island Child Health Study; THSD7A, thrombospondin type 1 domain containing 7A.

Figure 2.

Forest plot of associations between log-Cd and gene expression for CPA4, GRB10, and ILK within study- and sex-specific strata (, ), along with meta-analysis estimates for female- and male-specific expression associations (data for study and sex-specific associations are available in Excel Table S3, while inverse variance–weighted fixed-effects estimates for males and females are available in Excel Tables S4 and S5, respectively). Note: Cd, cadmium; Coef., regression coefficient; CPA4, carboxypeptidase A4; GRB10, growth factor receptor bound protein 10; ILK, integrin-linked kinase; NHBCS, New Hampshire Birth Cohort Study; RICHS, Rhode Island Child Health Study.

Figure 3.

Forest plot of associations between log-Cd and gene expression for DLX5, H19, NDN, IGF2-AS, IGF2, GRB10, and THSD7A within study- and sex-specific strata (, ), along with meta-analysis estimates across all four strata (black diamonds) (data for study and sex-specific associations are available in Excel Table S3, while inverse variance–weighted fixed-effects estimates are available in Excel Table S6). Note: Cd, cadmium; CI, confidence interval; Coef., regression coefficient; DLX5, distal-less homeobox 5; GRB10, growth factor receptor bound protein 10; H19, h19 imprinted maternally expressed transcript; IGF2, insulin-like growth factor 2; IGF2-AS, insulin-like growth factor 2 antisense 1; NDN, necdin, MAGE family member; NHBCS, New Hampshire Birth Cohort Study; RICHS, Rhode Island Child Health Study; THSD7A, thrombospondin type 1 domain containing 7A.

Sensitivity Analysis: Reproducibility of Results among Nonsmokers and Interactions with Maternal Race

We performed a sensitivity analysis to assess whether our results were robust to the exclusion of samples from those that reported any MSDP. The coefficients across the 74 genes and all four models were very highly correlated (; ; all ) (Figure S1). The estimated overall and sex-specific associations yielded almost identical results after excluding smokers, and all seven genes that yielded FDR-significant associations with Cd in the original models again yielded FDR-significant associations when restricting the analysis to nonsmokers (Excel Table S8). For instance, the overall estimate of association between log-Cd and DLX5 expression among nonsmokers (; 95% CI: 0.184, 0.595; ) was strikingly similar to the association from the original model in which smokers were included (; 95% CI: 0.207, 0.566; ). The estimated associations became modestly stronger among mothers who did not smoke for H19 (; 95% CI: , ; ), for NDN (; 95% CI: , ; ), and for IGF2-AS (; 95% CI: , ; ). Among the sex-specific associations that we identified, the female-specific relationships between log-Cd with ILK and CPA4 were modestly attenuated among nonsmokers, although these association were still substantially larger among female placenta vs. male placenta, and these female-specific associations were still nominally significant ().

We also tested whether the relationships between log-Cd and gene expression may have differed by race in the RICHS sample and whether potential interactions between Cd and race may have influenced our findings. We found that three genes had nominally significant interactions (; Excel Table S9): cysteine-rich angiogenic inducer 61 (CYR61; also known as cellular communication network factor 1; CCN1), cyclin-dependent kinase inhibitor 1C (CDKN1C), and COPG2. We additionally tested for three-way interactions between log-Cd, race, and sex, and found the SGCE and CDKN1C exhibited statistically significant interactions (; Excel Table S10).

Test of Enrichment for Cd-Associated Differential Methylation

We explored whether the genes identified in this analysis yielded associations between placental Cd and DNA methylation in our previously published EWAS (Everson et al. 2018). From that prior analysis, 442 epigenetic loci were identified that were of the start and end coordinates for DLX5, H19, IGF2, IGF2-AS, NDN, CPA4, GRB10, ILK, and THSD7A, which had been assessed for associations between placental Cd and DNA methylation. None of these loci were identified as being differentially methylated at the FDR or suggestive significance threshold as part of that EWAS, although 9% (42 loci) of these loci did yield nominally significant associations (Excel Table S11). All but three of these 42 nominally significant results were annotated to four genes, and these three genes were significantly enriched, via Fisher’s exact test, for Cd-associated differential methylation: DLX5 (), IGF2-AS (), IGF2 (), and GRB10 (). For DLX5, the majority of Cd-associated DNA methylation occurred upstream of the transcription start site, while for IGF2, IGF2-AS, and GRB10, these loci were not as spatially correlated and resided across multiple different potential regulatory elements, including the transcription start sites, gene body, as well as the 3′ and 5′ untranslated regions.

Associations between Imprinted Genes and Birth Metrics and Placental Dimensions

Finally, we tested whether the imprinted genes that associated with placental Cd concentrations also associated with z-scores for birth weight, birth length, and head circumference (Table 4). Most notable, the top three Cd-associated genes, DLX5, H19, and NDN, were associated with multiple birth size metrics. Higher expression of DLX5 was associated with smaller birth weight (; 95% CI: , ; ) and smaller head circumference (; 95% CI: , ; ), while higher expression of H19 was associated with larger birth size (; 95% CI: 0.014, 0.082; ) and longer birth length (; 95% CI: 0.013, 0.087; ). Higher expression of NDN was also associated with larger birth weight (; 95% CI: 0.041, 0.125; ) and longer birth length (; 95% CI: 0.027, 0.117; ). Interestingly, many of the imprinted genes that appeared to have female-specific Cd-associated differential expression did not tend to yield female-specific associations with birth size. In fact, higher expression levels of both CPA4 and ILK were associated with larger z-scores for birth weight, birth length, and head circumference, but only among male placenta. We also tested for associations between imprinted gene expression and placental weight (grams), while adjusting for gestational age, fetal sex, maternal education, maternal age, and maternal smoking with linear models in the NHBCS () (Excel Table S12). We only observed a significant association () between CPA4 with placental weight among males.

Table 4

Parameter estimates of overall and sex-specific associations between placental gene expression levels with z-scores for birth weight, birth length, and head circumference that were adjusted for maternal age, maternal educational attainment, and maternal smoking during pregnancy in all models and additionally adjusted for maternal race/ethnicity in RICHS; these analyses were performed within the strata in which these genes yielded FDR-significant associations (FDR ) with log-Cd.

Gene	Strata	BW z-scores		BL z-scores		HC z-scores
		Beta (95% CI)	p-Value	Beta (95% CI)	p-Value	Beta (95% CI)	p-Value
DLX5	Overall	−0.050 (−0.086, −0.014)	0.0070	−0.030 (−0.069, 0.0081)	0.12	−0.045 (−0.079, −0.011)	0.0089
H19	Overall	0.048 (0.014, 0.082)	0.0063	0.050 (0.013, 0.087)	0.0083	0.025 (−0.0094, 0.060)	0.15
NDN	Overall	0.083 (0.041, 0.13)	0.000097	0.072 (0.027, 0.12)	0.0017	0.015 (−0.024, 0.054)	0.46
IGF2-AS	Overall	0.021 (−0.24, 0.066)	0.36	0.016 (−0.029, 0.062)	0.49	0.014 (−0.025, 0.053)	0.47
IGF2	Overall	0.037 (−0.0075, 0.081)	0.10	0.044 (−0.0042, 0.092)	0.074	0.019 (−0.025, 0.063)	0.39
THSD7A	Overall	0.0049 (−0.046, 0.055)	0.85	−0.015 (−0.068, 0.038)	0.59	−0.0057 (−0.052, 0.040)	0.81
CPA4	Female	−0.023 (−0.095, 0.048)	0.53	−0.063 (−0.14, 0.010)	0.092	−0.0080 (−0.060, 0.044)	0.76
	Male	0.078 (0.019, 0.14)	0.0095	0.062 (0.00093, 0.12)	0.047	0.080 (0.013, 0.15)	0.019
GRB10	Female	−0.0099 (−0.055, 0.035)	0.66	0.0024 (−0.045, 0.050)	0.92	−0.00051 (−0.039, 0.038)	0.98
	Male	0.017 (−0.025, 0.058)	0.43	0.044 (−0.0014, 0.090)	0.057	0.0061 (−0.043, 0.055)	0.81
ILK	Female	0.0075 (−0.020, 0.035)	0.59	0.0074 (−0.021, 0.036)	0.62	0.0034 (−0.018, 0.025)	0.75
	Male	0.042 (0.014, 0.069)	0.0027	0.044 (0.014, 0.074)	0.0036	0.041 (0.0094, 0.073)	0.011

Note: —, no data; BL, birth length; BW, birth weight; Cd, cadmium; CI, confidence interval; CPA4, carboxypeptidase A4; DLX5, distal-less homeobox 5; FDR, false discovery rate; GRB10, growth factor receptor bound protein 10; HC, head circumference; H19, h19 imprinted maternally expressed transcript; IGF2, insulin-like growth factor 2; IGF2-AS, insulin-like growth factor 2 antisense 1; ILK, integrin-linked kinase; NDN, necdin, MAGE family member; RICHS, Rhode Island Child Health Study; THSD7A, thrombospondin type 1 domain containing 7A.

Discussion

We identified Cd-associated variations in the expression of six imprinted genes, DLX5, H19, NDN, IGF2-AS, IGF2, and THSD7A across studies and within both sexes, as well as sex-specific associations for CPA4, GRB10, and ILK. Among our top three hits, higher Cd concentrations were associated with higher expression of DLX5, lower expression of H19, and lower expression of NDN. Additionally, higher expression of DLX5 was associated with smaller birth weight and smaller head circumference, while lower expression of H19 was associated with smaller birth weight and shorter birth length, and lower expression of NDN was associated with smaller birth weight and smaller birth length. The placental expression of some of the genes identified in this study have also been associated with early life growth and cognitive outcomes in prior research from our group but had previously not been investigated for associations with environmental exposures. For instance, NDN, H19, and IGF2 were observed to be more highly expressed from placental tissues of babies born LGA compared with those born AGA or SGA, and in a factor analysis, these three genes all loaded onto the same factor that was also associated with birth weight, which could be an indication that they are involved in a coexpression pattern that may regulate fetal growth (Kappil et al. 2015). Additionally, the placental expression of DLX5 loaded onto a different factor that was also strongly associated with birth size (Kappil et al. 2015), and DLX5 expression was selected as the third most important imprinted gene for classifying infant neurobehavioral profiles (Green et al. 2015). These prior findings provide important context to our results, since prenatal Cd exposure associates with decreased anthropometric measures at birth (Al-Saleh et al. 2014, 2015; Wang et al. 2016) and/or fetal growth restriction (Llanos and Ronco 2009; Wang et al. 2018), as well as with impaired cognition and or neurobehavioral outcomes in childhood (Gustin et al. 2018; Kippler et al. 2016).

Our findings of Cd-associated variations in placental expression may represent important intermediates between prenatal Cd exposure and developmental outcomes. Experimental and epidemiologic studies have demonstrated the critical functions that some these genes play in growth and development. We observed inverse associations between Cd and H19, IGF2, and IGF2-AS, an mRNA that is expressed antisense to IGF2 that tends to be coregulated with IGF2 (Duart-Garcia and Braunschweig 2014). H19 and IGF2 are some of the most well-studied imprinted genes due to their coregulation and competing roles in growth and development. The IGF2 gene produces a growth factor that promotes both placental and fetal development, and is involved in nutrient/waste transfer between the mother and the fetus, while H19, which is transcribed but not translated into a protein, functions as an inhibitor of IGF2 and thus represses placental and fetal growth (Nordin et al. 2014). Variations in fetal genotype, DNA methylation, and expression of H19/IGF2 have been observed to be strongly associated with fetal growth (St-Pierre et al. 2012; Su et al. 2016). On the other hand, the functions of IGF2-AS have not been as thoroughly studied, and its roles in growth and development are not well established. Experimental studies have found that lack of functional IGF2-AS may contribute to fetal growth restriction due to its status as a paternally expressed transcript (Duart-Garcia and Braunschweig 2014), and inhibition of IGF2-AS has been suggested to promote angiogenesis (Zadora et al. 2017) and protect neuronal cells against apoptosis (Song et al. 2017) via up-regulation of IGF2 expression. Our study observed that H19, IGF2-AS, and IGF2 may have lower expression levels in placenta that have higher Cd concentrations, which may indicate Cd-associated disrupted activity of these coregulated genes or an adaptive placental response to Cd exposure. Other epidemiologic studies that have examined the relationships between maternal Cd exposure and fetal epigenetic responses at candidate ICR did not observe Cd-associated differential methylation at the H19/IGF2 ICR (Cowley et al. 2018; Vidal et al. 2015), which is not consistent with our findings. However, these studies examined DNA methylation rather than expression, and their molecular measures were obtained from maternal and cord blood rather than placental tissues. It is possible that Cd-associated responses differ across fetal tissues and that epigenomic and transcriptomic responses vary.

While the coregulation of H19 and IGF2 is a well-recognized phenomenon, other imprinted genes can also be coregulated via a shared ICR or potentially as part of a larger trans-regulatory mechanism, similar to the imprinted gene network (IGN), which has been thoroughly characterized in mice (Patten et al. 2016). Early studies of the IGN found that altering the expression of a few imprinted genes, ZAC1 or H19, could influence the expression patterns of a number of other imprinted genes (Gabory et al. 2009; Varrault et al. 2006). It has been suggested that the IGN acts as a compensatory mechanism to regulate appropriate fetal growth, which is supported by experimental mouse models in which the majority of the IGN in the placenta was up-regulated in response to assisted reproductive technologies, but the pregnancies produced phenotypically normal embryos (Fauque et al. 2010). In human placenta, the transcription factor PLAGL1 has been shown to be coexpressed with a number of other imprinted genes, most notably H19 and IGF2, suggesting that there may also be an IGN for human placenta and that perturbations to the human placental IGN is associated with growth restriction (Iglesias-Platas et al. 2014). Although the expression of PLAGL1 was not associated with Cd in our study, both H19 and IGF2 were. We also observed numerous moderate to strong correlations between the imprinted genes that were associated with Cd. Thus, it is possible that our findings are related to larger perturbations to the placental IGN rather than the Cd-associated responses of individual imprinted genes.

We observed the most statistically significant and most consistent relationships between Cd and expression of the imprinted DLX5 gene. Of note, although this gene was not identified as being differentially methylated in response to Cd in our prior EWAS (Everson et al. 2018), approximately 20% of the CpGs annotated to DLX5 were associated with Cd at a nominal significance level (), providing further evidence for the relationship between its expression control and Cd exposure. DLX5 is a transcription factor that is primarily recognized for its involvement in bone growth and repair. DLX5 expression from placental cells appears to decrease with increasing gestational time (Novakovic et al. 2017), and loss of imprinting at DLX5 leads to increased expression, which has been shown to be up-regulated in the placentas of preeclamptic pregnancies (Zadora et al. 2017). Additionally, expression of DLX5 plays a critical role in neurogenesis (Perera et al. 2004) and the development of the olfactory and GnRH systems (Garaffo et al. 2015). Interestingly, DLX5 may work in concert with necdin, produced from the NDN gene, to regulate the differentiation of neuronal cells (Kuwajima et al. 2006). The NDN gene, which was differentially expressed in association with Cd in our study, plays important roles in repressing the cell cycle and inhibiting cellular growth, and is primarily expressed in post-mitotic neurons (Taniguchi et al. 2000). Deletion of NDN leads to numerous defects in the axonal migration, arborization, and growth of some neuronal cells (Pagliardini et al. 2005), and NDN is a candidate gene for Prader-Willi syndrome, a neurobehavioral condition characterized by poor growth, feeding issues, developmental delay, respiratory problems, and learning disabilities (Cheon 2016).

Overall, we found that female placentae tended to be more likely to exhibit Cd-associated differential expression compared with male placentae. The CPA4 gene, which is located within a carboxypeptidase gene cluster on chromosome 7, was inversely expressed in association with log-Cd only among female infants. This gene is primarily recognized for its role in carcinogenesis and cancer progression, appearing to be overexpressed in several types of cancer (Sun et al. 2016a, 2016b). Although the functional role of CPA4 in placental tissues is largely unknown, it has been observed to be primarily expressed from the maternal allele in fetal tissues and suggested as a potential candidate gene for Silver-Russell syndrome (Bentley et al. 2003). In other tissues, CPA4 has been implicated to be involved in inhibiting adipogenesis and modulating in insulin sensitivity (He et al. 2016). Expression of the integrin-linked kinase (ILK) gene was positively correlated with placental Cd concentrations among females. Placental expression of the ILK gene promotes trophoblast syncytialization (Butler et al. 2017), while deletion of ILK leads to multiple vascular pathologies and placental insufficiency (Friedrich et al. 2004), demonstrating the critical role of this gene in placental vascular development. We also found the GRB10 was more highly expressed in female placenta with higher Cd concentrations. Interestingly, a GRB10 ICR has been observed to be differentially methylated in response to Cd exposure, although this association was observed in maternal blood rather than fetal tissues (Cowley et al. 2018). The expression of placental GRB10 has been shown to inhibit placental growth and reduce placental efficiency in experimental models (Charalambous et al. 2010). This gene has also been shown to regulate energy homeostasis (Liu et al. 2014) and thus has implications for pathologies related to lipid metabolism, thermogenesis, and adipogenesis.

Additionally, we observed an overall association with THSD7A, whose expression was up-regulated in association with Cd. Although this gene did not yield a statistically significant interaction term, the magnitudes of association were more pronounced among female placenta in our study, and this gene had previously been identified to be differentially methylated in the cord blood of female infants in association with maternal Cd exposure, but not male infants (Kippler et al. 2013). Given that we did observe more pronounced Cd associations among female placenta in both RICHS and NHBCS, it is possible that we were merely underpowered to detect a statistically significant interaction for this gene. The Cd-associated differences in the activity and regulation of THSD7A, generated from two independent studies, using different molecular markers (DNA methylation and gene expression) that were measured in different tissues (cord blood and placenta), are quite striking and deserve additional study to elucidate the role that this gene may play in Cd-associated developmental and reproductive toxicity, especially for female infants. THSD7A has been shown to be highly expressed in placental vasculature, particularly at the leading end of human umbilical vein endothelial cells (Wang et al. 2010), and thus could play an important role in the appropriate vascularization of the placenta. Additionally, zebrafish models have demonstrated that knockdown of THSD7A can cause angiogenic defects (Wang et al. 2011) and that it plays critical roles in promoting the development of the nervous and vascular systems (Liu et al. 2016). THSD7A has also been shown to be expressed in human trophoblast subtypes and associated with trophoblast invasion, while it appears to be down-regulated among placentae that are complicated by severe preeclampsia and down-regulated with hypoxia (Luo et al. 2016).

Some epidemiologic studies have found that female infants appear to be more susceptible to Cd-associated reductions in anthropometric measures (Kippler et al. 2012a; Taylor et al. 2016). Although our study focused on the relationships between Cd and imprinted expression patterns rather than anthropometry, we did observe a greater number of nominally significant and FDR-significant associations between log-Cd and imprinted expression among the placenta of female infants. Thus, the Cd-associated molecular and functional responses in the placenta may be more pronounced among female infants, which may, in part, explain their increased susceptibility to growth restriction.

Cigarette smoke is recognized to be a primary source of cadmium exposure among smokers and has been shown to alter placental morphology (Jauniaux and Burton 2007) and associate with perturbations to the placental epigenome and/or transcriptome (Bruchova et al. 2010; Morales et al. 2016). Although we adjusted for self-reported smoking during pregnancy in all models, we also performed a sensitivity analysis by excluding all samples from which the mothers reported any smoking during pregnancy. Among nonsmokers, we observed the same associations between Cd and imprinted gene expression, with very high correlations between regression coefficients before and after excluding smokers. Thus, our findings are likely independent of the potential confounding effects of smoking during pregnancy. We also tested whether race may act as an effect modifier, rather than just as a potential confounder of the associations between imprinted gene expression and placental Cd. None of the genes that yielded FDR-significant associations with Cd (DLX5, H19, NDN, IGF2-AS, IGF2, THSD7A, CPA4, GRB10, or ILK) had statistically significant interactions with maternal race. However, a subset of imprinted genes did exhibit associations with placental Cd that may differ among racial subgroups in the RICHS study: CYR61, CDKN1C, COPG2, and SGCE. These analyses could only compare the associations between white vs. nonwhite mothers due to small sample sizes within the specific racial subgroups. Additionally, these models could not be tested in a meta-analytic framework due to NHBCS being racially homogenous. Thus, further study with larger sample sizes and more heterogenous populations is needed to better understand whether and how placental imprinted genes respond to Cd within racial and ethnic subgroups.

Some of the observed associations yielded more homogeneous associations with Cd, while others were fairly heterogenous within strata of sex and study site. Among the genes that yielded overall associations with Cd in, only DLX5 yielded homogenous results across all four strata. Additionally, the associations between DLX5 expression and birth weight z-scores were the most homogenous among the analyses of birth size. We were unable to disentangle why the RICHS female subgroup produced null results, while all three of the other strata yielded relatively strong inverse relationships with H19, IGF2-AS, and IGF2. This same pattern of associations (null for RICHS females and strong for all other strata) was also observed for NDN. Additionally, THSD7A did not meet our criteria for a sex-specific association with Cd, although it did have an interaction p-value () close to our threshold and produced consistently weak or null associations among males and substantially stronger positive associations among females. Because of this heterogeneity, the associations for these genes needs to be interpreted with caution and explored further for both sex-specific and overall relationships with prenatal Cd exposure. We have also previously studied the associations between placental Cd and placental DNA methylation (Everson et al. 2018). Although the imprinted genes from the current study were not differentially methylated at an FDR threshold for statistical significance in that study, we showed that four genes (DLX5, IGF2, IGF2-AS, and GRB10) were enriched for nominal levels of differential methylation, providing an additional layer of evidence that the regulation of these genes may be influenced by placental Cd concentrations.

In summary, we identified Cd-associated gene expression in placental tissue in nine genes. These associations appear to be specific to female placenta for CPA4, GRB10, and ILK, while the associations with DLX5, H19, NDN, IGF2-AS, IGF2, and THSD7A were not sex specific. These findings were independent of MSDP and are consistent with prior work showing an apparent greater susceptibility to the effects of Cd among female fetuses. We found that the expression of some of these Cd-associated genes were also associated with birth weight, birth length, and head circumference. Higher placental expression of DLX5, which was the top hit from the Cd association models, was observed to be inversely associated with birth weight and birth length. These findings provide mechanistic insights into how Cd may elicit some of its toxic effects on growth and development via perturbed expression of placental imprinted genes.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.