Pediatric Psychiatric Emergency Department Utilization and Fine Particulate Matter: A Case-Crossover Study

Outdoor air pollution adversely affects respiratory, cardiovascular, perinatal, and other health outcomes. Globally, such pollution is estimated to be responsible for 3.3 million premature deaths annually (Dockery et al. 1993; Dockery 2009; Dominici et al. 2006; Pope and Dockery 2006; Zanobetti et al. 2009). Within the complex mixture of particles, gases, trace metals, and other toxic compounds that comprise outdoor air pollution, particulate matter with an aerodynamic diameter $<2.5\mathrm{\mu }\mathrm{m}$ (${\mathrm{PM}}_{2.5}$) is of particular health relevance due to its size, composition, and ability to induce inflammation and oxidative stress (Allen et al. 2017; Bell et al. 2007, 2010; Franklin et al. 2008; Kelly 2003; Kelly and Fussell 2012; Wei et al. 2009). Importantly, ${\mathrm{PM}}_{2.5}$ may affect the central nervous system (Block and Calderón-Garcidueñas 2009; Block et al. 2012; Calderón-Garcidueñas et al. 2002; Costa et al. 2014), possibly resulting in secondary neuroinflammation (Block and Calderón-Garcidueñas 2009; Block et al. 2012; Brockmeyer and D’Angiulli 2016; Sunyer 2008), which may, in turn, influence affective processing, anxiety, cognition, and behavior (Khandaker et al. 2015). Previous epidemiologic studies have reported associations between short-term (days to weeks) ${\mathrm{PM}}_{2.5}$ exposure and anxiety symptoms in older adults (Pun et al. 2017), suicide completion (Bakian et al. 2015), suicide attempts (Szyszkowicz et al. 2010), and emergency department (ED) visits for exacerbation of psychiatric disorders (Szyszkowicz et al. 2009, 2016; Oudin et al. 2018; Kim et al. 2018). However, null studies exist (Fernández-Niño et al. 2018) and the evidence on the effects of ${\mathrm{PM}}_{2.5}$ exposures on psychiatric disorders in adults is scarce.

Children may be particularly susceptible to neurotoxic effects of air pollution, but most epidemiologic studies to date have focused on the effects of long-term exposures and cognitive or behavioral outcomes. For instance, significant associations have been observed between increased prenatal and childhood exposure to air pollution and both lower intelligence quotient scores (Suglia et al. 2007; Porta et al. 2016; Perera et al. 2009) and increased hyperactivity during later childhood (Newman et al. 2013) as well as chronic exposure to air pollution at age 12 and diagnosis of major depressive disorder by age 18 (Roberts et al. 2019). However, few studies have examined the impact of air pollution exposure on depressive and anxiety disorders in youth. In fact, to our knowledge, there are no studies of short-term air pollution exposure and acute psychiatric outcomes in children and adolescents. Inflammation is involved in nearly all of the negative health effects of ${\mathrm{PM}}_{2.5}$ (Feng et al. 2016) and given that acute systemic inflammation is found in humans after exposure to ${\mathrm{PM}}_{2.5}$ (Zhao et al. 2013; Niu et al. 2013; Schneider et al. 2010), we hypothesized that acutely elevated ${\mathrm{PM}}_{2.5}$ concentrations could significantly contribute to exacerbations of psychiatric disorders, which have previously been linked to inflammation and microglia activation (Réus et al. 2015; Calcia et al. 2016). Our hypothesis parallels that of asthma, wherein an individual with an underlying inflammatory disorder is prone to acute exacerbations after experiencing acute increases in air pollution exposure. Indeed, high rates of comorbidity between asthma and psychiatric disorders have been noted and are both hypothesized to be mediated through inflammation of the central nervous system (Kewalramani et al. 2008). With this in mind, we sought to investigate the relationship between short-term exposure to ${\mathrm{PM}}_{2.5}$ in children and adolescents and the risk for pediatric psychiatric ED visits. Using a satellite-based model to estimate residential-level ${\mathrm{PM}}_{2.5}$ exposures and a case-crossover study design, we hypothesized that elevated ${\mathrm{PM}}_{2.5}$ exposure would be associated with an increased risk for pediatric psychiatric ED visits. Because neighborhood-level factors, including socioeconomic conditions, may modify the association between environmental exposures and health outcomes, we also hypothesized that community deprivation would modify the association between air pollution exposure and pediatric psychiatric ED visits.

