Study population. The NBDPS recruits cases from population-based, active surveillance congenital anomaly registries in nine U.S. states and includes live births and stillbirths > 20 weeks gestation or at least 500 g, as well as elective terminations of prenatally diagnosed defects when available (Yoon et al. 2001 (
link)). Arkansas, Iowa, and Massachusetts ascertain cases statewide, whereas California, Georgia, New York, North Carolina, Texas, and Utah ascertain cases in select counties. Cases are reviewed by clinical geneticists using standardized study protocols to determine study eligibility and classification, and cases with chromosomal/microdeletion disorders and disorders of known single-gene deletion causation are excluded. Controls are unaffected livebirths who are randomly selected from vital records or hospital records, depending upon study center. The NBDPS has been approved by the institutional review boards (IRBs) of all participating centers, and all participants provided written or oral informed consent before participation. These analyses were reviewed and approved by the University of North Carolina IRB.
For this analysis, the study population consisted of all controls and eligible cases with a simple, isolated CHD (i.e., a single CHD with no extra-cardiac birth defects present) and an estimated date of delivery (i.e., due date) from 1 October 1997 through 31 December 2006. During this time period, the participation response was 69% among all cases and 65% for controls. Within the NBDPS, a team of clinicians with expertise in pediatric cardiology reviewed information abstracted from the medical record and centrally assigned a single, detailed cardiac phenotype to each case whose diagnosis was confirmed by echocardiography, cardiac catheterization, surgery, or autopsy and documented in the medical record. Phenotypes were then aggregated into individual CHDs and defect groupings (Botto et al. 2007 (
link)). The following additional groups were created because of limited sample size of individual defects:
a) other conotruncal defects, which included common truncus, interrupted aortic arch–type B (IAA-type B), interrupted aortic arch–not otherwise specified (IAA-NOS), double outlet right ventricle associated with transposition of the great arteries (DORV-TGA) and not associated with TGA (DORV-other), and conoventricular septal defects (VSD-conoventricular); and
b) atresias that included both pulmonary and tricuspid atresia. Simple, isolated CHD cases represented 64% (
n = 12,383) of the total CHD cases. We restricted the analysis to offspring with a single CHD to create more etiologically homogeneous case groups, although this limits the generalizability of our findings. Women who reported having pregestational diabetes (types 1 and 2) during their pregnancy were excluded owing to the established association with CHD (Correa et al. 2008 ). Women living > 50 km from a pollutant-specific air monitor were excluded from that analysis.
Exposure assignment. Each woman reported the due date that was provided by her clinician during pregnancy to obtain the gestational age of the infant at birth. Using the gestational age to estimate the date of conception, we assigned calendar dates to each week of pregnancy. Women’s residential addresses during pregnancy were centrally geocoded to ensure consistency across study centers. Each geocoded address during weeks 2–8 of pregnancy was matched to the closest air monitor for each pollutant, with > 50% of the data available using ArcGISv10 (ESRI, Redlands, CA) and monitor locations obtained from U.S. EPA’s Air Quality System (U.S. EPA 2013 ). Participants from 1996–1998 were excluded from the analysis of PM
2.5 because monitoring began in 1999.
We used the daily maximum hourly measurement for CO, NO
2, and SO
2, the daily maximum 8-hr average for O
3, and 24-hr measurements of PM
10 and PM
2.5 to assign exposure. We averaged over the daily maximum or 24-hr measurements for weeks 2–8 of pregnancy to assign a 7-week and also 1-week averages of the daily values. We included week 2 in addition to the standard window of cardiac development, because of the potential for lag effects of air pollution (van den Hooven et al. 2012 (
link)). If only a single measurement was taken during a given week, it was assigned as the weekly exposure. Ambient levels of each pollutant except O
3 were categorized into the following categories, using the distribution of pollutant concentration among controls: less than the 10th centile (referent), 10th centile to less than the median, the median to less than the 90th centile, and greater than or equal to the 90th centile. These categories captured the departure from linearity observed in initial, exploratory analyses (data not shown). For similar reasons, O
3 was categorized into quartiles. Centiles were calculated separately for the 7-week and 1-week measures of exposure.
Statistical analysis. The following variables obtained from the maternal interview were identified as potential confounders through directed acyclic graph analysis (Greenland et al. 1999 (
link)) and included in the final adjustment set: maternal age, race/ethnicity, educational attainment, household income, tobacco smoking in the first month of pregnancy, alcohol consumption during the first trimester, and maternal nativity. Maternal age was represented as a single, continuous term, measured at the time of conception. Race/ethnicity was self-reported and categorized into the following groups: white non-Latino, black non-Latino, Latino, Asian or Pacific Islander, and other. Educational attainment was collapsed into six categories: 0–6 years of education, 7–11 years, high school graduate or equivalency, 1–3 years of college or trade school, 4 years of college or completion of a bachelor’s degree, and an advanced degree. Household income was self-reported as < $10,000 annually, > $50,000 annually, or in-between. We adjusted for any tobacco use in the first month of pregnancy and differentiated between some alcohol consumption (less than four drinks) and binge drinking (four or more drinks) during the first trimester. Maternal nativity was defined as self-report of being born outside the United States.
To account for potential differences in case ascertainment by study center, models were also adjusted for the center-specific ratio of septal defects to total CHDs. Identifying septal defects often depends on method of case ascertainment (Martin et al. 1989 (
link)). All potential confounders, as well as distance to major roadway, prepregnancy body mass index (BMI), and maternal occupation status during pregnancy were assessed for effect measure modification by constructing logistic regression models with and without interaction terms and conducting likelihood ratio tests using an
a priori alpha level of 0.1. Distance to the closest major road—defined as an interstate, U.S. highway, state, or larger county highway—was constructed using ArcGISv10 and then dichotomized at 50 m. Prepregnancy BMI was defined using self-reported maternal height and weight and categorized according to National Institutes of Health (1998) guidelines into underweight (BMI < 18.5), normal weight (18.5 ≤ BMI < 25), overweight (25 ≤ BMI < 30), and obese (BMI ≥ 30). Maternal occupation status was defined as ever working outside the home during any time during pregnancy.
For each pollutant, models were constructed to explore individual defects and defect-groupings. If a woman did not have at least one monitoring value for each week of exposure, she was excluded from the weekly analysis. We explored the relationships between all weeks and all defects because of uncertainty in pregnancy dating when using an estimated date of conception and the lack of clearly elucidated mechanisms by which cardiac development could be disrupted by exposure to air pollution. Animal research suggests that exposures outside the typical period of development for an individual heart structure could also be etiologically relevant (Morgan et al. 2008 (
link)).
Because we simultaneously assessed multiple weeks of exposure and multiple defects/groupings, we constructed two-stage hierarchical regression models to account for the correlation between estimates and partially address multiple inference (Greenland 1992 (
link); Witte et al. 1998 (
link)). The first-stage, represented in
Equation 1, was an unconditional, polytomous logistic regression model of individual CHDs on exposure (
x) defined as either all 1-week averages of maximum or 24-hr pollutant values or the single 7-week average, and the full adjustment set (
w) detailed above.
bd is the vector of regression coefficients corresponding to pollutant exposure for an individual CHD (
d),
cd is the vector of regression coefficients corresponding to the covariates for a given defect, and
m is the total number of individual types of CHDs. The second-stage model, which defines how the first-stage betas are associated, is given in Equation 2:
β
i =
Zir + δ
i, [2]
where
Zi is a row in the design matrix that includes an intercept term and then indicator variables for type of defect, broader defect grouping, and exposure week/level for the
ith β,
r is the vector of coefficients corresponding to the variables included in the design matrix, and δ
i are independent normal random variables with a mean of zero and a variance of τ
2 that describe the residual variation in β
i. The obtained second-stage coefficients,
r, are used to estimate values toward which the first-stage coefficients will be shrunk, with the magnitude of the shrinkage depending on the precision of the maximum-likelihood estimate obtained in stage 1 and the value of the second-stage variance, τ
2 (Greenland 1992 (
link); Witte et al. 1998 (
link)). We fixed τ
2 at 0.5, corresponding to a prior belief with 95% certainty that the residual odds ratio (OR) will fall within a 16-fold span.
To assess whether our results were robust to changes in model specification, we conducted sensitivity analyses by setting the value of τ
2 to 0.25, corresponding to a 7-fold OR span, as well as to a value of 1, corresponding to a 50-fold span. We also explored different specifications for the design matrix, in turn defining the prior value as a common mean for all defects, a common mean for each defect, or a common mean for each exposure week/level, across defects. Individual defects with > 10 but < 100 cases were excluded from hierarchical models and explored using Firth’s penalized maximum-likelihood method to address the quasi-complete separation that occurred due to small sample size (Heinze and Schemper 2002 (
link)). These defects included the individual defects collapsed into the other conotruncals and atresia categories described above; Ebstein’s anomaly, which was part of the right ventricular outflow tract obstruction (RVOTO) defect grouping; and muscular ventricular septal defects (VSD
muscular), which was part of the septal defect-grouping. IAA-type A and partial anomalous pulmonary venous return had < 10 cases each and were excluded from all individual analyses, but were included in the left ventricular outflow tract obstruction (LVOTO) and anomalous pulmonary venous return (APVR) defect groupings, respectively. To assess whether pollutant–defect relationships conformed to a monotonic dose response, we reanalyzed the data using incremental coding which compares each category of exposure to its immediate predecessor. If the incremental ORs are all above (or below) 1, the relationship conforms to a monotonic dose response (Maclure and Greenland 1992 (
link)).
To explore associations with CHDs within a multipollutant context, a principal component analysis (PCA) was conducted among participants who lived within 50 km of each type of monitor. PCA is used to reduce the number of correlated variables into a smaller number of artificial variables that capture most of the variance of the original variables while being uncorrelated with each other (Hatcher 1994 ). This allows the resulting factors to be included within the same model, reducing issues of multicollinearity. Applying PCA, we retained components that accounted for at least the same or more variance than one of the original pollutant variables. We then applied a varimax rotation and calculated factor scores for each participant. These factor scores were categorized using the 10th, 50th, and 90th centiles and used to assign exposure in hierarchical models.
Stingone J.A., Luben T.J., Daniels J.L., Fuentes M., Richardson D.B., Aylsworth A.S., Herring A.H., Anderka M., Botto L., Correa A., Gilboa S.M., Langlois P.H., Mosley B., Shaw G.M., Siffel C, & Olshan A.F. (2014). Maternal Exposure to Criteria Air Pollutants and Congenital Heart Defects in Offspring: Results from the National Birth Defects Prevention Study. Environmental Health Perspectives, 122(8), 863-872.