Hamilton County has about 190,000 children living across 222 urban, suburban, and rural census tracts, and Cincinnati Children’s Hospital Medical Center (CCHMC) is a 629-bed pediatric academic health center that sees 99% of all hospital admissions and 81% of all ED visits among 0- to 14-y-olds and 78% of all hospital admissions and 47% of all ED visits among 15- to 18-y-olds within the county (Beck et al. 2018). Of all patients hospitalized for psychiatric care at CCHMC, 90% are admitted from the ED. ED visits to CCHMC for psychiatric care occurring between 1 January 2011 and 31 December 2015 were identified based on primary diagnosis of International Statistical Classification of Diseases and Related Health Problems, Tenth Revision (ICD-10) (WHO 2016) codes and having a ZIP code likely to be within Hamilton County. Primary diagnosis codes within the Mental Illness classification, as defined by the Clinical Classifications Software Healthcare Cost and Utilization Project (HCUP), were grouped to parallel Diagnostic and Statistical Manual of Mental Disorders (DSM)-5 (APA 2013) nosology and categorization (e.g., depressive disorders, anxiety disorders). A board-certified child and adolescent psychiatrist (JRS) familiar with the institutional practices of ICD coding/documentation reviewed the excluded or altered codes based on both a) current clinical practice and b) current child and adolescent psychiatry diagnostic formulation. Some ICD-10 codes were not assigned to any category (e.g., alcohol-related disorders, cognitive disorders) because we did not believe that they represented an acute exacerbation that could have been possibly caused by an acute chemical insult. Table S1 presents the 277 ICD-10 codes that were considered and the 15 categories that were generated. For each visit, the date of ED visit, residential address, parent-reported race (defined as white or Caucasian, black or African American, or other), insurance payor status (defined as public, private, or self-pay), age, sex, and ICD-10 code for primary diagnosis were extracted. Institutional addresses, for example, the address for CCHMC or Family and Child Welfare Services, and post office boxes were excluded and all residential locations of study participants at the time of each encounter were geocoded using our stand-alone and validated geocoder (Brokamp et al. 2018b). Of 15,489 total encounters, 13,812 (89.2%) were retained for analysis because they were within Hamilton County and had a residential address that could be geocoded with sufficient precision (range or street level and score $>0.5$) to assign ${\mathrm{PM}}_{2.5}$ concentrations and community deprivation. If a child had another visit for the same reason within 7 d, then only the initial visit was considered for analysis. Table S2 describes the inclusion criteria as well as the number of cases meeting each criterion. All data were obtained following review and approval by the CCHMC institutional review board. Informed consent was not required to use patient data because there was no contact with the patients or their families and the risk for loss of confidentiality was deemed to be minimal.

We conducted a time-stratified case-crossover study (Maclure 1991; Janes et al. 2005) to evaluate the potential association between ${\mathrm{PM}}_{2.5}$ exposure and psychiatric pediatric ED visits. A case-crossover study design is appropriate when both the exposure of interest and the outcome are transient, as in the case of elevated ${\mathrm{PM}}_{2.5}$ concentrations and ED visits. A key advantage to the design is that each case serves as its own control and is, therefore, matched for fixed individual characteristics including age, sex, and socioeconomic status. We chose controls by stratifying on a 45-d period and day of the week to control for temporal trends by design (Levy et al. 2001) and to ensure global exchangeability (Vines and Farrington 2001). Most time-stratified case-crossover studies stratify on month, year, and day of the week, but unlike most of these studies, our exposure time-series was not fixed, for example, each subject had a different estimation of ${\mathrm{PM}}_{2.5}$ exposures even if selected on the same day. Thus, we chose our stratification periods to maximize our statistical power (both the number of controls and the within-individual variability of ${\mathrm{PM}}_{2.5}$ exposure estimates) while minimizing the possibility for temporal confounding by long-term changes. Based on our time-stratified design, case windows could have between five and six control windows. However, because estimated exposures to ${\mathrm{PM}}_{2.5}$ were only available through 2015, some of the time periods were truncated, resulting in 26 (1.4%) cases having between two and four control windows. To assess the sensitivity of our results to our choice of using a 45-d period for temporal stratification, we also repeated our main analyses using a 30- and 60-d period for temporal stratification. The odds ratios (ORs) and 95% confidence intervals (CIs) were extremely similar among all three lengths of the temporal stratification (see Table S3), which suggests that our results were not affected by this choice.

Daily ambient concentrations of ${\mathrm{PM}}_{2.5}$ (as 24-h averages) were estimated at the residential location on the dates specific to both cases and controls using a previously developed and validated spatiotemporal model (Brokamp et al. 2018a). Briefly, our ${\mathrm{PM}}_{2.5}$ model is based on satellite-derived measurements of aerosol optical depth, a measure of the scattering of electromagnetic radiation due to aerosols in the atmosphere, calibrated with ground-based ${\mathrm{PM}}_{2.5}$ monitoring data. Additional spatiotemporal predictors in the model include meteorological measurements from the North American Regional Reanalysis data set: visibility (in meters), planetary boundary height (in meters), air temperature at $2\mathrm{m}$ (in Kelvin), relative humidity at $2\mathrm{m}$, precipitation rate (in kilograms per meter squared per second), accumulated total precipitation (in kilograms per meter squared), pressure (in pascals), wind speed and direction at $10\mathrm{m}$ (in meters per second) as well as land use (in percent developed imperviousness; National Land Cover Database), total length of major roadways, greenspace (mean normalized differential vegetation index), and indicators for grid, day of year, and year. After harmonization of the spatiotemporal data sets to a $1\times 1$ km grid, we trained machine learning models, specifically random forest models, to predict ${\mathrm{PM}}_{2.5}$ concentrations. Cross validation of the model using empirical data obtained from air monitoring stations indicated that the ${\mathrm{PM}}_{2.5}$ predictions were excellent across the study area with a cross-validated median absolute error of 0.95 ${\mathrm{\mu }\mathrm{g}/\mathrm{m}}^3$ and a cross-validated $R^2$ of 0.91. We also considered lagged effects of ${\mathrm{PM}}_{2.5}$ by estimating concentrations at residential locations 1, 2, and 3 d prior to the case and control dates.

Community deprivation was assigned using a previously validated census tract–based deprivation index (Brokamp et al. 2019) wherein a principal components analysis was used to create an index from 0 to 1, with higher scores indicative of greater deprivation. The index is based on six different census tract–level variables derived from the 2015 5-y U.S. Census’ American Community Survey: a) fraction of households with income below the poverty level, b) median household income, c) fraction of population 25 y and older with educational attainment of at least high school, d) fraction of population with no health insurance coverage, e) fraction of households receiving public assistance income or enrolled in the Supplemental Nutritional Assistance Program, and f) fraction of housing units that were vacant. The census tract–level community deprivation index was assigned to each child based on the census tract to which their residential address was geocoded to at the time of their encounter. For modeling, the deprivation index was dichotomized at the median of all 222 census tracts in Hamilton County (0.40), which was equivalent to the median deprivation index of all cases included in our analysis.

Conditional logistic regression models (Breslow and Day 1980) were used to estimate the OR of a psychiatric ED encounter from 0 to 3 d after a $10{\text{-}\mathrm{\mu }\mathrm{g}/\mathrm{m}}^3$ increase in ${\mathrm{PM}}_{2.5}$ exposure. Separate models were created for each psychiatric ED encounter type and lag days, defined as the time between the exposure and outcome. Models were estimated using a stratified Cox model with each encounter assigned to its own stratum (Gail et al. 1981). Covariates in each model included cubic splines of daily temperature and relative humidity, as well as two indicator variables: one for the presence of a holiday on the observation date and one for the presence of a holiday 1–2 d preceding the observation date. Temperature and humidity were obtained from the North American Regional Reanalysis data set (Mesinger et al. 2006).

To evaluate potential differences in associations between children living in high- and low-deprivation communities, we repeated models for each outcome and lag after adding an interaction term between ${\mathrm{PM}}_{2.5}$ and an indicator for high- or low-deprivation index, and tested for a significant improvement in model fit (chi-square $p<0.05$ for the difference in the log-likelihood). In accord with our prespecified model selection strategy, we estimated separate associations for high and low community deprivation only for psychiatric ED encounter types with a significant interaction for at least one exposure lag.

All statistical and geospatial computing was done in R (version 3.4.3) using the sf (Pebesma 2018) and survival packages.

During the 5-y study period, 13,176 psychiatric ED encounters meeting our inclusion criteria (see Table S2) occurred, representing 6,812 unique children. Of all children utilizing the ED for psychiatric reasons, 4,340 (63.7%) had only 1 visit, 1,219 (17.9%) had 2 visits, and 1,253 (18.4%) had at least 3 visits. The median [interquartile range (IQR)] age at encounter was 14.4 y (11.7, 16.1). Among all encounters, patients were 50% ($n=\mathrm{6,643}$) female, 44% ($n=\mathrm{5,756}$) African-American, and 2% ($n=264$) Hispanic; they were 66% ($n=\mathrm{8,740}$) publicly insured, 32% ($n=\mathrm{4,190}$) privately insured, and 2% ($n=245$) self-pay. The median ${\mathrm{PM}}_{2.5}$ (IQR) exposure estimated on the date of each case was $10.5{\mathrm{\mu }\mathrm{g}/\mathrm{m}}^3$ (8.1, 14.5). Based on the encounters, we generated 71,481 control windows, meaning that, on average, each case window had 5.4 control windows. Each psychiatric encounter was grouped into 1 of 12 categories and the number of encounters as well as demographic summary information by psychiatric grouping is presented in Table 1. The most frequent categories for psychiatric ED encounters included depressive disorders ($n=\mathrm{3,847}$), externalizing disorders ($n=\mathrm{1,850}$), other mood disorders ($n=\mathrm{1,903}$), and impulse control disorders ($n=\mathrm{1,755}$). Developmental disorders and personality disorders were excluded from further analysis because of the low number of ED visits (88 and 142, respectively).

ORs from conditional logistic regression models adjusted for temperature, humidity, and holidays indicated significant associations between a $10{\text{-}\mathrm{\mu }\mathrm{g}/\mathrm{m}}^3$ increase in ${\mathrm{PM}}_{2.5}$ and any psychiatric ED visit 1 d [$\text{OR}=$1.07 (95% CI: 1.02, 1.12)] and 2 d later [$\text{OR}=$1.05 (95% CI: 1.00, 1.10)] (Figure 1; see also Table S3). When grouped by psychiatric encounter type, we found that increased ${\mathrm{PM}}_{2.5}$ was significantly associated with ED visits related to schizophrenia on the same day [$\text{OR}=$1.25 (95% CI: 1.00, 1.57)], adjustment disorder 1 d [$\text{OR}=$1.24 (95% CI: 1.02, 1.52)] and 2 d later [$\text{OR}=$1.24 (95% CI: 1.02, 1.51)], other mood disorders 2 d later [$\text{OR}=$1.15 (95% CI: 1.02, 1.30)], and suicidality 1 d later [$\text{OR}=$1.44 (95% CI: 1.03, 2.02)].

p-Values for modification of associations between ${\mathrm{PM}}_{2.5}$ and ED visits by community deprivation (high vs. low) were significant ($p<0.05$) for at least one lag period for ED visits related to anxiety (lag 0,1), adjustment disorder (lag 0,1), and suicidality (lag 0) (see Table S4). Therefore, in accord with our a priori selection strategy, we report estimates stratified by community deprivation for these outcomes only (Figure 2; see also Table S5). Associations between a $10{\text{-}\mathrm{\mu }\mathrm{g}/\mathrm{m}}^3$ increase in ${\mathrm{PM}}_{2.5}$ and psychiatric ED visits were stronger for children living in high- versus low-deprivation communities for anxiety-related ED visits on the same day [OR 1.39 (95% CI: 0.97, 2.00) vs. 0.77 (95% CI: 0.57, 1.05)] and 1 d later [OR 1.39 (95% CI: 0.96, 2.01) vs. 0.85 (0.62, 1.17)], and for ED visits related to suicidality on the same day [OR 1.98 (95% CI: 1.22, 3.23) vs. 0.93 (0.60, 1.45)]. Conversely, associations with adjustment disorder-related ED visits were weaker for children living in high- versus low-deprivation areas, for visits on the same day [OR 0.82 (95% CI: 0.62, 1.08) vs. 1.25 (95% CI: 0.96, 1.62)] and 1 d later [OR 1.00 (95% CI: 0.76, 1.33) vs. 1.50 (95% CI: 1.16, 1.93)].

To our knowledge, this is the first report of associations between short-term ${\mathrm{PM}}_{2.5}$ exposures and psychiatric ED visits in children and adolescents. We estimated significant associations for all visits combined and for ED visits with primary diagnoses of suicidality and adjustment disorder. In addition, we found evidence that associations were modified by community deprivation, with stronger associations between ${\mathrm{PM}}_{2.5}$ and ED visits for suicidality and anxiety disorders for children in high-deprivation census tracts, and stronger associations with ED visits for adjustment disorder for children in low-deprivation census tracts. Notably, all daily exposures within our study domain were below the National Ambient Air Quality Standards set by the U.S. Environmental Protection Agency (EPA; i.e., $35{\mathrm{\mu }\mathrm{g}/\mathrm{m}}^3$ for ${\mathrm{PM}}_{2.5}$ measured over 24 h).

The neuroendocrine sympathetic–adrenal–medullary and hypothalamic–pituitary–adrenal axes mediate a wide array of health effects associated with both chemical and nonchemical stressors, and it has been suggested that more study on the interactive influence of air pollution, lifestyle, and environmental factors on physiologic stress responses is required (Snow et al. 2018). Animal studies have demonstrated that psychosocial stressors such as chronic restraint, social isolation, and repeated social defeat result in increased microglial activation (Schmeer and Yoon 2016), and it is known that childhood poverty contributes to stress and inflammation and harms normal brain and immune system development (Wadsworth and Rienks 2012; McEwen 2007; Shonkoff et al. 2009, 2012; Evans and Kim 2012). Moreover, epidemiologic studies suggest that children exposed to chronic stress may be more susceptible to air pollution-related respiratory deficits than other children (Clougherty et al. 2007; Islam et al. 2011).

The possibility that material community deprivation strengthens the relationship between ${\mathrm{PM}}_{2.5}$ exposure and ED presentation for anxiety disorders and suicidality may have clinical implications. Cross-sectional and longitudinal studies have suggested that the risk of developing an anxiety disorder is increased among those who experience factors associated with deprivation, including low self-esteem, childhood sexual abuse, educational attainment, number of traumatic experiences by 21 years of age, childhood separation, and a “disturbed family environment” (Beesdo et al. 2010; McCulloch 2001; Meltzer et al. 2007; Kingsbury et al. 2015). Similarly, in cross-sectional epidemiologic studies of youth, the risk of suicidality relates to several factors that relate to deprivation (e.g., poor family environment, low parental monitoring) (King et al. 2001) and is greater in urban settings (Husky et al. 2012). Although these adverse childhood experiences may not work exclusively through stress pathways, they are frequently co-occurring with deprivation factors and are considered one of multiple dimensions of early adversity in general (McLaughlin et al. 2014). Our results raise the possibility that primary or secondary prevention strategies directed either at reducing ${\mathrm{PM}}_{2.5}$ exposure or decreasing the vulnerability to ${\mathrm{PM}}_{2.5}$ exposure, for example, interventions that attenuate the impact of deprivation in youth, could abort the development or forestall the exacerbation of this psychopathology.

Two advantages of our study were the high spatiotemporal resolution of ${\mathrm{PM}}_{2.5}$-exposure assessment and the large availability of the vast majority of psychiatric hospitalizations within Hamilton County, which allowed us to conduct the individual-level case-crossover analysis at a near population level. This increased precision allowed us to examine susceptibility characteristics using individual-level data, an advantage that is not feasible in time-series and ecologic studies. The case-crossover study design prevents confounding by covariates that do not differ within 45 d, such as race and individual-level socioeconomic status; however, like all observational studies, unobserved confounding is possible, and temporal confounders that we did not consider could have biased our results.

Although this is the first examination of the impact of ${\mathrm{PM}}_{2.5}$ exposure and community deprivation on psychiatric outcomes in children and adolescents, there are important limitations that warrant additional discussion. These data were abstracted from medical records and the diagnoses were not established with the use of structured, clinician-administered assessments that apply strict diagnostic criteria. For example, suicidality reflects a heterogeneous clinical phenomenon and includes suicidal ideation, preparatory suicidal behavior and planning, aborted suicide attempts, and so on. The ability to resolve these aspects of suicidality in the present sample is limited; however, regardless of its severity, suicidality represents a phenomenon that necessitates clinical evaluation and, generally, intervention. Furthermore, we examined psychiatric disorder-related ED visits, ostensibly reflecting either a severe symptom burden or severe exacerbation; however, it is possible that exacerbations of other disorders that did not result in an ED visit may have been present. In addition, outcomes were classified based on the primary diagnosis entered into the medical record, without considering co-occurring disorders or secondary indications for the ED visit. In general, the use of data extracted from the electronic health record and a crossover study design allowed us to study a rare outcome at a population level, which might otherwise have been possible only by assembling a prospective cohort or ascertaining cases and matching controls.

In conclusion, we found evidence of an association between short-term exposure to ambient ${\mathrm{PM}}_{2.5}$ and the exacerbation of psychiatric disorders in children as indicated by increases in psychiatric ED visits. These findings need to be confirmed in other populations, but they strongly support the need for further research on the potential influence of ambient air pollution on child and adolescent mental health.

The authors did not receive any grants that supported this work.

The authors declare they have no actual or potential competing interests.