1 

# Children’s Health Project: Linking Asthma to PM 10 

# in Central Phoenix – a report to the Arizona Department of Environmental Quality

H.J.S. Fernando 1 

### , R. Dimitrova 1


### , G. Runger 2


### , N. Lurponglukana 2


### , P. Hyde 1 

### , B. Hedquist 1 

### , J. Anderson 1



### Arizona State University’s 1 

### Center for Environmental Fluid Dynamics and
2 

### Center for Health Information and Research

February 16, 20092 TABLE OF CONTENTS

List of Figures..............................................................................3

List of Tables...............................................................................6

Executive Summary........................................................................7

Chapter 1: Introduction....................................................................7

Chapter 2: PM Monitoring Networks....................................................8

Chapter 3: Interpolation Methods......................................................13

Chapter 4: Neutral Network.............................................................17

Chapter 5: Air Quality Modeling.......................................................19

Chapter 6: Neutral Ne tworks vs. CMAQ..............................................37

Chapter 7: Linking PM 10 Concentrations.............................................40

Chapter 8: Asthma warning systems....................................................54

Chapter 9: Conclusions..................................................................55

References.................................................................................56

Appendices................................................................................59

A. PM monitoring sites.............................................................59

B. Neural Network..................................................................65

C. MM5 - modeling results and discussion......................................67

D. Definitions of statistics..........................................................73

E. Mathematics of the case-crossover method..................................74

Glossary.....................................................................................773 LIST OF FIGURES

Figure 2.1 PM 10 monitoring sites with elevations for the Children’s Health Project 10

Figure 2.2 Spatial distribution of annual PM 10 concentrations –central Phoenix for 2006 12

Figure 3.1 Comparison between IDW and Ordina ry Kriging methods (Note: SS1 = SS = Phoenix Supersite) 15

Figure 3.2 Percentage differences between IDW and

Ordinary Kriging methods – December 19, 8am 15

Figure 3.3 Mean Absolute Errors (MAE) for December between IDW and Ordinary Kriging methods – based on daily maps 16

Figure 3.4 Root Mean Square Error (RMSE) for December between IDW and Ordinary Kriging methods – based on daily maps 16

Figure 4.1 24-hour average PM 10 predictions from the ne ural network versus observations at the Central Phoenix site for 2006 18

Figure 5.1 Central Phoenix study domain with monitoring locations 20
Figure 5.2 24-hour averaged PM 10 observations for 2005-2006 21
Figure 5.3 Daily emergency department visits and admissions to hospitals of children ages 0-17 in Maricopa Coun ty with an asthma diagnosis for 2005 – 2006 22

Figure 5.4 PM 10 concentrations and asthma events (2005-2006) 23
Figure 5.5 Predicted and observed hourly PM 10 concentrations in central Phoenix for periods of high pollution (left) and low po llution (right) – from all five monitors 24

Figure 5.6 Predicted and observed hourly PM 10 concentrations [μg/m 3 ] in central Phoenix for the two sets of monitors shown separately – high pollution period (November 2005) 24

Figure 5.7 The "soccer goal" plot measures both bias and error for central Phoenix 27

Figure 5.8 Daily Mean Bias for high polluti on conditions (top) and low pollution 28
Figure 5.9 CMAQ predictions versus observations of PM 10 concentrations at five sites for the high pollution period of November 2005. The lines are the median values of the 24 hour ly predictions for each day; the box plots are 25 th to 75 th percentiles of the 24 predictions 29 100 * 10 10 10

IDW KR IDW PM PM PM

−4 Figure 5.10 CMAQ predictions versus observations of PM 10 concentrations at five sites for the low pollution pe riod of March 11 – April 9, 2006. The lines are the median values of the 24 hourly predictions for each day; the box plots are 25 th to 75 th percentiles of the 24 predictions 30

Figure 5.11 Hourly variation of the five-sit e averages, calculated and observed PM 10 concentrations for November 2005 31

Figure 5.12 Contributions of different processes to PM 10 concentrations for 3- days periods in November 2005 and April 2006 32

Figure 5.13 The spatial distribution of the Index of Agreement between CMAQ and Inverse Distance Weighting (IDW) interpolated surfaces 33

Figure 5.14 The spatial distribution of Mean Bias between CMAQ and Inverse Distance Weighting (IDW) interpolated surfaces 33

Figure 5.15 Comparison of the daily, “modeled” PM 10 concentrations with the observations at Central Phoenix: vi olet is straight CMAQ; blue is “CMAQ-new” -- CMAQ combined with interpolation 34

Figure 5.16 Comparison of the daily, “modeled” PM 10 concentrations with the observations at West Phoenix (WP): vi olet is straight CMAQ; blue is “CMAQ-new” --CMAQ combined with interpolation 34

Figure 5.17 Comparison of the daily, “modeled” PM 10 concentrations with the observations at Supersite (SS): viol et is straight CMAQ; blue is “CMAQ-new” -- CMAQ combined with interpolation 34

Figure 5.18 Comparison of the daily, “modeled” PM 10 concentrations with the observations at West 43 rd Avenue (WF): violet is straight CMAQ; blue is “CMAQ-new” -- CMAQ co mbined with interpolation 35

Figure 5.19 Comparison of the daily, “modeled” PM 10 concentrations with the observations at Durango Complex (D C): violet is straight CMAQ; blue is “CMAQ-new” -- CMAQ co mbined with interpolation 35

Figure 6.1 Time series of 24-hour average PM 10 concentrations predicted by CMAQ and EnviNNet with observations at the Central Phoenix site 37

Figure 6.2 Scatter plots of hourly PM 10 concentrations predicted by CMAQ (left) and EnviNNet w ith observations [μg/m 3 ]: November 2005, Central Phoenix site 38

Figure 6.3 Comparison of hourly PM 10 concentrations predicted by Models 3 and EnviNNet with the observati ons of November 5-9, 2005, Central Phoenix site 39

Figure 6.4 Comparison of hourly PM 10 concentrations predicted by CMAQ and EnviNNet with the observations of November 25-29, 2005, Central Phoenix site 395 Figure 7.1 Washout periods for asthma incidences within two and five miles of a continuous PM10 monitor 41

Figure 7.2 Asthma incidences by different characteristic within two miles of PM 10 monitoring sites 42

Figure 7.3 Asthma incidences by different characteristic within five miles of PM 10 monitoring sites 43

Figure 7.4 Asthma incidences by zip code within two and five miles of PM 10 monitoring sites 44

Figure 7.5 Comparison of AZHQ and ADHS asthma incidences within two and five miles of PM 10 monitoring sites 45- 46
Figure 7.6 Correlation matrix among Dailymean PM 10 and lag variables 50
Figure 7.7 Adjusted odd ratios with 95% lower and upper bound confidence limits for the four asthma populations 53

Figure 7.8 Natural logarithm of the adjusted odds ratio (Log Odds) of asthma incidents versus quartiles of mean PM 10 concentrations 54

Figure A1.1 West 43 rd Avenue PM 10 monitoring site 61
Figure A1.2 Greenwood PM 10 monitoring site 61
Figure A1.3 Durango Complex PM 10 monitoring site 62
Figure A1.4 Central Phoenix PM 10 monitoring site 62
Figure A1.5 West Phoenix PM 10 monitoring site 63
Figure A1.6 Supersite PM 10 monitoring site 63
Figure A1.7 South Phoenix monitoring site 64
Figure C-1 Modeling domains used with weather predicting model MM5 66
Figure C-2 Scatter plots of wind speed and temperature for all stations for high pollution episode November 2005 68

Figure C-3 Scatter plots of wind speed and temperature for all stations for the low pollution period of March 11 – April 9, 2006 68

Figure C-4 Comparison between predicted a nd observed wind speeds and directions at different locations – with strong synoptic winds 69

Figure C-5 Comparison between predicted a nd observed wind speeds and directions at different locations – with calm synoptic conditions 70

Figure C-6 Simulated horizontal winds at 10m and surface temperature fields during nocturnal (left) and daytime conditions, November 24, 2005 716 LIST OF TABLES

Table 2.1 Central Phoenix continuous PM 10 monitoring sites 9 Table 2.2 Annual average (2006) PM 10 concentrations 11
Table 2.3 December 2007 PM 10 monthly average concentrations 13
Table 5.1 Statistical measures in comparing model-predicted hourly PM 10 concentrations to observations for bo th study periods at different sites 27

Table 6.1 Statistics of CMAQ and EnviNNet predictions of hourly PM 10 concentrations at the Central Phoenix site 38

Table 7.1 Referent selection based on tim e-stratified scheme within a 28-day
stratum 47

Table 7.2 Model statistics fo r each main effect (pre school age group omitted) 50
Table 7.3 Results for all asthma incidents 51 Table 7.4 Results for asthma incidents excluding the preschool group 51 Table 7.5 Result for ambulatory events only 52 Table 7.6 Result for ambulatory events without the preschool group 52
Table 7.7 Case-crossover statistics relating asthma incidents to PM 10 concentrations in central Phoenix 53

Table A-1 Characteristics and images of continuous PM 10 monitoring sites 59 Table A-2 Traffic near PM monitoring sites 60
Table C-1 Grid definitions for MM5 model 65 Table C-2 MM5 configuration for different domains 65 Table C-3 RMSE for wind sp eed at different sites 71 Table C-4 Summary of statistical measures , in comparing obse rvations to model predictions of wind speed at different sites 727 Executive summary

The effects of elevated PM 10 concentrations on childhood asthma in cidents (emergency room visits & hospital admissions) in central Phoenix, AZ we re determined. A major portion of this work involved the application of a deterministic, numerical air quality model to produce gridded concentration fields, the development of a ne ural-network based pr ediction system for PM 10 , and the interpolation of PM 10 concentrations from five monitoring sites to produce censu s-tract specific concentrations for the statistic al health effects analyses. Through the case-crossover method, PM 10 concentrations were shown to have statistically significant associations with asthma incidents in central Phoenix at the 95% confidence level. For child ren ages 5-17 an effect much stronger than in previous studies was determined, while children ag es 0-4 exhibited virtually no effects. Previous studies have reported the risk of adverse health effects as a function of the cha nge of daily mean PM 10 from the 25th to 75th percentile; in this study that change was 36 μ g/m 3 . This change in PM 10 is associated with a 13% increase in the probabi lity of asthma incidents among children ages 5-17. Within this group, neither age nor gender was signi ficant. The lack of an effect on preschool children was assumed to be due to the difficult y of making accurate asthma diagnoses in the younger ages and/or more time spent indoors.

Chapter 1 Introduction
Voluminous studies on the health effects of pa rticulate matter (PM) ha ve been carried out. Atmospheric scientists deal w ith the conditions conducive to high PM episodes; the exposure community delves into the duration and severity of ambient PM concentrations to which various populations are subject; toxicologists and bi ochemists study the effects of the hazardous constituents of PM on human and surrogate organisms and in in vitro systems; physiologists and physicians conduct exposure-effect studies on humans in chambers, epidemiologists relate health outcomes to ambient concentrations and exposures, and modelers (including risk assessors) utilize various data to estimate exposure-response relations hips for determining safe levels of exposure. Integration of such disparate disciplines is difficu lt, given the complexity inherent in each of the approaches, though the U.S. Environmental Prot ection Agency (EPA) and the World Health Organization (WHO) Air Quality Criteria Documents, the EPA Staff Papers, and the EPA Office of Air do make such an attempt. For issuing health wa rning system is sufficient to establish, if there be any, a reliable association betw een PM concentrations and health effects such as admittances due to respiratory (or cardiovascular) ailments to allo w prediction of such even ts. The principal goal of the “Children’s Health Project” was to clarif y the association of children’s health with PM 10 levels to develop a better health warn ing system in Phoenix. We brought together a partnership among Arizona Department of Environmental Quality ( ADEQ), Arizona Department of Health Services (ADHS), Ross & Associates, and two groups of resear chers at Arizona State University: the Center for Health Information and Research (CHIR), and the Center for Environmental Fluid Dynamics (EFD). We discussed our findings with Univers ity of Arizona medical personnel and the Asthma Coalition. Work elements of the project

• With the central portion of P hoenix as the study area, PM 10 concentrations and asthma records were analyzed for each day in 2005 – 2006.

• PM 10 concentrations were estimated with bot h deterministic and st atistical models.8 • The deterministic models were the EPA’s Community Multi-scale Air Quality Modeling System (CMAQ) and a meteorological model cal led “MM5.” This grid-based weather and pollution simulation system relied on the most recent Western Regional Air Partnership PM inventories. Design days incl uded periods of low, intermed iate, and high pollution. The hourly and daily simulated concentrations were compared with data from five permanent monitoring sites with continuous PM 10 instruments.

• The first type of statistical model employed was interpolation . With the interpolation technique called “inverse distance weighting” (IDW), spatial concentration fields were built that provided a PM 10 concentration for each of the 168 ce nsus tracts from measurements at the five central-city monitoring s ites. The efficacy of this sc heme was verified by deploying four temporary monitors for f our months and by using hourly PM 10 concentrations from two Maricopa County permanent sites whose continuous PM 10 monitors began operating after 2005. These six additional sites enabled the use of a more sophisticated interpolation technique called “ordinary krigin g”, which produced results sim ilar to those of the IDW.

• The second type of statistical model was the neural network , developed by Italian scientists and applied to a single permanent monitoring si te -- Maricopa County’ s “Central Phoenix” site. This model is a statistical system to predict the next day’s PM 10 concentrations without the labor and emissions-inventory intensive a lternative of CMAQ. In the present work, however, its use was limited to estimating PM 10 concentrations already recorded in 2006.

• About 6,000 ADHS asthma events in children ages 0-18 – in a unified area within five miles of the five monitoring sites – comprise d 75% of the Maricopa County total. These ADHS records of emergency room treatment and hospital admissions were cross-checked against an independent data se t at CHIR and were found to be consistent once duplicate records had been eliminated.

• The census-tract-based health data were linked with the PM 10 concentrations derived from the IDW interpolations by a statistical method called “ case-crossover ”. Asthma has a strong seasonal component and the case-crossover analys is associates a reference time period as a control for each case. Consequently, each patien t was treated as a matched case-referent pair with control exposures obtained from the sa me patient at different times. With an appropriate referent period, th e case-crossover analysis cont rolled for long-term trends, seasonal effects, epidemic, and other covariat es that change slowly with time such as behavior and diet. In addition to PM 10 on the “event” day, 24-hour averages from one to seven days before the event were also examined (this delay is referred to as “lag time”). The patients were also stratified by age group, ge nder, and race/ethnicity, and interactions between air quality and these variables were explored.

Chapter 2 PM monitoring networks

Permanent PM 10 network in the central city Considering “metropolitan Phoenix” to include th e city and all of the adjoining suburbs (but excluding Casa Grande and environs ), five governmental agencies conduct air quality monitoring:

• Maricopa County Department of Air Quality (MCDAQ)

• Arizona Department of Environmental Quality – Air Qualit y Division (ADEQ)

• Pinal County Air Quality Control District9 • Salt River Pima Maricopa Indian Community, and

• Gila River Indian Community.

Together, these agencies operate about 50 different sites with a variety of instruments and several different objectives. What the U.S EPA terms “criteria pollutants”, such as particulate matter, ozone, carbon monoxide, and others, may be the prime focus; but ma ny other aspects of pollution in the planetary boundary layer are also monitore d. For example, this monitoring includes meteorological parameters such as wind speed and wi nd direction, visibility characteristics such as light extinction and photographic views, and the ch emical composition of particulates and gaseous pollutant mixtures (MCDAQ, 2007; ADEQ, 2007).

The present work, based on the PM 10 network and the asthma cases in central Phoenix, relies on the “continuous” meteorological and air pollutant data from the first two agencies in the list above. “Continuous” means that the instru ments run round the clock with their meteorological parameters and air pollutant con centrations, including PM 10 , reported as hourly averages . As described later in this chapter, Arizona State University (ASU ) and ADEQ staff set up and operated temporary PM 10 continuous monitors to supplement the perm anent network (Table 2.1 and Figure 2-1).

Table 2.1 Central Phoenix continuous PM 10 monitoring sites Symbol Name Major streets/city 1 Duration Agency WF West 43 rd Avenue 43 rd Avenue & Broadway Rd. DC Durango Complex 27 th Avenue & Durango Street WP West Phoenix 39 th Avenue & Thomas Road CP Central Phoenix 16 th Street & Roosevelt Street SS 2 Supersite 15 th Avenue & Camelback Rd. SP 3 South Phoenix Central Ave & Broadway Rd. GR 3 Greenwood 27 th Avenue & I-10 Permanent MCDAQ or ADEQ

MRV Maryvale 59 th Ave & Indian School Rd. VGC Valley Garden Center 15 th Avenue & McDowell Rd.

WVR Weaver’s Auto Service 29 th Street & Thomas Road

CSA Community Service of AZ 59 th Ave. & Glendale Avenue, Glendale, AZ Temporary: Nov. 2007 – March 2008 ASU/ADEQ

1. All except CSA are in Phoenix, Arizona. 2. SS, Supersite, is ADEQ; other permanent sites are MCDAQ. 3. Used only in the ordinary kriging interpolation, along with the other permanent and all temporary sites; not used in the inverse distance weighting interpolation for 2005 – 2006, because continuous PM 10 monitoring did not begin until January 2006.
Each of these sites has its unique combination of local emission sources and its own susceptibility to urban transport that togeth er, along with the desert backgr ound, result in the measured concentrations. Each has its own degree of represen tativeness: i.e. how fa r away from the monitor do the measured concentrations ac tually prevail? The permanent county sites, fully described in Maricopa 2007, and the ADEQ Supersite, are also ch aracterized in Appendix A. What follows is a discussion of how their average PM 10 concentrations fit into the spatial distribution of PM 10 concentrations throughout metropolitan Phoenix.10 SP GR WVR VGS CSA

MRV

SP GR WVR VGS CSA

MRV

Figure 2.1. PM 10 monitoring sites with elevations for the Children’s Health Project. The permanent monitors are black, the temporary – blue.

The PM 10 variation among all of the metropolitan sites can be seen in their 2006 annual averages, given in Table 2.2 and also expressed as a percenta ge of the lowest urban concentration site, Mesa. The central city sites, from WF to CP, have con centrations considerably higher than Mesa’s (from 38 to 162% higher), varying widely among themselves. The PM 10 gradients are steepest here, tapering off in the more outlying po rtions of the cen tral city.

Shown in Figure 2.2, these concentrations have thei r highest values and steepest gradients in southcentral Phoenix, with almost all of the highest on es being south of I-10. As one moves west of 75 th Ave, north of Bethany Home Rd, east of 36 th St., and south of Dobbi ns Rd, the prevailing concentrations generally fall to the mid to low 30’ s, or even lower in the case of the Phoenix Mountain Preserve of the South Mountains. Because measured PM 10 concentrations exhibit high sensitivities to localized emissions ; however, elevated concentrations can occur in the vicinity of such sources as major earthmoving projects or intensively tilled agricultural lands. Two such examples of these PM 10 discontinuities in metropolitan Phoe nix are Higley in the southeast and Buckeye in the far west. The former’s concentrati ons are elevated above the values expected with its distance from the city cen ter by emissions from major earthmoving activities from road construction, commercial building, and residential housing developments, as well as from continued farming. The latter’s con centrations are elevated from the same source categories, with the agricultural activities being the most influential. Th e greatest uncertainty in this concentration field lies to the west, but even this has been dimini shed by including measurements from a discontinued11 county site called “Maryvale”, at 61 st Avenue and Encanto Boulevard, normalized to those of Central Phoenix (the former is 1.7% higher than the latter).

Table 2.2. Annual average (2006) PM 10 concentrations (μg/m 3 ) (bold sites used in the present work)

Symbol Site name Location or type [PM 10 ] % from Mesa WF West 43 rd Avenue Central 79.9 162 DC Durango Complex Central 69.2 127 HI Higley Suburban 60.6 99 SP South Phoenix Central 55.0 80 BU Buckeye Urban fringe 53.2 74 GR Greenwood Central 51.7 70 WP West Phoenix Central 49.8 63 MV 1 Maryvale Central 42.7 40 CP Central Phoenix Central 42.0 38 GL Glendale Suburban 36.3 15 SSC South Scottsdale Suburban 35.0 19 NP North Phoenix Suburban 34.4 13 SS 2 Supersite Central 32.9 8 DY Dysart Urban fringe 32.3 6 ME Mesa Suburban 30.5 0 ES 3 Estrella Park (1995-2004) Urban fringe 30.2 -1 AJ 4 Apache Junction Urban fringe 23.6 -23 PV 5 Palo Verde (1997-2004) Urban fringe 21.9 -28 AO 6 Ajo (1995 – 2003) Remote small town 20.2 -34 OP 7 Organ Pipe (1995-2004) Background 10.6 -65 1. MV, Maryvale, a discontinued county site at 61 st Avenue and Encanto Blvd. The 2006 annual average is estimated by the MV/CP long-term average. 2. SS, Supersite, 2006 average was estim ated by the SS/CP long-term average. 3. ES, Estrella Park, a discontinued ADEQ site, near Estr ella Parkway and the Gila River, 18 mi (29 km) westsouthwest of downtown (Van Buren, Grand, 7 th Avenue intersection) 4. AJ, Apache Junction, a Pinal County site, 32 miles (51 km) east of downtown 5. PV, Palo Verde, a discontinued ADEQ site, 43 mi (69 km) west of downtown 6. AO, Ajo, Arizona, an ADEQ site, 83 mi (133 km) south-southwest (SSW) of downtown 7. OP, Organ Pipe National Monument, a National Park Service site, 105 mi (169 km) SSW of downtown

Although the two highest con centration sites, West 43 rd Avenue (WF) and the Durango Complex (DC) garner most of the regulatory attention, the Greenwood site deserves sp ecial mention, as well. Unlike the top two, whose elevated concentrations can in part be attributed to the somewhat localized windblown and m echanically generated PM 10 emissions along the Salt River, Greenwood, named for the cemetery across the street, is situat ed in a residential neig hborhood with little in the way of fugitive dust emissions. Instead, this mon itor is affected by heavy traffic on nearby paved roads. (Traffic on paved roads generates PM 10 emissions by resuspending soil particles from the road surface and by exhaust emissions, with the fo rmer having ten times the weight of the latter.) From the monitor inlet it is 10 mete rs east to the bus y arterial of 27 th Avenue (19,000 vehicles on an average weekday (ADT)), 905 m east to I-17 (1 12,000), and 85 m north to Interstate 10 (229,000 ADT). PM 10 concentrations from this monitor, with suitable adjustments for wind direction and traffic volume, could be assigned to all of the freeway corridors in the central city.12 xxxxxxxv

Figure 2.2. Spatial distribution of annual PM 10 concentrations –central Phoenix for 2006

Temporary PM 10 sites for this project Given the five permanent monitoring sites with complete records in 2005 – 2006, and considering that the Greenwood and South Phoenix continuous PM 10 monitors began operation after 2005, additional continuous PM 10 monitors were thought desirable fo r this project. ASU and ADEQ staff established four temporary sites and operate d them from November 2007 through March 2008. These additional concentrations, with the Greenw ood and South Phoenix sites included as well, provided the 11 points requi red for the more sophisticated inte rpolation technique called “ordinary kriging” (the simpler technique is inverse distan ce weighting, or IDW). Th e temporary sites were chosen to fill the gaps in the ge ographic coverage of the permanent PM 10 network and to be placed in areas with moderate to high numbers of asthma cases. These four tem porary sites are shown as three-initial labels in Table 2.1 and Figure 2.1, with additional details given in Appendix A. Even though the PM 10 concentration fields built with the 11 monitors did not differ much from those based on the original five, the efforts involved were considered well worth while. Without these six additional sites, the fidelity of the interpolated fields based on only the five permanent monitors would have remained too uncertain for comfort.

The PM 10 concentrations at these four temporary si tes were generally lower than those at the permanent sites (Table 2-3).13 Table 2.3. December 2007 PM 10 monthly average concentrations (μg/m 3 )

Chapter 3 Interpolation methods

Introduction Interpolating air pollutant measurem ents to produce spatial concentratio n fields has a rich history on all scales -- urban, regional, and global. What dis tinguishes the present work is the nature of the pollutant. Instead of gaseous air pollutants or fine particles that mimic gases (particles 2.5 microns and smaller, “PM 2.5 ”), this project concerns particle s 10 microns and smaller, or “PM 10 ”. In metropolitan Phoenix the coarse fraction of the aer osol (2.5 – 10 microns) comprises two thirds of the mass of PM 10 and consists of pa rticles that are 90% geological. Even the fine fraction has an appreciable geological component (20%), in a ddition to the usual secondary and combustion particles. Speciated aerosol measur ements show that somewhat less than one fourth of the regulated PM 10 comes from the better understood and more homoge neously distributed emission sources such as fuel combustion from motor ve hicles, industrial combustion, a nd food cooking, as well as the secondary components of nitrate and sulfate. Most of the remaining three quarters, virtually all geological, comes from paved road dust, constr uction dust, windblown dust, and unpaved road dust. It is precisely this geological, or “soil”, or “cru stal” component that defies efforts to build precise emission inventories, and, with a few notable exceptions, includi ng the present work, makes it difficult to simulate with numerical models on an urban scale. Because the continuous PM 2.5 monitoring network in 2005-2006 wa s limited to two stations, PM 10 was the choice for this project.

Practical constraints make it impossi ble to obtain data at every desi red point. Thus, interpolation is important and fundamental to gr aphing, analyzing and understanding the spatial distribution of air pollution throughout a metropolitan area. Interpolation methods selected in this study were inverse distance weighting and ordina ry kriging, the former with the five permanent PM 10 sites in 2005- 2006; the latter, with the expanded networ k from December 2007 – February 2008. Hourly concentrations were converted to daily 24-hour averag es for the daily spatial concentration fields. In addition, select hours were also ex amined on a case by cas e basis. Some consideration was given to using emission or traffic estimates to adjust the interpolated values, but this was rejected. Inverse distance weighting Inverse distance weighting (IDW) is one of the simplest interpolation methods. A neighborhood in the vicinity of the interpolated point is identified and a weighted average is taken of the observation Label Name Type [PM10] DC Durango Complex 49.6 WF West 43 rd Avenue 46.1 GR Greenwood 42.9 SP South Phoenix 39.6 WP West Phoenix Permanent

38.2 CSA CSA of Glendale Temporary 33.7 CP Central Phoenix Permanent 32.8 VGC Valley Garden Center 31.4 MRV Maryvale Temporary 28.5 SS Supersite Permanent 28.3 WVR Weaver’s Auto Service Temporary 21.814 values within it. The weights are a decreasing func tion of distance. The user has control over the mathematical form of the weighting function, the size of the neighborhood (e xpressed as a radius or a number of points), in addition to other options. In this project, IDW was used to interpolate the PM 10 concentrations of 2005-2006 from the five permanent monitoring sites with an automated function in Geographic Information Systems (GIS ). Within ArcGIS 9.2 software, the Spatial Analyst extension tool was used to create IDW interpolations with in the study area. Within Spatial Analyst, a pull-down menu is found in which the user can input the specific information pertaining to the dataset utilized. In the pu ll-down menu, default input parameters were selected. This included selecting PM 10 for the z-value, a power value of two, a va riable search radius with 12 points, andthe search radius distance value unspecified. Once the proper data files were selected in the pulldown menu, the GIS software then automatically cr eated a temporary inter polation raster, based on the PM 10 value at each of the five mon itor sites on the sp ecific day or hour selected. (A “raster” is a scanning pattern of parallel lines that form the display of an imag e projected on a display screen.) Next, the temporary interpolation ra ster files were made permanent by saving each hourly file in the table of contents section of ArcM ap, the principal graphical user in terface to create maps in ArcGIS. Finally, the interpolated surface value tables were a ggregated by census tracts in Excel. Census tract shapefiles, containing the new information from these interpolations, were then finalized as graphics in ArcMap. Ordinary kriging Ordinary kriging is a geostatistical approach to modeling. Instead of wei ghting nearby data points by some power of their inverted di stance, it relies on the spatial co rrelation structure of the data to determine the weighting values. This is a more ri gorous approach than IDW, as correlation between data points determines the estimated value at an unsampled point. Ordinary kriging makes the assumption of normality among the data points, as well as an unknown constant trend. This is unlike simple or universal kriging, which a ssume a known constant trend and/or model.

The word “kriging”is synonymous w ith “optimal prediction”. Kr iging methods, including ordinary kriging, utilize a variogram model to express the spatial variation of data points, and then minimize their error.

For this project, the Spatial Analyst extension in ArcGIS (the same used for IDW analysis) was chosen as the optimal automated met hod to create kriged surfaces of PM 10 . After selecting an hourly file containing PM 10 for the z-value, the spherical semiva riogram model was selected in the pull-down menu as the best way to create a kriged surface with the least amount of error. The advanced parameters tab, which allows for input on the lag size, major range , partial sill, and nugget size, was not used. Default options, including a va riable search radius w ith 12 neighbors, and an unspecified search distance, were used to create tem porary interpolated raster files. As was the case with IDW, interpolated files were first made perman ent, and were then aggregated into census tract shapefiles in Excel and ArcMap, to create maps of kriging values within study area census tracts. Comparison between IDW and ordinary kriging Each interpolation method has its advantages and disadvantages. Although ordi nary kriging is more rigorous and better at predicting su rfaces at a distance from known da ta points, IDW can be utilized with fewer data points. For example, ordinary krig ing performs best with a large number of data points, ideally above 50, although as few as 10 can be used. IDW on the other hand, can create a prediction surface, though based mainly on nearest ne ighbors, with as few as 3-4 points within a given area. Figure 3.1 illustrates resu lts from the two methods. Five PM 10 station values were used in IDW on the left; 11 stations were used in ordinary kriging on the right. Figure 3.2 shows the15 percentage differences betw een the two methods (PM10 IDW -PM10 KR )/PM10 IDW . In Figures 3.3 and 3.4, statistical comparisons between the two interpolation methods re veal better agreement close to the central part of the st udy area, in close proximity to the five permanent PM 10 monitors. Greater disagreement between the two models occurs in the outer edges of the study area, in census tracts close to the temporary PM 10 monitors.

Figure 3.1. Comparison between IDW and Ordinary Kriging methods

Figure 3.2. Percentage differen ces between IDW and Ordinary Kriging methods – December 19, 2007, 8am 100 * 10 PM 10 PM 10 PM

IDW KR IDW −16 Figure 3.3. Mean Absolute Errors (MAE) for December between IDW and Ordinary Kriging methods – based on daily maps

Figure 3.4. Root Mean Square Error (RMSE) fo r December between IDW and Ordinary Kriging methods – based on daily maps (see Appendix D for math ematical definitions of MAE and RMSE)17 Chapter 4 Neural network

Introduction Traditionally, the term “neural network” referred to a circuit of biologi cal neurons. The modern usage of the term often refers to networks made up of artificial neurons or “nodes.” These networks consist of interconnecting artifi cial neurons such as programm ing constructs that mimic the properties of biological neurons. Ne ural networks have been develo ped and used widely in both applied mathematics and applied science, includi ng predictive modeling base d on a time series of measurements – its application in the present work. Other application areas include system identification and control (vehic le control, process control), game-playing and decision making (backgammon, chess, racing), pattern recogniti on (radar systems, face identification, object recognition), sequence recognition (gesture, sp eech, handwritten text recognition), medical diagnosis, financial applications, data mining, and e-mail spam filtering. Thes e networks are always imbued with the capability of “learning” through “training,” w ith a variety of different methods employed. Although this may appear somewhat metaphys ical, it really boils dow n to a set of data and predictive equations, which, exercised daily, acquires “knowle dge” and an improved predictive ability through the simple acquisition and prediction of additional data. In this project the neural network’s application was limited to the prediction of historical PM 10 concentrations. An operational neural network, on the other hand, would make the daily predictions of tomorrow’s air pollution. This network would acquire more “ski ll” with time as additional daily relationships between PM 10 concentrations and meteor ological variables were encountered, processed, and predicted.

In the present work, an artificial neural netw ork called “EnviNNet” was de veloped and applied to the PM 10 monitoring record at the Central Phoenix (CP) site, with the goal of constructing an air pollution predictive system that might complement the deterministic system of the grid-based air quality model, discussed in Chapter 5. Four steps were necessary:

1. Assemble hourly records of PM 10 and meteorological variable s from one monitoring site,

preferably for three years or more (only one year was examined in this project).

2. Derive and calculate statistical relationships between the PM 10 concentrations and the

meteorological measurements.

3. Run a meteorological predictive model such as MM5 (see Chapter 5) for each day of the

historical record examined, determining quant itatively how well its predictions matched the
measurements.
4. Build the artificial neural network and th e meteorological modeling system such that

predictions of PM 10 are forthcoming for each day of the historical record. Compare these predictions with the measurements.

Should an organization wish to apply the neural network in a genuine air pollution predictive system, then the requisite software would need to be installed, the data gathering mechanisms built, the meteorological model would be set up for daily operational fore casts, and considerable staff training would have to take place. The remainder of this chapter discusses some aspects of these artificial neural networks and desc ribes how the EnviNNet was built.

Considering its wide usage in atmospheric applications (Gardne r and Dorling 1998), a type of neural network called “Multi La yer Perceptron (MLP)” was empl oyed. The MLP structure consists of an interconnected system of nodes (neurons ) within a hidden layer that employs non-linear18 continuous transfer functions c onnecting input and output vector s. (Appendix B provides more details.) Phoenix case study As the first step, PM 10 and meteorological data from the Ce ntral Phoenix (CP) monitoring station were assembled. EnviNNet was trained using a select ed set of appropriately pre-processed historic time series data at CP, whereupon it was used to predict PM 10 in the hindcasting mode. The input vector is made of the following variables:

1. Pollution data - These constitute the endogenous component of the NN (A component or

variable is called endogenous if it is explained within the model in which it appears). Time lags are selected by optimizing the data reconstr uction efficiency of the NN, so as to achieve the best trade-off between sensitivity and specif icity. Since these are multi-step operations, the measured pollutant levels are used in a rolli ng mode to generate the reconstructed values.

2. Meteorological measurements - These comprise the exogenous component of the NN, and

have a modulating effect on PM 10 concentrations ( exogenous refers to an action or object coming from outside a modeling system). The tim e lags were selected as before, and the balance between data reconstruction efficacy and network efficiency was optimized.

3. Statistical and descriptive indicators - These represent the behaviour of the “Central

Phoenix” neighbourhood via classes of input parameters that account for aggregated homogeneous spatio-temporal band s of microclimatic factors and air pollutant concentrations. The latter are dependent on the wi nd velocity, wind direction, and emissions specified as either weekday or weekend/holidays.

The meteorological and emissi ons characteristics of low, moderate, and elevated PM 10 concentrations at CP have to be accounted for in developing the NN, which makes the training strategy critical for reliable predictions. Since daily PM 10 concentrations vary widely throughout the month and year, specific input models had to be built by combining time section data of different years to make the NN robust and adaptive to the local climate. The data for training came from selected four-to-six-month windows, as well as from time periods that show noteworthy patterns such as exceptional events.

Input variables selected for training are:

• Historical series of PM 10 and meteorological data collect ed at CP during 2005-2006, prior to the dates selected for hindcasting;

• A selection of days based on analysis and corr elations of historical series representing PM 10 trends in the time window selected;

• Noteworthy patterns observed outside the time window of training (e.g. significant PM 10 peaks representative of exceptional events).

In this phase, the intrinsic parameters of the NN were estimated by numeric optimization until a satisfactory convergence was obtaine d. The weights and time lags so identified naturally depend on the site (e.g. loca l weather and neighbourhood emission character istics). The most important aspect of this particular traini ng is to obtain a non-lin ear dynamic relationship that connects input data with PM 10 output so that the NN can be ope rated in the hindcasting mode.19 Results The neural-network system gave dependable results for the Central Phoenix monitoring site. These results, compared with those of CMAQ in Chapter 6, are illustrated in Figure 4.1. PM 10 predictions in arid climates are difficult a nd a correlation coefficient of 0. 38 is better than average.

0 50 100 150 200 250 Observed PM10 concentrations [ μ g/m 3 ] 0 50 100 150 200 250Modeled PM10 concentrations [μg/m3] November, 2005

R 2 =1

R 2 =0.38

Figure 4.1. 24-hour average PM 10 predictions from the neural netw ork versus observations at the Central Phoenix site for 2006

Chapter 5 Air Quality Modeling

Introduction This chapter describes the deterministic, grid-b ased numerical models employed in simulating the central Phoenix PM 10 concentrations. The three kinds of models are mete orological, emissions, and air quality. In addition, near the ch apter’s end, the air qual ity model results are discussed in light of the interpolated concentrati on fields of Chapter 3.

The MM5/SMOKE/CMAQ modeling system was employed to simulate PM 10 concentrations in central Phoenix. This widely used, sophisticated platform has been a dopted as a regulatory tool. The predictive system consists of three integrated models: the Pennsylvani a State University and National Center for Atmospheric Research (NCAR) Mesoscale Me teorological Model - MM5 v3.7 (Grell and Dudia, 1994) for simula ting the weather, the Sparse Matr ix Operator Kernel Emissions – SMOKE v2.2 model for emissions processing (Car olina Environmental Program, 2005). and the Community Multiscale Air Qu ality Model – CMAQ v.4.5 (Byun a nd Ching, 1999) for simulating pollutant concentrations.

MM5 is a terrain-following model that simulates mesoscale atmos pheric circulation. It has been applied to a broad spectrum of studies, including land-sea breezes , mountain-valley circulation, and real-time weather forecasting. The data for mo del initialization and la teral boundary conditions were obtained from the operational NCEP/Eta m odel (40km grid resolution) and NCEP global surface and upper air observations every three and si x hours, respectively. (“NCEP” is the National20 Weather Service’s National Center for Environmental Prediction an d “Eta” is one of the global weather models.) Three-dimensional analysis and surface nudging with observed data were used to generate grid-based meteorological fields. MM5 provided the meteorol ogical input to CMAQ. Appendix C contains further modeling details.

The SMOKE model produces gridded, hourly emissi ons for CMAQ, using meteorological fields from MM5 and a host of emission estimating techniques. Tempora lly, spatially, and chemically resolved model-ready emissions were prepared using the recent Western Regional Air Partnership (WRAP) emissions inventory (including the BRAVO- 1999 study for the Mexico emissions inventory).

Developed by EPA in the 1990s, CMAQ is a multi- pollutant, multi-scale air quality model that contains techniques for simulating all atmospheric and land processe s that affect the transport, chemical reactions, and deposition of atmosphe ric pollutants and/or their precursors on both regional and urban scales. It is designed as a sc ience-based modeling tool for handling all the major pollutants such as ozone and particulate matter. This “Eulerian” or grid -based model produces hourly concentrations of a number of pollutants in the frame of a fixed, thre e-dimensional grid with uniformly sized horizontal grid cells and variable vertical layer thicknesses. This model has undergone numerous applications and performance evaluations, in Phoenix, in Arizona (Yuma and Douglas), in the Western United States, and throughout the world. Modeling domains The modeling domain was based on a La mbert Projection centered at (97 ° W, 40 ° N), with three nested domains having 36, 12 and 4km grids to predict meteorology, emissions, and air quality, respectively. In the nested simulations, the resu lts obtained from the outer domain were used as initial and boundary conditions for the inner doma in (Figure 5.1). Vertically 29 levels were specified with the layer closest to the ground being seven meters thick to capture boundary-layer processes. The inner domain with 4x4km grids wa s centered in the Salt River Valley in central Phoenix.

Supersite

West Phoenix

Durango Complex

West 43 rd Avenue Central Phoenix Supersite

West Phoenix

Durango Complex

West 43 rd Avenue Central Phoenix

Figure 5.1. Modeling domains for MM5/SMOKE/CMAQ simulations21 Modeling periods The modeling periods selected were constrained by the available air quality and health data in 2004- 2006. As was described in Chapter 2, there are five monitors in the study area with hourly PM 10 concentrations, but only two with hourly PM 2.5 concentrations. Therefore, the study examined only the PM 10 concentrations. One part of the project was to develop a nd test different interpolation techniques to build spatial concentration fields. These interpolation methods require from five to ten monitoring sites, another reason that the PM 2.5 data could not be used. In 2004 the sites of West Phoenix and West 43 rd Avenue had incomplete records, re stricting the analysis to 2005 and 2006. The 24-hour averages of PM 10 are shown for the five monito rs in Figure 5.2, whose names and designations are given below: CP – Central Phoenix DC – Durango Complex WF – West 43 rd Avenue WP – West Phoenix SS – Supersite

Different periods with high and low polluti on levels were selected for analysis.

• High pollution (above the National Ambient Air Quality Standard of 150μg/m 3 for a 24-

hour average); November 2005; December 12- 19, 2005; January 10-20, 2006; November 2-16, 2006

• Low pollution (less than 75μg/m 3 ); September 3-15, 2005; August 2005; March 11 - April

9, 2006; August, 2006

Figure 5.2. 24-hour averaged PM 10 observations for 2005-2006

The daily counts of emergency depa rtment visits and hospital admissi ons with an asthma diagnosis in Maricopa County, shown in Figure 5.3, we re divided into four categories:

• Worst : greater than the sum of the aver age and the standard deviation;

• Bad : between the annual average and the sum of the average and the st andard deviation;

• Moderate : between the average and the difference between this average and the standard deviation;

• Good : less than the difference between the average and the standard deviation. 0 75 150 225

1 366 Julian dayPM10 concentrations CP DC WF WP SS

2005 2006 0 75 150 225

1 366 Julian dayPM10 concentrations CP DC WF WP SS

2005 2006 μg/m 322 01/01/05 06/30/05 12/27/05 06/25/06 Day of the year 0 10 20 30 40Asthma events

Figure 5.3. Daily emergency department visits and ad missions to hospitals of children ages 0-17 in Maricopa County with an asthma diagnosis for 2005 – 2006

The average number of daily asthma visits is 9 w ith a range from 1 to 31. The seasonal pattern is virtually identical for both years. The good category occurs almost exclusively in July; the worst category is mostly limited to a period from the la st week of October until March 1, with a hiatus from December 20 to January 20, when the level of admissions stays in the moderate range.

Simple statistical analyses of the 24-hour PM 10 concentrations and daily asthma visits were carried out to help select the modeling periods. Curves in Figure 5.4 illustrate the variations throughout the two years, with a seasonal periodicity in both PM 10 and asthma being quite evident. The periods when the two data sets diverge are well displa yed. The highest number of asthma incidents and elevated PM 10 concentrations occurred together in November 2005 through January 2006. A decrease in asthma incidents near the end of December, while PM 10 concentrations remain elevated, may be due to behavioral changes during the Ch ristmas holidays. November 2005 had the highest number of asthma incidents; this month was select ed to represent “high po llution” conditions. It was more difficult to select the period of low polluti on because of incoherence in both data sets. The lowest PM 10 concentrations were recorded in August, while the lowest numbers of asthma incidents came in July. One period with both low PM 10 levels and low numbers of asthma cases was found in March-April 2006. This period was chosen to represent “low pollution” conditions.

The modeling was conducted for the two selected periods:

• November 2005, with high pollution, and

• March 11 - April 9, 2006, with low pollution.

Modeling runs were carried out with 31 hours of spin-up time for the weather predictions and 24 hours for the air quality modeling. “Spin-up” time m eans that the model run is begun before the day (or period) of interest. This allows its calculati ons to begin before the period of interest, reducing the influence of the initial c onditions on the model results. average 

### } standard deviation23 01/01/05 06/30/05 12/27/05 06/25/06 Day of the year 0 10 20 30 40Asthma events

0 50 100 150 200 250PM10 concentration [μg/m3]

Figure 5.4. PM 10 concentrations and asthma events (2005-2006)

Modeling results and discussion
Numerical simulations were carried out with the MM5/SMOKE/CMAQ modeling system for both high and low pollution periods and were compared with observations from the five permanent air monitoring stations. This model validation was n ecessary to estimate the model’s performance in the study region. Air pollutant concentrations de pend on meteorological conditions such as wind and rainfall, as well as on the strength and patter ns of air pollutant emissi ons. Consequently, in any air quality model’s performance evaluation, one must take into account the uncertainties in both the emission and meteorological models.

Model runs for two time period s: November 2005 and March 11 - April 9, 2006, show reasonable results in comparison with the observations . The scatter plots of simulated hourly PM 10 concentrations versus the observa tions at all five locations in central Phoenix are shown for both periods with their corresponding regr ession coefficients (Fig. 5.5).

CMAQ generally overestimated the low (<50μg/m 3 ) and underestimated the elevated (>150μg/m 3 ) observed PM 10 concentrations for both periods, with better correspondence between observed and modeled data for the high pollution period. Elevated PM 10 levels were a frequent occurrence at the Durango Complex (DC) and West 43 th Avenue (WF) sites, well a bove the other three sites in winter. For this reason the model results were also compared with the measurements in two groups:

• CP, WP, and SS, the three lowe r concentration sites, and

• DC and WF, the two higher con centration sites (Figure 5.6). November, 2005 March 11 - April 9, 200624 0 100 200 300 400 500 Observed PM10 concentrations [ μ g/m 3 ] 0 50 100 150 200 250Modeled PM10 concentrations [μg/m3] November, 2005

R 2 =1

R 2 =0.16

0 100 200 300 400 500 Observed PM10 concentrations [ μ g/m 3 ] 0 50 100 150 200 250Modeled PM10 concentrations [μg/m3] March 11 - April 9, 2006

R 2 =1

R 2 =0.12

Figure 5.5. Simulated and observed hourly PM 10 concentrations in centra l Phoenix for periods of high pollution (left) and low pollution (right) – from all five monitors

The coefficients of determination, which quantif y the goodness of the fit and are 1.0 for perfection, were less than 0.2 for all five sites within the stu dy domain with the 4km grids. Perhaps this level of performance can be attri buted to various inadequacies of the emissions inventory. Even a perfect description of the boundary laye r air flow from the meteorol ogical model would not produce excellent modeled air pollutant concen trations as long as the estimated PM 10 emissions were inaccurate in magnitude, space, and time. As th e WRAP emission inventories came in 36 and 12km grids, the 4km emissions for the central Phoenix m odel runs were obtained by setting them equal to one ninth of the 12km emissions (there are nine 4km grids within a 12x12km square).

PM10 - November, 2005

R 2

= 0.2078

R 2

= 1

R 2

= 0.33

R 2

= 0.2544

0 100 200

0 100 200 300 400 500 Observed [ug/m3]Predicted [ug/m3] CP

WP

SS PM10 - November, 2005

R 2

= 0.2056 R 2 = 1

R 2 = 0.1625

0 100 200

0 100 200 300 400 500 Observed [ug/m3]Predicted [ug/m3] DC

WF

Figure 5.6. Predicted and observed hourly PM 10 concentrations in central Phoenix for the two sets of monitors shown separately – high pollution period (November 2005)25 One explanation for the elevated PM 10 concentrations measured at th e DC and WF monitors is that weekday “industrial” and roadway emissions add to the concentration load ing already high from urban transport. The air quali ty model’s inability to simula te these higher concentrations, demonstrated by the near-horizontal angle of the DC and WF regression lines in the right figure above, may lie with an insufficient emission inve ntory. One should keep in mind that while CMAQ has produced excellent simulations for ozone and chemically speciated fine particulates, whose atmospheric behavior mimics that of gase ous pollutants its ability to simulate PM 10 has not been so good. The reasons for this lack-luste r performance are many, and include:

• a single, spatially uniform emi ssion value throughout a single gr id – at complete odds with

the spatial variability of real-world emissions;

• deposition phenomena that depend on the partic le size distribution, as different emission

sources produce particles with different distri butions; The average particle size of the total mass of PM 10 in central Phoenix is from five to seven microns, large e nough to be deposited somewhat rapidly to the ground, in contrast to the smaller particles comprising the fine fraction (smaller than 2.5 microns), which remain suspended considerably longer.

• the sporadic and unknown temporal distributio n of most fugitive s ources that cannot be

captured by any emissions inventory; and

• the heterogeneous spatial distribution of ambient PM 10 , exemplified by monitors influenced

by localized emissions that vary widely in thei r distance from the monitor, in their timing, in their release height, a nd in their magnitude.

Despite its shortcomings, a deterministic model su ch as MM5/CMAQ offers three advantages over its statistically based competitors, providing:

1. nearly complete geographic coverage, albeit li mited to the size of the grid, and adequate temporal coverage in the fo rm of hourly averag es, with simulated values of both meteorological variables and air pollutant concentrations calcul ated for each and every grid cell;

2. air quality predictions years or decades into the future, necessarily dependent on the reliability of a future-year emissions inventory; and

3. a way to quantify the benefits of specific ai r pollution control strategies, accomplished by modifying a future, base-case emissions inve ntory to reflect the anticipated emission reductions from the controls.

While only the first advantage is of relevance to the Children’s Health project, the other two should be considered in any general assessment of the utility of deterministic models.

A disadvantage of CMAQ is that, as a determinis tic model, it has none of the interpolation techniques that could better fit its output to the observed data. Furthermore, the emissions data and meteorological data that drive CMAQ each have c onsiderable uncertainties , the former to a far greater degree than the latter. Thes e uncertainties, as well as those i nherent in its embedded science, cause it to predict better for longer time scales (monthly, annual) than for shorter times (1-hour or 24-hour); and to be the most error prone with smaller grid cells. The "soccer goal" plot, recommended by USEPA (Guidance on the Use of M odels, 2007), is a convenient way to visualize26 model performance, as both bias and errors are shown on a single pl ot. Two statistics -- the Mean Fractional Error and the Mean Fractional Bias -- for both the high and low pollution periods are plotted for the different sites in Figure 5.7. As bias and error approach 0, the points are plotted closer to or within the "goal" represented here by the dashed boxes. (See Appendix D for the mathematical definitions of these and other statistics.) The plot shows reasonable model performance, with only a single site be ing out of the goal area -- West 43 th Avenue for the high pollution period.

0 25 50 75 100 125

-90 -60 -30 0 30 60 90
Mean BiasMean Absolute Error

0 25 50 75 100 125

High pollution Low pollution goal area

WF DC

WP SS CP WF DC WP SS CP

0 25 50 75 100 125

-90 -60 -30 0 30 60 90
Mean BiasMean Absolute Error

0 25 50 75 100 125

High pollution Low pollution goal area

WF DC

WP SS CP WF DC WP SS CP goal area

WF DC

WP SS CP WF DC WP SS CP

Figure 5.7. The "soccer goal" plot measures both bias and error for central Phoenix

In addition to this visual estimate of performan ce, the different statistics that were calculated to estimate the model performance for the study peri ods also show reasonable model performance (Table 5.1).

Of all these statistics, perhaps th e most telling one is the Index of Agreement (IA): its values range from 0.46 to 0.72, with perfect agreement being 1.0. The best IA was achieved at WP during the winter and at DC during the spring. The better know n Coefficients of Determination (RSQ), also known as the linear “correlation coefficient” abbreviated as “R 2 , ranges from 0.10 to 0.33, indicative of a slightly positive correlation. This weak correlation is often found in matching numerically modeled concentrations with the measured ones in air pollution studies, especially for PM 10
. The Mean Absolute Errors (MAE) are all less than 30 except for DC and WF in November 2005. The local emissions not taken into account with sufficient spatial and temporal accuracy in the emissions inventory contribute to this disagreement. CMAQ overestimated the PM 10 concentrations at CP and SS: the Modified Coeffici ent of Determination (RSQ*) is greater than 1; but it underestimated the observations at WF and DC , where RSQ* is less than 1 for both high and low pollution periods. At West Phoenix (WP) the model under-estimated the higher concentrations measured in November 2005 but over-estimated the lower ones in March-April 2006. The Root Mean Square Errors (RMSE) ranged from 26 to 42, except for the two problematic sites of WF and DC in November 2005, which had an omalously high values about 90.27 Table 5.1. Statistical measures in co mparing model-predicted hourly PM 10 concentrations to observations for both study pe riods at different sites



### Monitor MAE RMSE RS Q RS Q * IA Perfect agreement 

### 0 n a 1.00 n a 1

Central Phoenix 26.1 2 34.39 0.21 1.45 0.68 West Phoenu x 27.81 41.99 0.33 0.60 0.7 2 S u p er Site 23.48 30.0 6 0.25 3.30 0.6 2 West 43rd Avenue 62.21 89.85 0.16 0.0 6 0.51 D uran g

o C p

m

p lex 58.49 85.80 0.21 0.68 0.55

Central Phoenix 25.00 33.0 6 0.10 4.28 0.4 6 West Phoenux 22.78 30.25 0.16 2.39 0.59 S u p er Site 22.19 27.15 0.12 5.01 0.47 West 43rd Avenue 25.31 26.47 0.17 0.4 2 0.60 Durango Cpmplex 23.94 35.86 0.25 0.54 0.68 November 2005

March 11- April 9, 2006

Abreviation: MAE - Mean Absolute Error; RMSE - Root Mean Square Error; RSQ - Coefficient of Determination (aka R 2 ); RSQ* - Modified Coefficient of Determination 2 ; IA - Index of Agreement
2 Modified Coefficient of Determination: values greater than 1 mean the model ov erestimated (orange); values less than 1 mean the model underestimated (green) The health data are available only on a daily basis, and consequently, the prediction of accurate 24- hour averages of PM 10 was of greater importa nce than the hourly valu es. Figure 3 shows the 24- hour averaged Mean Bias for both high and low pollution periods . The difference between daily calculated and observed concentrati ons can be substantia l. While the average MB for 30 days was in the range of -15 to 15, except at WF and DC in November (see Fig. 5.7), the 24-hour averaged MB doubled for some days and is almost zero for the others, as can be seen from Figure 5.8.

The model underestimations (the negative MB) were most pronounced at Durango Complex and West 43rd Avenue, and a weekly periodicity can be found for the MB. The weekday industrial and roadway emissions add to already-high concentra tion loading from urban transport. These two regulatory problematic sites, which represent the worst PM 10 air quality in metropolitan Phoenix, have relatively high frequencies of violations of the NAAQS, in part because of somewhat dense localized emissions from the extractive and material handling i ndustries along the Salt River. Some of these sporadic and unknown tem porally distributed fugitive s ources cannot be captured by an emissions inventory, and they contribute to the heterogeneous spatial distribution of ambient PM 10 and hence to prediction errors. The averag e particle size of the total mass of PM 10 in central Phoenix is from five to seven microns, la rge enough to be deposited on the ground somewhat rapidly, in contrast to the smaller particles comprising the fi ne fraction (smaller than 2.5 microns) that remains suspended for considerably longer ti mes. The aerosol module in CMAQ considers both PM 2.5 and PM 10 , but represents the particle size distribut ion as the superposition of three lognormal sub-distributions called modes. The model takes into account on ly one mode of the "coarse particles" (2.5-10 μ m) and assigns the sa me deposition velocity.28 -150 -125 -100

-75 -50 -25

0 25 50 Mean Bias

1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 Day November, 2005

CP DC WF WP SS

-150 -125 -100

-75 -50 -25

0 25 50 Mean Bias

1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 Day March 11 - April 9, 2006

CP DC WF WP SS
Figure 5.8. Daily Mean Bias for high po llution conditions (top) and low pollution

Figures 5.9 and 5.10 show the CMAQ predictions co mpared with 24-hour averaged observations for both high and low pollution periods. The calculated data were extracted from the grid cell where the corresponding monitor is located. A good corresponde nce between observed and calculated data for WP and CP sites can be found in November 200 5 (Fig. 5.9). The model overestimates the 24-hour averaged observed values at SS and strongly unde restimates the values at WF and DC sites. The same behavior can be found for the low pollu tion case (Fig. 5.10). CMAQ overestimates the 24- hour averaged measured values at SS, CP, WP and underestimates th e values at WF and DC sites for several days. Nonetheless, th e model's predicted values generally follow the trends of the observations. In spite of MM5’s inability to simu late the higher wind speeds for several days, the model did predict some higher wind speeds, reflected by the low PM 10 concentrations calculated by CMAQ for 11, 26-27 of November 2005, March 30, and 2, 5-6 of April, 2006.29 0 10 20 30 0 50 100 150 CP

Da yPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 WP

Da yPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 DC

DayPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 WF

Da yPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 SS

Da yPM10 concentration (μg/m3) Observation CMAQ Observation CMAQ

Figure 5.9. CMAQ predictions versus observations of PM 10 concentrations at five sites for the high pollution period of November 2005. The lines are the median values of the 24 hourly predictions for each day; the box plots are 25 th to 75 th percentiles of the 24 predictions.30 0 10 20 30 0 50 100 150 CP

Da yPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 WP

Da yPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 DC

Da yPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 WP

Da yPM10 concentration (μg/m3)

0 10 20 30 0 50 100 150 SS

Da yPM10 concentration (μg/m3) Observation CMAQ Observation CMAQ

Figure 5.10. CMAQ predictions versus observations of PM 10 concentrations at five sites for the low pollution period of March 11 – April 9, 2006. The lines are the median values of the 24 hourly predictions for each day; the box plots are 25 th to 75 th percentiles of the 24 predictions.31 The PM 10 concentrations were averaged for the five sites for the entire period for every hour. Figure 5.11 presents the averaged con centrations from the five monito rs in central Phoenix with the spatially-averaged PM 10 concentrations from the numerical modeling.

0 50 100 150

0 2 4 6 8 10 12 14 16 18 20 22 Time [h] Observed Calculated

Figure 5.11. Hourly variation of the five-s ite averages, calculated and observed PM 10 concentrations for November 2005

The model is able to produce the morning and even ing maxima at the right times, in agreement with the observations, but it underestimates the magn itude of the peaks. Th ere are at least two explanations for this lack of agreement. First, consider the weather-pr edicting MM5 model, which is in need of some modifications itself. Its pred icted wind fields have speeds considerably higher than the measured values during the morning a nd evening transitions (for more details see Appendix C). Provided that the wind speeds remain below the dust resuspen sion threshold - roughly 15 to 20 miles per hour (7 to 10 me ters per second) - higher wind speeds generally lead to lower PM 10 concentrations; hence, the underestimation. S econd, consider that MM5 does not differentiate among the various types of urban land surface. Th e treatment of urban land in modeling grids results in too uniform a surface roughness and other properties, which are not representative of the local-scale variation. Sub-categorie s are needed for the so called "u rban" land characteristics. These unrealistic uniformities can cause incorrect transport of pollutants, as well as too little or too much vertical mixing, both of which affect their resultant concentrations. Especially for stable conditions, which can last for several days, improvements in the application of this model to metropolitan Phoenix would be helpful.

CMAQ contains “accurate” and “state-of-the-scien ce” parameterizations of atmospheric processes affecting transport, transformation, and deposition of po llutants. Process analysis is a useful tool for understanding the contributions of differe nt physical processes to the total PM 10 concentrations. Through process analysis the contribution th at each phenomenon makes to the resultant PM 10 concentration can be quantifie d (Fig. 5.12). The emissions (EMI S), dry deposition (DDEP), and vertical diffusion (VDIF) give a jo int contribution ten times greater than the combined contributions of advection (horizontal HADV + vertical VADV), horizontal diffusion (HDIF), atmospheric chemistry (CHEM), aerosol composition (AER O), and cloud chemistry (CLDS) combined. Inadequate parameterizations of both the diffusion and the depositional processes in the planetary boundary layer, combined with uncertain emission i nventories, may be responsible for the model’s performance. μg/m 332 Process Analysis - November 7-9, 2005

-800 -600 -400 -200 0 200 400 600 800

0 12 24 36 48 60 72 time [h]Concentration rate [ug/m3]

-60 -40 -20 0 20 40 60

DDEP EMIS VDIF AERO CLDS
HADV VADV HDIF CHEM PM10

Process Analysis - April 1-3, 2006

-800 -600 -400 -200 0 200 400 600 800

0 12 24 36 48 60 72 time [h]Concentration rate [ug/m3]

-60 -40 -20 0 20 40 60

DDEP EMIS VDIF AERO CLDS
HADV VADV HDIF CHEM PM10

Figure 5.12. Contributions of different processes to PM 10 concentrations for 3-days periods in November 2005 and April 2006: the va lues are averaged from several grid cells in central Phoenix, with EMIS, DDEP, and VDIF plotted on the left Y-axis, with all others on the right

One reason for the differences between the modelled and observed PM 10 concentrations is that measured point values cannot represent the la rger grid-based volumes of 36, 12, or 4 km surrounding the monitoring sites, given the large spatial inhomogeneities already noted. The spatial distribution of the Index of Agreement (IA) between the in terpolated and modeled PM10 concentrations is shown in Figure 5.13 for November 2005. The inter polated surfaces were constructed from the observati ons and mapped into census tract concentrations, from which the average value for each modeling grid cell was cal culated. The map shows the IA between the interpolated PM10 concentrations from the IDW a nd the CMAQ estimates for each 4 km grid cell. A very good correspondence can be found with the IA greater than 0.5 for the whole domain. The modeled data fit very well to the interpolated surfaces (IA between 0.6 - 0.7) in the northeast part (WP, CP, SS) and some disagreement can be seen at the south-west (l ess than 0.5). The CMAQ model overestimates the PM10 concentrations in the northeast around CP and SS, unlike in the EMIS
, DDEP
, VDIF

EMIS
, DDEP
, VDIF33 southwest around WF and DC, where it generally unde restimates the observat ions, as can be seen from the map of Mean Bias for November 2005 (Figure 5.14).

Figure 5.13. The spatial distribu tion of the Index of Agreem ent between CMAQ and Inverse Distance Weighting (IDW) interpolated surfaces

Figure 5.14. The spatial di stribution of Mean Bias between CM AQ and Inverse Distance Weighting (IDW) interpolated surfaces34 One way to improve the numerical predictions for the central part of the study area would be to combine CMAQ and IDW surfaces. Application of the Daily Mean Biases to CMAQ predictions produces a more realistic pattern of PM 10 concentrations and reduces the disagreement between calculated and observed data. Examples of this t echnique are shown in Figures 5.14-5.18 for all five permanent sites for November 2005.

ErrBar=25%

0 20 40 60 80 100 120

0 5 10 15 20 25 30 Day of November, 2005PM10 concentration CP-obs CMAQ CMAQ-new

Figure 5.15. Comparison of the daily, “simulated” PM 10 concentrations with the observations at Central Phoenix: violet is straight CMAQ; blue is “CMAQ-new” - CMAQ combined with interpolation

ErrBar=25%

0 50 100 150

0 5 10 15 20 25 30 Day of November, 2005PM10 concentration WP-obs CMAQ CMAQ-new

Figure 5.16. Comparison of the daily, “simulated” PM 10 concentrations with the observations at West Phoenix (WP): violet is straight CMAQ; blue is “CMAQ-new” - CMAQ combined with interpolation

ErrBar=25%

0 20 40 60 80 100

0 5 10 15 20 25 30 Day of November, 2005PM10 concentration SS-obs CMAQ CMAQ-new

Figure 5.17. Comparison of the daily, “simulated” PM 10 concentrations with the observations at Supersite (SS): violet is straight CMAQ; bl ue is “CMAQ-new” - CMAQ combined with interpolation.35 ErrBar=25%

0 50 100 150 200 250 300

0 5 10 15 20 25 30 Day of November, 2005PM10 concentration WF-obs CMAQ CMAQ-new

Figure 5.18. Comparison of the daily, “simulated” PM 10 concentrations with the observations at West 43 rd Avenue (WF): violet is straight CMAQ; blue is “CMAQ-new” - CMAQ combined with interpolation

ErrBar=25%

0 50 100 150 200 250 300

0 5 10 15 20 25 30 Day of November, 2005PM10 concentration DC-obs CMAQ CMAQ-new

Figure 5.19. Comparison of the daily, “simulated” PM 10 concentrations with the observations at Durango Complex (DC): violet is straight CMAQ; blue is “CMAQ-new” - CMAQ combined with interpolation

Although this "nudging" technique can be applied to historical periods, it cannot be of help in forecasting, as the observations that would provide th e interpolated surfaces are yet to be made. Any improved health warning system would have to rely either on a statistical-based prediction such as the neural-network EnviNNet or on CM AQ or another deterministic model prediction driven by the meteorological model’s pred icted weather for the next day.

Conclusions For the meteorological fields :

• MM5 fails to predict high wind speeds, perhap s due to the inconsis tency in global model

predictions used as i nput and boundary conditions.

• MM5 slightly overestimates the low wind speed s during the generally stable conditions of

the night and early morning; furthermore, it can not faithfully simulate rapid changes of wind direction.

• MM5 does accurately simulate the temperature fi eld, albeit with slight underestimates of

nocturnal temperatures.36 For the PM 10 concentration fields :

• CMAQ generally underestimates the highe r observed PM10 concentrations and

overestimates the lower ones.

• CMAQ captures the diurnal variation of PM10 concentrations.

• CMAQ especially underestimates the ately simula tes accurate such as the neural-network

based EnviNNet or on a CMAQ or CMAQ-type deterministic model driven byelevated PM10 concentrations at West 43rd Avenue and the Durango Complex in November 2005. Dust emissions in the vicinity of the monitors, but inadequately accounted for in the emissions inventory, are the likely explanation.

• CMAQ adequately simulates the surface PM10 di stribution in the centr al Phoenix study area

for the needs of this project.

Three difficulties plague the appl ication of a deterministic system such as MM5/SMOKE/CMAQ to air quality surveillanc e and management:

• First, these models need pr ecise initial conditions, boundary conditions, and emissions to

provide reasonable predictions;

• Second, they are computationa lly demanding and somewhat labor-intensive to operate,

troubleshoot, update, and maintain; and

• Third, the meteorological model MM5 is actually designed to predict larger-scale, regional

phenomena, not the more local-scale characteris tics that prevail in and around specific urban air pollution monitoring sites.

As done in this project, the importance of the initial conditions can be reduced through spin-up times. As also done in this project, the impor tance of boundary conditions can be diminished through choosing larger modeling domains. The impor tance of spatially and temporally accurate hourly emissions, however, cannot be overlooked: th eir contribution to the performance of air quality model predictions is paramount . In particular, airsheds in ar id and semiarid areas, including most of the southwestern U.S. a nd much of its intermountain West, ar e subject to pervasive fugitive emissions of dust under all meteorological condi tions and to windblown dust under more turbulent and stormy conditions. Locally suspended dust from widespread fugitive emissions is a major component of everyday particulate matter concentr ations in this regio n. Production of soil dust aerosols depends on the nature of the mechanical disturbance, the wind ener gy, and the soil surface properties. Improvements in quan tifying both human-caused, mechanically generated emissions, as well as windblown dust emissions, would improve deterministic model performance in the study

area.
The computational demands of deterministic models cannot be casually dism issed, in spite of the recent, substantial improvements in speed and power stemming from pa rallelization and other techniques. Accurate weather for ecasting on the synoptic scale of one to five days is now a facet of our daily lives; still, the comput ational time needed for more precise predictions on the local scale remains about equal to the typica l real-world time of the local phenomena themselves. These phenomena include thermally driven upslope and valley winds, stagnation, extremely low wind speeds, stratification, and canaliza tion effects – all of which last for hours and require equal hours of simulation. Considering human e xposure, in spite of the general lack of short-term health37 standards (let’s say from one to six hours), the hi gher concentrations occurring in a few to several hours may be as important or even more so than the exposures of 24 hours or longer.

Meteorology’s critical control of air pollutant concen trations is well know n: the emission and chemical transformation of particles, their transp ort and dispersion – both hor izontal and vertical -- and their loss through dry and wet deposition are al l consequences of the weather. During the meteorological simulations of this project, it ha s become evident that the MM5 model would benefit from a number of modifications. For example, disagreements between a measured wind vector (speed and direction) and the MM5 volume-averag ed predictions for the 36, 12, or 4 km grids are commonly encountered. MM5 is not scientifically formulated to reproduce the exact wind speed and direction at each measuremen t point; rather, it has been designe d to generate su b-regional and regional wind, temperature, and moistu re fields that are internally c onsistent and strictly satisfy the laws of mass, momentum, and energy conserva tion. Conversely, local wind measurements are seldom sited in such a manner as to be representative of a la rge volume of space surrounding the observation site. In the complex terrain of me tropolitan Phoenix and its surrounding mountains, because of terrain forcing and orographic eff ects, the MM5 volume-averaged meteorological estimates have greater uncertainties than in flat terrain.

Chapter 6 Neutral Network versus CMAQ Predictions
For these comparisons November 2005 was selecte d, considering that the late fall and winter months tend to have the highest PM 10 concentrations and the most hospital visits for respiratory illnesses. EnviNNet predictions were limited to a single site, Central Phoenix. Furthermore, this November has both high and low PM 10 days, including two periods of elevated PM 10 – November 10 and November 22-23, for which statistics and gr aphics have been produced. The daily averaged values are presented first, given that the health data for comparison with PM 10 are available on a daily basis and that the National Ambient Air Quality Standard is expressed as a 24-hour average, midnight to midnight. The data are from the Mari copa County Air Quality Department, whose staff has verified their accuracy thr ough quality control procedures and has reported them to EPA’s AIRS system. Figure 6.1 shows the time series of these daily predictions by EnviNNet and CMAQ in comparison with the daily observations; Fi gure 6.2 presents comparisons of the hourly observations from November 2005 with the two models.

0 20 40 60 80 100

1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 DayPM10 concentration[ug/m3] EnviNNet CMAQ Observations

Figure 6.1. Time series of 24-hour average PM 10 concentrations predicted by CMAQ and EnviNNet with observations at the Central Phoenix site38 0 50 100 150 200 250 Observed PM10 concentrations [ μ g/m 3 ] 0 50 100 150 200 250Modeled PM10 concentrations [μg/m3] November, 2005

R 2 =1

R 2 =0.21

0 50 100 150 200 250 Observed PM10 concentrations [ μ g/m 3 ] 0 50 100 150 200 250Modeled PM10 concentrations [μg/m3] November, 2005

R 2 =1

R 2 =0.38

Figure 6.2. Scatter plots of hourly PM 10 concentrations predicted by CMAQ (left) and EnviNNet with observations: November 2005, Central Phoenix site

Predicting the higher PM 10 concentrations well is more important than doing well on the low ones. To this end the model’s predictions for elevated PM 10 concentrations were evaluated for the two high pollution periods of Novemb er 6-10 and November 22-26. The results are shown in Figures 6.3 and 6.4; the performance measures for the thr ee time periods and the two models are in Table 6.1, are defined in Appendix D, and are discussed below.

The calculations were made for one month and fo r the two five-day periods with elevated PM 10 concentrations. The Index of Ag reement (IA) is greater than 0.6, which shows good correspondence between the calculated and observed data. Generall y, EnviNNet gives better IA in comparison with CMAQ for all periods considered here. The Mean Absolute Errors (MAE) are less than 32 for CMAQ and less than 24 for EnviNNet. The Root Mean Square Errors (RMSE) are in the range of 25 – 40 for different periods. The statistical mode l EnviNnet yields smaller errors than CMAQ.

Table 6.1. Statistics of CMAQ and EnviNNet predictions of hourly PM 10 concentrations at the Central Phoenix site



### Model MAE RMSE RS Q RS Q * IA Period CMA Q 26.1 2 34.40 0.21 1.45 0.68 November EnviNNet 19.01 25.0 2 0.39 0.82 0.77 November

CMA Q 32.24 40.18 0.14 1.05 0.61 6-10 Nov. EnviNNet 20.39 28.04 0.42 0.75 0.79 6-10 Nov.

CMA Q 29.33 39.17 0.10 1.19 0.60 22-26 Nov

. EnviNNet 24.2 6 34.01 0.13 0.78 0.61 22-26 Nov

. Abreviation: MAE - Mean Absolute Error; RMSE - Root Mean Square Error; RSQ - Coefficient of Determination (aka R 2 ); RSQ* - Modified Coefficient of Determination 2 ; IA - Index of Agreement
2 Modified Coefficient of Determination: values greater than 1 mean the model ov erestimated (orange); values less than 1 mean the model underestimated (green)39 November 6-10, 2005

0 50 100 150 200 250

0 12 24 36 48 60 72 84 96 108 120 Time [h]PM10 concentration[ug/m3]

Observed CMAQ EnvirNNet Figure 6.3. Comparison of hourly PM 10 concentrations predicted by Models 3 and EnviNNet with the observations of November 6-10, 2005, Central Phoenix site

November 22-26, 2005

0 50 100 150 200 250

0 12 24 36 48 60 72 84 96 108 120 Time [h]PM10 concentration[ug/m3]

Observed CMAQ EnvirNNet Figure 6.4. Comparison of hourly PM 10 concentrations predicted by CMAQ and EnviNNet with the observations of November 22- 26, 2005, Central Phoenix site

The neural network satisf actorily predicts the PM 10 peaks in contrast to CMAQ, which has problems capturing them satisfactorily. A lthough the deterministic model may have given much better predictions under ideal initial, boundary and pollution inventory c onditions with all s cales of motion fully resolved, the present status of the model is far from this state. In addition, CMAQ cannot resolve scales smaller than the size of the grid, which in this case is 4x4 km.

Conclusions As a workable alternative to the deterministic MM5/CMAQ, the neural network proved to be somewhat better than CMAQ in predicting the moderate to high PM 10 concentrations. Its principal advantage in an asthma warning system is that it would not require an em issions inventory or the daily computational and staff demands of the grid -based, deterministic modeling system. The neural network is much easier and quicker to use than CMAQ and could be partially automated for issuing health warnings. Its shortcoming lies in its lim ited geographical coverage, especially for pollutants with steep concentration gradients such as PM 10 . In addition, only one site was examined, suggesting that similar work at additional sites w ould be advisable before any claims of area-wide viability could be supported.40 Chapter 7 Linking asthma with PM 10 concentrations

Introduction The goal of the health effects analysis is to determine whether asthma incidents (primarily emergency department (ER) visits and hospital admissions with a diagnosis of asthma) can be quantitatively linked to the air quality measurem ents obtained from Maricopa County, ADEQ, and ASU monitoring sites. The primary health data are the asthma incident reports from the Arizona Department of Health Services (ADHS). This data includes hospita lization admissions and emergency room visits from throughout the state as well as patient demogra phics e.g., gender, age,
etc.

The health effects analysis was conducted by th e Center for Health Information & Research (CHIR), part of the School of Computing and Info rmatics in Arizona Stat e University’s Ira A. Fulton School of Engineering. CHIR is a multidisci plinary team whose interests are health care, clinical quality, the health care workforce, occupational illness and injury, medical malpractice, health care economics, and disab ility. CHIR is home of the Ar izona Health Query (AZHQ), a community/university database effo rt. More than fifty organizations voluntarily share their data with ASU to create this integrated database. Each partnership is governed by contracts that comply with the federal Health Insurance Portability and Accountability Act (HIPAA), and stringent security is maintained to avoid unauthorized di sclosures. AZHQ contains health information on millions of Arizona adults and children and pr ovides the unique capability of tracking patients across health care systems and over time. As the da ta system continues to expand to include more services, CHIR efforts are focused on building addi tional partnerships with large employers, third party administrators, and physician groups. The ADH S routinely reports its data to the AZHQ.

AZHQ data has already been utili zed to study the respiratory eff ects of agricultural burning in northern Mexico and southwest Arizona on Yuma County, AZ children study funded by the Southwest Center for Environmental Research and Policy (SCERP). Using AZHQ, we examined the frequency of health care visits, emergency de partment utilization’ and hospitalizations for respiratory illnesses for children ages 0-17 during th e burning and compared this frequency with the prevalence of respiratory illness before and after exposure to agricultural burning. A separate study was conducted that linked asthma inci dents to air quality in the Casa Grande area of Pinal County, AZ.

Data In the present work, the ADHS and CHIR conducted a high-level scan of air quality monitoring data for PM 10 and asthma incidents from hospital discharg e data [and ER admission data] to determine the relationship, if any, between poor air quality and higher rates of asthma admissions. Preliminary analyses also helped researchers define the de mographic groups, sample si ze, and data cleaning strategies that were subsequently applied in th e final analyses of the work being reported here. Based on these preliminary results and discussions, research focused on ADHS asthma incident data occurring within the 168 census tract areas of central Phoenix from January 1, 2005 through December 31, 2006, a period for which complete air m onitoring records were available. All asthma events within five miles of a continuous PM10 monitor were included.

In addition to the standard cleaning applied by CH IR to correct errors, additional cleaning steps were taken. For example, there were cases where th e same patient went to an emergency room and was later admitted to the hospital. This single patient had multiple records on the same day with his41 status changing from emergency to inpatient. Su ch redundant records were deleted to leave but a single asthma incident. Also, those subsequent medical encounters occurring sufficiently long after the initial visit were considered to be a distinct incident. This time period is commonly referred to as a “washout” period. We explored several wa shout periods in this study. Figure 7.1 shows the change in the number of reporte d incidents for selected washout periods. As the incident changes are small for washout periods of se ven days or greater, we selected a washout period of seven days in this work. If the same person had multiple reco rds within the seven days , only the initial record was kept.

Figure 7.1. Washout periods for asthma incidents within two and five miles of a continuous PM 10 monitor

Asthma incidents were described by age group, day of the week, month, and place of service (emergency department and inpatient admission), resulting in a dataset containing about 2,000 events for the two-mile radius and 6,000 events for the five-mile radius for children aged 0-17 years. See Figures 7.2 - 7.5. Gap (days)Number of cases

28 21 14 7 1 11000

10000

9000

8000

7000

6000

5000

4000

3000 Variable 2-mile 5-mile ADHS - Number of Records42 (a) Asthma incidents by age group (b) Asthma incidents by day of week

(c) Asthma incidents by month (d) Asthma incidents by place of service

(e) Asthma incidents by age group and time Figure 7.2. Asthma incidents by different characteristic within two miles of PM 10 monitoring sites43 (a) Asthma incidents by age group (b) Asthma incidents by day of week

(c) Asthma incidents by month (d) Asthma incidents by place of service

(e) Asthma incidents by age group and time Figure 7.3. Asthma incidents by different ch aracteristic within five miles of PM 10 monitoring sites44 (a) Asthma incidents by zip code

(b) Asthma inci dents by zip code and place of service

(c) Asthma incidents by ag e group and place of service Figure 7.4. Asthma incidents by zip code within two and five miles of PM 10 monitoring sites Age Ambulatory Emergency Inpatient Missing Other Group AZH Q ADH S AZH Q ADH S A ZH

Q A DH S A ZH Q A DHS A ZH Q ADHS 0 - 4 810 - 184 685 95 155 1 - 4 - 5 - 9 868 - 185 580 57 117 7 - 3 - 10 - 14 574 - 107 338 29 64 3 - 5 - 15 - 17 184 - 47 115 6 11 2 - 3 - 18 - 24 173 - 110 222 25 33 - - 1 - 25 - 34 192 - 113 216 27 32 4 - 2 - 35 - 44 262 - 120 227 34 45 4 - 4 - 45 - 54 211

- 73 160 4 8 6 9 1

- 5

- 55 - 64 176 - 51 94 35 85 1 - 5 - 65 - 74 99 - 22 57 17 41 - - - -

75

+ 183 - 24 23 21 31 - - 3 - Total 3732 - 1036 2717 394 683 23 - 35 - ZIP Am bulatory Em ergency Inpatient Missing Other A

ZH

Q ADH S AZH Q A DH

S A ZH

Q A DH

S A ZH

Q A DH

S A ZH

Q A DHS 85003 108 - 47 93 23 25 - - 2 - 85004 54 - 15 61 13 22 - - 1 - 85006 441 - 145 428 54 86 - - 2 - 85009 938 - 257 594 79 128 3 - 4 - 85019 347 - 100 326 32 89 - - 7 - 85031 546 - 131 336 55 83 5 - 11 - 85041 753 - 218 553 79 165 12 - 2 - 85043 267 - 6 9 21 8 3 2 47 2

- 5

- 85339 278 - 54 108 27 38 1 - 1 - Total 3732 - 1036 2717 394 683 23 - 35 - 0 200 400 600 800 1000 1200 1400

123456789 ADHS
AZHQ N o. Z IP A D

H

S A ZH

Q 1 850 03 118 180 2 850 04 83 83 3 850 06 514 642 4 850 09 722 1 281 5 850 19 415 486 6 850 31 419 748 7 850 41 718 1 064 8 850 43 265 375 9 853 39 146 361 T otal 3400 5 22045 (a) Asthma incidents by zip code

(b) Asthma incidents by zip code and place of service

No. ZIP ADHS AZHQ

1 85003 118 180 2 85004 83 83 3 85006 514 642 4 85007 331 461 5 85008 1049 1205 6 85009 722 1281 7 85012 102 323 8 85013 231 342 9 85014 325 423 10 85015 596 993 11 85016 421 434 12 85017 578 736 13 85018 310 389 14 85019 415 486 15 85031 419 748 16 85033 821 1213 17 85034 224 226 18 85035 633 1018 19 85040 624 1022 20 85041 718 1064 21 85043 265 375 22 85301 788 1509 23 85339 146 361 Total 10433 15514 0 200 400 600 800 1000 1200 1400 1600

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 ADHS
AZHQ

ZIP Am bulatory Em ergency Inpatient Missing Other A

ZHQ ADHS AZH Q A DHS A ZHQ A DHS A ZH Q A DHS AZHQ A DH S 85003 108 - 47 93 23 25 - - 2 - 8500 4 54 - 15 61 13 22 - - 1 - 85006 441 - 145 428 54 86 - - 2 - 8500 7 302 - 112 269 41 62 1 - 5 - 85008 806 - 265 861 114 188 10 - 10 - 8500 9 938 - 257 594 79 128 3 - 4 - 85012 257 - 42 72 23 30 - - 1 - 85013 19 5 - 102 19 5 4 2 36 - - 3 - 8501 4 251 - 105 274 63 51 1 - 3 - 85015 647 - 246 514 86 82 2 - 12 - 8501 6 308 - 78 332 38 89 2 - 8 - 85017 494 - 175 473 59 105 1 - 7 - 8501 8 262 - 67 235 56 75 1 - 3 - 85019 347 - 100 326 32 89 - - 7 - 85031 546 - 131 336 55 83 5 - 11 - 85033 891 - 240 663 64 158 7 - 11 - 85034 161 - 49 184 13 40 - - 3 - 8503 5 757 - 191 534 58 99 3 - 9 - 85040 713 - 207 508 86 116 7 - 9 - 85041 753 - 218 553 79 165 12 - 2 - 85043 267 - 69 218 32 47 2 - 5 - 85301 1031 - 343 602 115 186 7 - 13 - 85339 278 - 54 108 27 38 1 - 1 - Total 10807 - 3258 843 3 1252 2000 6 5 - 132 -46 (c) Asthma incidents by age group and place of service Figure 7.5. Comparison of AZHQ and ADHS asthma incidents within two and five miles of PM 10 monitoring sites

Analysis Objectives The processed asthma incidence data by census tract were linked with the PM 10 concentrations interpolated from the five continuous monitors. In addition to the daily PM 10 concentrations on the asthma “event” day, 24-hour averages of PM 10 from one to seven days be fore the event were also examined in the statistical analysis (these variables are referred to as “lags”). To detect the effects of PM 10 for each asthma incident while accounting fo r the confounding variables of seasonality, day of the week, patient-level covariates (e.g., age, gender, ethnicity), and the presence of other air pollutants, a case-crossover analys is was conducted. This type of analysis associates a reference time period as a referent (or control) for each case. Consequently, each patient was treated as a matched case-referent pair with control exposures obtained from the same patient in different time periods. The time-stratified design used a 28-day period or “stratum size” with three referents selected from the same day of the week as the cas e. With an appropriate time-stratified period, the case-crossover analysis controlled for long-term tr ends, seasonal effects, and various epidemiologic covariates that change slow ly with time by design (e.g., lifestyle behaviors, diet).

Analysis Methods Background The case-crossover design, similar to a crossove r and matched-pair cas e-control studies, was developed to study transient short-term exposure eff ects on the risk of rare acute events (Maclure, 1991). The case is a person with the event of interest at a certain time called the “hazard period”, while the control (called the referent) is the sa me person at a different time called the “control period.” The key feature of this design is that each case serves as its own control. The exposure information for each subject during the hazard period is compared to the exposure information during the control period, de fined as the referent period for that event. With an appropriate referent A

ge A m

bulatory E m

ergency Inpatient M issing O ther G

roup A ZH Q A D

H

S A ZH Q A D

HS A ZH Q A D

H

S A ZH Q A D

H

S A ZH Q A D

HS 0 - 4 2456 - 533 1966 301 446 2 - 6 - 5 - 9 2387 - 479 1006 157 167 12 - 14 - 10 - 14 1572 - 320 332 90 34 7 - 16 - 15 - 1 7 513 - 134 663 17 78 5 - 11 - 18 - 24 496 - 308 899 64 112 4 - 5 - 25 - 3 4 631 - 445 815 86 165 13 - 12 - 35 - 44 794 - 414 584 129 256 10 - 28 - 45 - 54 786 - 31 1 1674 15 1 323 7 - 1 5

- 55 - 6 4 554 - 149 278 122 202 3 - 12 - 65 - 74 314 - 86 139 68 113 1 - 4 -

75

+ 304 - 79 77 67 104 1 - 9 - Total 10807 - 3258 8433 1252 2000 65 - 132 -47 period, the case-crossover analysis controls for l ong term trends, seasonal effects, epidemic, and other covariates that change slowly with time by design. In the bidirectional approach Ba teson and Schwartz (1999, 2001) have examined methods to sample for the referent to reduce the bias by selecting from either side of the event Also, disjoint referent periods have been recommended for bias concerns (Levy et al., 2001) . With a control close in time to the event, these methods avoid much of the confounding due to subject differences and other long-term effects such as seasonality. However, a balance must be maintain ed between a control too close in time that generates auto correlation and a control too fa r apart which confounds the longterm effects. Previous work considered alternativ es within a few weeks of the event. Neas et al. (1999) studied daily mortality in Philadelphia where the case period was the 48 hour period ending at midnight on the day of death (Schwartz and Dock ery, 1992) and the control period was the day of week seven, 14, or 21 days before and after the case period. Medina-Ramon et al. (2006) used a matching scheme from Bateson and Schwartz (1999, 2001) and a time-s tratified partition by Lumley and Levy (2000) that chose control days onl y in the same month as the hospital admission for chronic obstructive pulmonary disease (COPD). Lin et al. (2005) calculated a one to seven day exposure average ending on the ad mission date as the exposure in the case period with control periods of two weeks before and after the admissi on date. Peel et al. ( 2005) studied the hospital admissions associated with ambient air pollution le vels and respiratory health effects. The case period of three-day moving averages of pollutant concentration was selected within two weeks of the case period. The average was the average of pollutant concentr ations on the same day as the visit, one day before, and two days before. Method of the present work We chose to select the reference period based on a time-stratified scheme because it can control confounding effects from seasonality and time trends. The referents are also selected within a 28- day stratum at seven, 14, and 21 days before or af ter the case, provided they occur in the same stratum shown in Figure 7.6. Each case has three controls.

Table 7.1. Referent selection based on time-st ratified scheme with in a 28-day stratum

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 Day Stratum S M T W TH F SA S M T W TH F SA S M T W TH F SA S M T W TH F SA

1 2 3

While the case-crossover method has been widely used to examine a tran sient effect of an intermittent exposure on an infre quent acute disease, the same pr oblem can be analyzed by other methods. Szyszkowicz (2006) applie d a generalized linear mixed mode ls (GLMM) and used days of the week to construct clusters to study the effects of air pollu tion on daily emergency department (ED) visits in Vancouver, Canada. Each cluster contained four or fi ve days of the same day of the week. We conducted preliminary analys es with this approach but did not find it as effective as the case-crossover method. Mathematic al details of the case-crossove r method employed in the present work can be found in Appendix B.

Results In the analysis, we considered many potential eff ects. The first primary effect was the 24-hour average PM 10 concentration from the Inverse Distance Weighting (IDW) interpolation of the observations on the day of the asth ma event. The second primary effect, called the “lag”, was the PM 10 concentration two to six days before the event, the “lag” being the delay in days between the Event Day = Case Referent Day = Control48 exposure and the health effect. For exam ple, a two-day lag is the average PM 10 of the previous two days, abbreviated as “lag 2 ”.

Secondary effects were the covari ates of age, gender, ethnicity, and place of service. Place of service is a categorical variable with two values: emergency and inpatient. Age is encoded into four categories: 0 - 4, 5 - 9, 10 - 14, and 15 -17 years. Ethnicity followed U. S. Census Bureau definitions. Preliminary analysis detected no diff erences between places of service; therefore, this variable was not considered further. Even t hough the case-crossover de sign does control for the main effects of these patient-level covariates, we nonetheless also considered the interactions of these covariates with the event-day PM 10 concentration and its variou s lags. For example, we considered interactions such as: age*event day PM 10 , age*meanlag 2 , age*meanlag 3 , age*meanlag 4 , age*meanlag 5 , age*meanlag 6 , gender*event day PM 10 , gender *meanlag 2 , gender *meanlag 3 , gender *meanlag 4 , gender *meanlag 5 , and gender*meanlag 6
. The analysis considered all subsets of these cova riates through conditional logistic regression. The event day PM 10 , also called “Dailymean”, a nd all lag variables were highl y correlated as shown in the scatter plots of Figure 7.7 an d the statistics of Table 7.1. Thes e and the other statistics that appear in the remainder of this chapter are e xplained below, as the relationships between an increased probability of an asthma incident and PM 10 :

• Variable

• Coefficient

• Chi Square

• p-value

• Standard Error (SE)

• Adjusted Odds Ratio (Adjusted OR*)

• 95% Confidence Interval (95% CI)

• Interquartile Range (IQR)

Variable : The mean PM 10 concentration on the day of the asthma inci dent, or the average concentration of the “lag” days before the incident.

Coefficient : The slope of a line obtained usin g linear least squares fitting is called the regression coefficient, or, for short, just “Coefficient”. In Table 7.1 the first five coeffici ents are similar, indicative of49 generally similar relationships between the “Dailymean PM 10 “ and the lagged PM 10 . The last variable, “Lag6”, has a coefficien t different from the others, sugge sting a different least squares relationship. In Tables 7.2 – 7.5 the “Coefficient” has a sim ilar meaning, although it comes from the more complex “conditional logistic regre ssions” of the case-crossover method. These coefficients do vary consider ably – from 0.0016 to 0.0055 – indicating somewhat different relationships between the incr eased probability of an asth ma event and the event-day PM 10 concentrations among the four different cases.

Chi Square : Generally speaking, the chi-square test examines whether two variables are independent: in the case of Table 7.1, the two variables are:
1. the PM 10 concentration on the day of the asthma incident, and 2. the average concentration of the “lag” days before the incident.

The “Chi Square” values of Table 7.1 exceed the five percent significance level of 3.84, meaning that with an uncertainty of five percent, the eff ects of the lags differ insignificantly from the effects of the PM 10 concentration on the day of the asthma incide nt (called “Dailymean” in the table).

p-value : In statistical hypothesis testing , the p-value is the probability of obtaining a result at least as extreme as the one that was actually observed, given that the null hypothesis is true. Generally, one rejects the null hypothesis if the p-value is smaller than or equal to the significance level . If the level is 0.05, then the results are only 5% likely to be as extraordinary as just seen, given that the null hypothesis is true. In Table 7.1 the p-values are all less than the si gnificance level of 0.05; thus, the null hypothesis that the lagged concentrations differ from the event-day concen trations is rejected. In Tables 7.2 – 7.5, all p-values ar e also less than the significance level of 0.05, indicating that the PM 10 concentration on the day of the asthma inci dent was significantly associated with it.

Standard Error (SE) : The standard error of a statistical analysis is the estimated standard deviation of the error in that analysis. Specifically, it estimates the standard devi ation of the difference between the measured or estimated values and the true values. Notice that th e true value of the standard deviation is usually unknown and the use of the term standard error carries with it the idea that an estimate of this unknown quantity is being used. It al so carries with it th e idea that it measures not the standard deviation of the estimate itself but the standard de viation of the error in the estimate, and these can be very different. Standard errors of the first th ree cases in Tables 7.2 – 7.4 are similar, but the one for the fourth case – the ambulatory events excl uding preschoolers – is r oughly twice that of the first three, indicating more uncertainty for this case.

Adjusted Odds Ratio (Adjusted OR*): While its mathematical form is in Appendix E, this odds ratio can be re adily interpreted as a measure of how likely or unlikely it is for an asthma incident to be associated with a Dailymean PM 10 concentration. Odds ratios greater than 1.0 m ean that the asthma event is more likely to be associated with the PM 10 concentration; less than 1.0, less likely; and 1.0, equally likely.

95% Confidence Interval (95% CI): In statistics , a confidence interval (CI) is an interval estimate of a population parameter . Instead of estimating the parameter by a single value, an inte rval likely to include the parameter is given.50 Thus, confidence intervals, such as the 95% confidence interval fr om the case-crossover statistical analysis of the present work, indicate the reliabil ity of a estimate: the wider the interval, the lower the confidence.

Interquartile Range (IQR): The difference between the 75 th and 25 th percentiles of the daily mean PM 10 concentrations in μg/m 3

.

Figure 7.6. Correlation matrix among Dailymean PM 10 and lag variables

Table 7.2. Model statistics fo r each main effect (p reschool age group omitted)

Case-crossover: final results The case-crossover statistics are now presented a nd discussed for four different samples of the childhood asthma population:

1. all events, ages 0 -17,

2. events excluding ages 0-4, the preschool group, Variable Coefficient Chi Square P-Value Dailymean 0.00331 11.7859 0.0006 Lag2 0.00350 9.3428 0.0022 Lag3 0.00303 5.7093 0.0169 Lag4 0.00312 5.0621 0.0245 Lag5 0.00380 6.4476 0.0111 Lag6 0.00465 8.4119 0.003751 3. all ambulatory events, ages 0 – 17, and

4. ambulatory events excluding ages 0-4, the preschool group.

All events, ages 0 -17 Because of the high positive correl ations between the Dailymean and its lags, further analysis was restricted to the Dailymean as the principal covariate of PM 10 . All subsets regression indicated that the Dailymean alone was a sufficient covariate for asthma incidents. A statistically significant relationship was detected between Dailymean and asthma incidents. The adjusted odds ratio and 95% confidence interval are shown in Table 7.2

Table 7.3. Results for all asthma incidents

Variable Coefficient Standard Error p-value Adjusted OR* 95% CI

Dailymean PM 10 0.00248 0.0007507 0.0010 1.093 (1.061, 1.153)
*Adjusted OR* is the odds ratio adjusted by the approximate interquartile range of 36: (Q 3 - Q 1 ) ≈ (71-35) ≈ 36 μ g/m 3

The change in daily average PM 10 from the 25th to 75th percentile is the interquartile range (IQR). For the data analyzed here the IQR = 36 μ g/m 3 . Every additional 36 μ g/m 3 of daily mean PM 10 increases the odds ratio of an asthma incide nt for a person under 18 years old by 9. The 95% confidence interval suggests that there could be as li ttle as 6% or as high as 15% increased risk of an asthma incident.

There was also evidence of an interaction between age and Dailymean. The research group considered an analysis to compare the 0-4 age gr oup to the others because asthma is more difficult to detect in the youngest age group. Further anal ysis detected a stronger relationship between Dailymean and asthma incidents when the 0-4 age group “preschool ers”) was excluded.

Events excluding ages 0-4, the preschool group The preschool (0-4 years old) group was exclude d and the analysis was repeated. All subsets regression again indicated that Dailymean alone was a sufficient covariate for asthma incidents. The adjusted odds ratio and 95% confiden ce interval are shown in Table 7.3

Table 7.4. Results for asthma inci dents excluding the preschool group Variable Coefficient SE p- value Adjusted OR* 95% CI Dailymean PM 10 0.00331 0.0009629 0.0006 1.127 (1.053, 1.206) *Adjusted OR* is the odds ratio adjusted by approximate interquartile range of 36: (Q 3 - Q 1 ) ≈ (71-35) ≈ 36 μ g/m 3

Every additional 36 μ g/m 3 of daily mean PM 10 increases the odds ratio of an asthma incident for a person between 5 to 17 years old by 13%. The 95% confidence interval suggests that there could be as little as 5% or as high as 21% increased risk of an asthma incident.

All ambulatory events, ages 0 – 17 The ADHS data base contains emergency and inpa tient records. There are also a large number of asthma incidents by ambulatory place of service in the AZHQ database. We conducted an analysis from the AZHQ database with onl y ambulatory data to compare to the previous ADHS results. The adjusted odds ratio and 95% confiden ce interval are shown in Table 7.452 Table 7.5. Results for ambulatory events only Variable Coefficient SE p- value Adjusted OR* 95% CI Dailymean PM 10 0.00163 0.0008124 0.0448 1.060 (1.001, 1.123)
*Adjusted OR* is the odds ratio adjusted by the approximate interquartile range of 36: (Q 3 - Q 1 ) ≈ (71-35) ≈ 36 μ g/m 3

The adjusted odds ratio is 1.060. The confidence inte rval suggests that there could be as little as 0.1% or as high as 12% increased risk of an asthma incident. The results from AZHQ and ADHS showed similar effects.

Ambulatory events excluding ages 0-4, the preschool group The analysis of the AZHQ ambulatory events was re peated excluding the pres chool (0-4 years old) group. The odds ratio and 95% confidence interval are shown in Table 7.5

Table 7.6. Results for ambulatory events without the preschool group

* Adjusted OR* is the odds ratio adjusted by the approximate interquartile range of 36: (Q 3 - Q 1 ) ≈ (71-35) ≈ 36 μ g/m 3

The adjusted odds ratio is 1.219. The confidence interv al suggests that there could be as little as 8% or as high as 38% increased risk of an asthma incident. The results fr om AZHQ and ADHS data were similar. Case-crossover summary Assembled in Table 7.6, the results from the f our childhood asthma populations exhibit certain differences. That the coefficients vary considerab ly indicate four somewhat different relationships between asthma and event-day PM 10 concentration. As the standard e rrors are all similar except for the fourth case, this case is c onsidered the most uncertain by th is criterion. Among the p-values, however, it is the third case that is the most uncertain, although stil l statistically significant at the five percent confidence level. The adjusted odd ratios show that in all four cases the risk of asthma increases with an increase in PM 10 concentration of 36 μ g/m 3 . The increased risk is highest for the case of ambulatory events excl uding preschoolers and lowest for all ambulatory events, but these are the two cases of greatest uncer tainty, based on the standard error of the former and p-value of the latter. Of greater importance in this work, the adjusted odd ratios of th e first two cases -- all events with and without preschool ers – are accompanied by low standa rd errors and by p-values far less than the significance level of 0.05. The increased risk of asthma from an incremental gain of 36 μ g/m 3 of daily-averaged PM 10 is 9% for all children 0-17 year s of age and 13% for children 5-17 years of age. Variable Coefficient SE p- value Adjusted OR* 95% CI Dailymean PM 10 0.00551 0.00172 0.0013 1.219 (1.080, 1.37753 Table 7.7. Case-crossover statistics relating asthma incidents to PM 10 concentrations in central Phoenix

Case Variable Coefficient Standard Error p-value Adjusted OR* 95% CI

All events Dailymean PM 10 0.0025 0.00075 0.0010 1.093 (1.061, 1.153)

All events w/o ages 0-4 Dailymean PM 10 0.0033 0.00096 0.0006 1.127 (1.053, 1.206)

Ambulatory
events Dailymean PM 10 0.0016 0.00081 0.0448 1.060 (1.001 1.123)

Ambulatory events w/o

ages 0-4 Dailymean PM 10 0.0055 0.00172 0.0013 1.219 (1.080, 1.377)

0 5 10 15 20 25 30 35 40

1234
CaseAdjusted odds ratios (OR*) as percentages

95% CI-lower OR* 95% CI-upper case 1: all asthma events case 2: all events without preschoolers case 3: ambulatory asthma events case 4: ambulatory events without preschoolers

Figure 7.7. Adjusted odd ratios expressed as percentages with 95% lower and upper bound confidence limits for the four asthma populations

Model Validation The regression analysis is based on a linear re lationship between the Dailymean exposure and response. A model without any linear or monot onic assumptions determined the relationship between asthma incidents and Dailymean PM 10 . Divided into quartiles, the PM 10 concentrations were linked to asthma events by the indicator va riables in the regression models with the first quartile as the baseline. The natura l logarithm of the adjusted odds ra tio of an asthma incident, “log odds” of the y-axis, plotted against the exposure quartile in Figure 7.8, illu strates the monotonic, nearly linear relationship of the log odds versus the quartiles of the Dailymean PM 10
.54 Log Odds VS Grouped PM10

0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18

Q1=35 Q2=55 Q3=71 Q4=273

PM10Log Odds

Figure 7.8. Natural logarithm of the adjusted odds ratio (Log Odds) of asthma incidents versus quartiles of mean PM 10 concentrations

Chapter 8 Asthma warning systems Predicting air pollutant concentrations has been one of the mainstays of the Air Quality Division of the Arizona Department of Environmental Quality since the early 1990s. These predictions have been for ozone in Phoenix and fo r particulate matter in Phoenix an d Yuma, with more specialized predictions for various towns and cities during the forest fire and prescribed burning seasons. In general, predictions of this kind can be done in three ways, presented here in hypothetically decreasing degrees of rigor, reliability, and labor:

1. Apply deterministic grid-based meteorologica l and air quality models, coupled with a complete, model-ready emissions inventory, to predict concentrations of several different pollutants throughout a particular urban area or region. One exampl e of this is a predictive system for the eastern United States for both oz one and fine particulat es. Another example is in Seattle, Washington, wher e predictions of gaseous, particul ate, and air toxic pollutants are made throughout the Seattle-Portland region. Predictions for both of these progressive systems are distributed to the news media and are made available through the internet. The predictions are a result of first, having an adequate emissi ons inventory that is ready for the air quality model; second, runni ng a prognostic meteorological model such as MM5 to predict tomorrow’s meteorological fields of temperature, winds, and so forth; and third, running an air quality model such as the Community Multiscale Air Quality system (CMAQ) to provide the grid-based concentra tions of pollutants. Among the products of the Seattle system are color-coded concentration maps for the next day’s pollutants that can be viewed through the internet on home (and business) computers.

2. Apply a statistical model that relates historical pollutant levels to measured meteorological variables at a specific air pollutant monitori ng site to predict to morrow’s air quality. Though neither computationally no r labor intensive as the first method, this system still requires full-scale meteorologi cal modeling with a grid-based prognostic model. This meteorological modeling is first done for the historical period to derive quantitative relationships between the measured and simula ted variables such as mixing height, wind55 speed, and temperature. In cons tructing the system as a whole, the next step is to determine from the historical record the statistical relationship between measured meteorological variables and pollutant concentrations at a speci fic monitoring site. With this relationship, and with the link between the meteorologica l model’s predictive performance and the measured variables, a daily operational system can then be put into place. It works as follows and entails considerable automation. First, the prognostic me teorological model is run for the next day, producing tomorrow’s me teorological fields of importance in air pollution. Second, these predicted fields are aligned with those from the historical record and their corresponding air pollutant concentr ations. Third, tomorrow’s air pollution is predicted statistically from the predicted weather and the correspondence between the historical weather and measured air pollution. Into this system an additional component can be incorporated that adds the daily pollutant -weather measurements and simulated weather to the historical records in such a way that the system en toto effectively improves it performance with each additiona l day’s prediction. This last wrinkle is subsumed into a mathematical-statistical construc t called a “neural network.” As described in chapters 4 and 6 of this report, the historical portion of this work has already been done for PM 10 for Maricopa County’s Central Phoenix monitoring site. Based on one year of measurements, the system’s performance for November 2005 PM 10 concentrations was slightly better than the deterministic CMAQ model’s.

3. The third kind of predictive system, and the one employed by Air Quality Division meteorologists for over ten years, can be ch aracterized as “manual forecasting.” In this method the meteorologist scrutinizes diffe rent National Weather Service and other numerical models that predict a host of variables influencing pollutant concentrations. This scrutiny extends to regional a nd national weather maps and ch arts that are of help in predicting long-distance transpor t. This staff person synthesizes the various predictions, in part based on how well they have performed. The next step is to relate the most recent two or three days’ pollutant concen trations to the measured weat her variables, and, in many cases, to consult the historical records to understand how these relationships have led to specific levels of an air pollutant. Often pollu tant concentrations are assembled from broad swaths of Arizona and southern California, as well as from greater Phoenix. Because the prediction does not have to be made until noon, the staff person will usually examine the morning’s hourly concentrations to determine how quickly (or slowly) the present day’s air pollution is developing. Finally, with all this information and with genuine seat-of-the-pants knowledge and intuition, the fo recaster issues tomorrow’s predictions for ozone or particulates.

Since asthma has been clearly linked with PM 10 concentrations, not only in the present work but also in numerous studies worldwide, an effectiv e asthma warning system would consist of two parts:

• one (or more than one) of the three predictive systems just described, and

• an automated, or semi-automated communicatio n network that would not only advise the

news media on tomorrow’s PM 10 outlook but would al so send this predic tion directly to interested citizens by phone, fax, text, and/or e-ma il. These citizens would include parents of asthmatic children, school nurses and other school officials, pub lic health nurses, and, anyone of any age or occupation who wished to receiv e the notices because of respiratory sensitivity in themselves, their immediate families, relatives, friends, or neighbors.56 Chapter 9 Conclusions

Major findings on asthma and PM 10

• PM 10 concentrations have a statistically signifi cant association with asthma incidents in

central Phoenix at the 95 % confidence level.

• For children ages 5-17 an increase in the daily mean PM 10 from the 25 th to the 75 th percentile

- 36 μ g/m 3 in this study - is associated with a 13% increase in the probability of asthma incidents, an effect much stronger than in previous studies.

• The effects of age and ge nder were insignificant.

Additional findings

• CMAQ generally underestimated the highe r PM10 concentrations and overestimated the

lower, with the better correlations for the higher concentrations.

• The better correlations in the higher conc entration regimes, despite their under predictions,

still bode well for the viability of the predictive tool to warn asthma populations in central Phoenix for these reasons:

o The under predictions were most pronoun ced at the two maximum concentration

sites of Durango Complex and West 43 rd Avenue, which represent the worst PM 10 air quality in metropolitan Phoenix.

o Any predictive system that comes close to predicting these peak concentrations

provides adequate protection for asthmatic children in th e rest of greater Phoenix.

• Despite some shortcomings, the meteorol ogical model MM5 performed adequately for the

needs of this project, even though it

o failed to simulate the higher wind speeds;

o overestimated the low wind speeds during stab le conditions at nigh t and in the early

morning;

o could not capture the rapid changes of wind direction during the morning and evening transitions.

• Because the uncertainty in emissions far exceeds th at of the meteorological fields or of the

dynamics and chemistry inherent in the model, the construction, maintenance, and periodic updating of an adequate model-ready PM 10 emissions inventory ought to receive a higher priority and more resources.

• The Neural Network, which proved to be so mewhat better than CMAQ in predicting PM 10

concentrations, might be a complementary syst em for the “manual forecasting” now in place.

• Because the PM 10 concentration fields from inverse di stance weighting (IDW) and the more

sophisticated ordinary kriging were similar, onl y the former were used to produce the censustract-specific concentrations for the health analysis.57 References Arizona Department of Envir onmental Quality (ADEQ), 2008. 2007 Annual Report

. [Online] Available at: http://www.azd eq.gov/environ/air/monitoring/monitor.html [Accessed 27 October 2008]

Bateson, T. and Schwartz, J., 1999. Control for seasonal variation and time trend in casecrossover study of acute effects of environmental exposures. Epidemiology , 10(5), pp.539–544.

Bateson, T. and Schwartz, J., 2001. Select ion bias and confounding in case-crossover analyses of environmental time-series data. Epidemiology , 12(6), pp.654–661.

Byun, D. W. and Ching, J. K. S., 1999. Science Algorithms the EPA Model-3 Community Multiscale Air Quality (CMAQ) Modeling System , Washington, DC: U.S. Environmental Protection Agency, Office of Resear ch and Development, EPA-600/R-99/030.

Carolina Environmental Program, 2005. SMOKE v2.2 User’s Manual , Chapel Hill, NC

Figueiras, A., Carracedo-Martinez, E., Saez, M., and Taracido, M., 2005. Analysis of casecrossover designs using longitudina l approaches: a simulation study. Epidemiology
, 16(2), pp.239–246.

Gardner, M.W. and Dorling, S.R., 1998. Ar tificial neural netw orks (the multilayer perceptron)--a review of applicati ons in the atmospheric sciences. Atmospheric Environment , 32(14-15), pp.2627-2736

Grell, G.A. and Dudia, D.R.S., 1994. A description of the fifth generation Penn State/NCAR Mesoscale Model (MM5), Boulder, Colorado: National Center for Atmospheric Research, NCAR/TN-398+STR.

Hosmer, D. and Lemeshow, S., 2000. Applied Logistic Regression .2 nd ed. John Wiley & Sons.

Kemball-Cook, S., Jia, Y., Emery,C., Morris, R., Wang, Z., Tonnesen, G., 2004. 2002 Annual MM5 36 km simulation to support WRAP CM AQ visibility modeling for the Section 308 SIP/TIP: Preliminary Report on the In itial 2002 36 km MM5 Si mulation Performed during 2003 , ENVIRON International Corporation a nd University of California Riverside.

Lee, S. M., Giori, W., Princevac, M., Ferna ndo, H. J. S., 2005. Implementation of a stable turbulence parameterization for the mesoscal e model MM5: nocturnal flow in complex terrain, Boundary-Layer Meteorology , 119, pp.109-134.

Levy, D., Lumley, T., Sheppard, L., and Ch eckoway, H., 2001. Referent selection in casecrossover analyses of acute he alth effects of air pollution. Epidemiology , 12(2), pp.186–192.

Lin, M., Steib, D., and Chen, Y., 2005. Coarse particulate matter and hospitalization for respiratory infections in ch ildren younger than 15 years in Toronto: a case-crossover58 analysis. Pediatrics , 116(2):e235–e240.

Lumley, T. and Levy, D., 2000. Bias in the case- crossover design: imp lications for studies of air pollution. Environmetrics , 11, pp.689-704.

Maclure, M., 1991. The case-crossover design: a method for studying tr ansient effects on the risk of acute events. American Journal of Epidemiology , 133(2), pp.144–153.

Maricopa County Department of Air Quality (MCDAQ), 2007. 2007 Air Monitoring Network Review. [Online] Available at: http://www.maricopa. gov/aq/divisions/mon itoring/network.aspx [Accessed 27 October 2008]

Medina-Ramon, M., Zanobetti, A., and Schwar tz, J., 2006. The effect of ozone and PM 10 on hospital admissions for pneumonia and ch ronic obstructive pulmonary disease: a national multicity study. American Journal of Epidemiology , 163(2), pp.579–588.

Morris, R.E., Koo, B., Lau, S., Tesche, T.W., McNally, D., Loomis, C., Stella, G., Tonnesen, G., Wang, Z., 2004. VISTAS Emissions and Air Quality Modeling – Task 4cd Report: Model Performance Evaluation and M odel Sensitivity Tests for Three Phase I Episode , ENVIRON International Corporation, Novato, California.

Neas, L., Schwartz, J., and Dockery, D., 1999. A case-crossover analysis of air pollution and mortality in Philadelphia. Environmental Health Perspective , 107(8), pp.629–631.

Peel, J., Tolbert, P., Klein, M., Metzger, K ., Flanders, W., Todd, K., Mulholland, J., Ryan, P., and Frumkin, H., 2005. Ambient air polluti on and respiratory em ergency department visits. Epidemiology , 16(2), pp.164–174.

Redelmeier, D. and Tibishirani, R., 1997. Asso ciation between cellular-telephone calls and motor vehicle collisions. The New England Journal of Medicine , 336(7), pp.453–458.

Schwartz, J. and Dockery, D.W., 1992. Incr eased mortality in Philadelphia associated with daily air polluti on concentrations. American Review of Respiratory Disease , 145, pp.600–604.

Szyszkowicz, M., 2006. Use of generalized linear mixed models to examine the association between air pollution and health outcomes. International Journal of Occupational Medicine and Environmental Health , 19(4), pp.224–227.

U.S. Environmental Protection Agency, Offi ce of Air Quality Planning and Standards, 2007. Guidance on the Use of Models and Other A nalyses for Demonstrating Attainment of Air Quality Goals for Oz one, PM2.5 and Regional Haze , Research Triangle Park, North Carolina: EPA-454/B-07-002.

Western Governors’ As sociation (WGA), 2006. Final Report for the Western Regional Air Partnership (WRAP) 2002 Visibility Model Performance Evaluation , WGA 30203.59 Appendix A PM monitoring sites
Table A-1. Characteristics and images of continuous PM 10 monitoring sites Symbol Name Objective 1 Scale 2 Height 3 Cover 4 Permanent sites (all MC, except SS, which is ADEQ)
Gravel
Paved WF West 43 rd
Avenue #1 #3 Middle/ Neighborhood 5 Desert, bare dirt
Paved
Bare dirt DC Durango Complex #1 #3 Middle/ Neighborhood 3.9 Paved, bare dirt, homes, offices Gravel Paved WP West Phoenix #2 Neighborhood 3.6 Lawns/Houses Paved Paved CP Central Phoenix #2 Neighborhood 11.3 Commercial & residential Gravel Lawns SS Supersite #2 Neighborhood 4.4 Homes, apartments Paved Paved, lawns GR Greenwood #2 Micro/Middle 4.4 Freeways Paved parking lot Pavement, dirt lot SP South Phoenix #2 #3 Neighborhood 4.0 Commercial & residential Temporary sites (ASU), November 2007 – March 2008
Lawn Golf Coarse MRV Maryvale #2 Neighborhood 4.4 Homes Lawns Homes VGC Valley Garden Center #2 Neighborhood 2.0 Garden Area Paved Lawns, homes WVR Weaver’s Auto Service #2 Neighborhood 2.0 Lawns, apartments Gravel Paved CSA Community Service AZ #2 Neighborhood 4.4 Homes
1. Objectives 1. Determine highest concentrations expected to occur in the area covered by the network. 2. Determine representative concentrations in areas of high population density. 3. Determine the impact on ambient pollution levels of significant sources or source categories. 2. Choices of spatial measurement scales (radius from site) Micro Scale 0 to 100 meters Middle Scale 100 to 500 meters Neighborhood Scale 0.5 to 4 kilometers Urban Scale 4 to 50 kilometers60 Regional Scale 10 to 100s of kilometers Note: In addition to the scales specified in the table, all sites have PM 10 concentrations with some (unknown) degree of urban contributions from transport and with 10 ug/m3 of desert background on an annual basis. 3. Height of instrument inlet above ground in meters 4. Predominant ground surface in (1) immediate vicinity of the monitor (0 – 25 meters); (2) 25 – 100 meters; and (3) 100 – 500 meters.

Table A-2. Traffic near PM monitoring sites Symbol Name Streets Distance 1 / Direction Traffic 2

WF West 43 rd Ave. Broadway Road 37, S 8,000 27 th Avenue 78, E 9,000 DC Durango Complex Durango Street 325, N 9,000 39 th Avenue 30, W 2,000 Earll Drive 25, N <1,000 WP West Phoenix

Thomas Road 360, S 32,000 16 th Street 91, W 24,000 Roosevelt Street 75, N 12,000 I-10 – E-W portion 450, N 276,000 CP Central Phoenix

I-10 – N-S portion 865, E 297,000 17 th Avenue 11, E <1,000 15 th Avenue 425, E 17,000 Camelback Road 635, N 40,000 19 th Avenue 370, W 17,000 SS Supersite

I-17 1525, W 208,000 27 th Avenue 10, E 19,000 I-10 85, N 229,000 GR Greenwood

I-17 850, E 101,000 Broadway Road 360, N 25,000 Central Avenue 200, E 16,000 SP South Phoenix

7 th Avenue 680, W 22,000 51 st Avenue 115, E 30,000 MRV Maryvale Indian School Road 300, S 34,000 15 th Avenue 30, W 13,000 VGC Valley Garden Center McDowell Road 267, S 21,000 29 th Street 33, W <500 WVR Weaver’s Auto Service Thomas Road 75, S 38,000 59 t h Avenue 36, E 31,000 CSA Community Service AZ Glendale Ave./Grand Ave/59 th Ave Intersection 375, NNE 73,000

1 Distance in meters and direction towards the roadway 2 Average daily traffic for weekdays: vehicles per day61 Figure A1.1. West 43 rd Avenue PM 10 monitoring site

Figure A1.2. Greenwood PM 10 monitoring site62 Figure A1.3. Durango Complex PM 10 monitoring site

Figure A1.4. Central Phoenix PM 10 monitoring site63 Figure A1.5. West Phoenix PM 10 monitoring site

Figure A1.6. Supersite PM 10 monitoring site64 Figure A1.7. South Phoenix monitoring site65 Appendix B

Neural Network This appendix goes with Chapter 4 and explains some of the detail s of the network’s construction.

Bearing in mind that the main purpose of EnviNNet is the prediction of PM 10 , the following design choices were made:

• architecture : three-layer MLP network, w ith the number of hidden node s selected to reliably rebuild data from a test data-set;

• topology : full connection between layers and no connection between neur ons in the same layer;

• transfer function : exponential, to ensure pos itive functions at the ou tput node, and hyperbolic tangent for hidden nodes, with the latter an ex cellent compromise for a non-linear function both globally and locally;

• information flow : feed-forward;

• representation of input variables : standardization of inputs to eliminate problems due to different measurement scales for different predictors;

• reconstruction of the output signal : through linear combination and exponentialization of the outputs of hidden nodes; and

• training method : parameters are learned or estimated through the conjugate-gradient method. It is preferred over the back-propa gation method as it exploits bo th first-order (gradient) and second-order (curvature) data during op timization of the objective function.

Mathematically, the three-layer neural network has the following form: 

# ( ) ) w , x ( f y φ

=

where x represents input data, w the coefficients (paramet ers estimated by learning), f the activation functions from layers 2 to 3 and φ the activation functions from layers 1 to 2. The choice of f and
φ determines the output; for example, if φ is hyperbolic tangent and f linear, an input of meteorological and pollution variables are transforme d to an output with both negative and positive values, but if f is exponential the output can have only positive values. As such, an exponential function for the output and a hyperbolic tangent activation function for th e hidden neurons were selected.66 Appendix C MM5 - modeling results and discussion This appendix goes with Chapter 5 and explains some of the more detailed aspects of the meteorological modeling.

The modeling domain was based on a La mbert Projection, centered at (97 ° W, 40 ° N), and three nested domains with 36, 12 and 4km grids were ut ilized. In the nested simulations, the results obtained for the outer domain were used as initial and boundary conditions for the inner domain. Vertically 29 levels were applie d with the layer closes t to the ground being seve n meters to capture the surface-layer processes. The parent domain for MM5 covered the entire North American continent with 36km grids, the inner domain with 4x4km grids, was centered in Phoenix and its surrounding mountains (Figure C-1).

Figure C-1. Modeling domains used with weather predicting model MM5

The number of rows (running west-east), columns (running south-north) and the definition of the x

y origin (i.e. the southwest corner) are listed in Table C-1. For each domain the values indicated as the x

- y origin are the distances in meters of the doma in's southwest corner from the center of the parent domain (97 ° W, 40 ° N).67 Table C-1. Grid definitions for MM5 model
South-west corner (in m ) Model Columns Dot (Cross) 1 Rows Dot (Cross) 1 x Origin Y Origin MM5 – 36km 130 (129) 94 (93) -2304000 -1656000 MM5 – 12km 77 (76) 77 (76) -1836000 -972000 MM5 – 4km 83 (82) 83 (82) -1648000 -908000

1 Dot nodes are defined at grid cell vertices; cross nodes – at grid cell centers;

A model validation was necessary to estimate the model performance in the study region. Without a doubt concentration fields depend on the meteorol ogical conditions (wind, temperature, rainfall, humidity, and so forth). So it is very important to have good performance for the meteorological flow field. Several tests were made for different physical com ponents to achieve better model predictions. The sensitivity tests (three days we re performed for each case - November 2005 and April 2006) were carried out specifically to so lve temperature and humidity problems in the southwestern US. It was found in earlier studies (WGA, 2006; Morris et al., 2004; Kemball-Cook et al., 2004) that the model underestimates the amplitude of the diurnal temperature cycle and greatly overestimates the humidity during the summer. The difficulties associated with the parameterization of turbulence in the stable nocturnal planetary boundary layer (PBL) have been overcome by the use of a modified Medium Range Forecast (MRF) sc heme (Lee et al., 2005). MRF is high resolution PBL based on non-local K-theory with an additional counter gradient term that incorporates the contribution of large-scale eddies to the total fl ux. A unique aspect of this parameterization is its stability-dependant turbulent Prandtl number that allows momentum to be transported by the internal waves, while heat diffu sion is impeded by the stratificati on. The new schemes showed an improvement in predictions, particularly for nocturnal near-surface te mperatures. Surface wind predictions also improved slightly, but not to the extent of the temperature predictions. The physical options ultimately chosen are shown into the Table C-2. Table C-2. MM5 configurati on for different domains *
Analysis FDDA

Grid resolution Radiation LSM PBL Cumulus Microphysics 3D Surface D1 – 36km RRTM 5 Layers MRF mod. Betts-Miller Mixed Phase W/T/H W/T/H D2 – 12km RRTM 5 Layers MRF mod. Grell Mixed Phase W/T/H W/T/H D3 – 4km RRTM 5 Layers MRF mod. None Mixed Phase W/T/H W/T/H * Abbreviations: LSM – land surface model; PBL – planetar y boundary layer; FDDA – four dimensional data assimilation; RRTM – rapid radiative transfer model; MRF mod. – Medium Range Forecast scheme - modified; W/T/H – wind/temperature/humidity.

Numerical calculations were completed for bot h 30-day periods – November 2005 and March 11-- April 9, 2006. The model evaluati ons have been made for all th ree domains, but outcomes only for the innermost one are presented in this report. The scatter plots based on hourly wind speeds and temperature are shown for two di fferent cases with high and low pollution (Figures C-2 and C-3) with their corresponding regression coefficients. The conclusion for both cases is that we achieve excellent model performance for temperature but only average predictions for wind fields. The model was not able to capture events with high synoptic wind speeds (Figure C-4), because of the inconsistency of the input data from the operational NCEP/ETA gl obal model. The model predicted wind speeds considerably higher than the measur ed values during the night and early morning periods (Figure C-5), when stagnation with low wind speeds occurs, especially in the northeast valley. The model was able to pr edict the wind direction for a st rong synoptic wind (Figure C-4), but the disagreement between predicted and observ ed directions is quite pronounced for the case68 with low synoptic wind conditions (Figure C-5). The change in the wind direction with 12 hours periodicity is well displayed within the model predictions (Figure C-5). Mesoscale thermal circulation induced by the diurnal heating-cooli ng cycle of valley and slope flows, with downslope/down-valley flows occurring at night and up-valley/up-slope flows during the day, was well simulated by the model, especially under low synoptic wind condi tions (Figure C-6).

0 2 4 6 8 10 12 14 16 18 20 Observed wind speed [m/s] 0 2 4 6 8 10Modeledwind speed[m/s] November, 2005

R 2 =1

R 2 =0.3

-5 0 5 10 15 20 25 30 35 40 Observed temperature [ ° C] -5 0 5 10 15 20 25 30 35 40Modeled temperature [°C] November, 2005

R 2 =1

R 2 =0.88

Figure C-2. Scatter plots of wind speed and temperature for all st ations for high pollution episode November 2005

0 2 4 6 8 10 12 14 16 18 20 Observed wind speed [m/s] 0 2 4 6 8 10Modeledwind speed[m/s] March 11 - April 9, 2006

R 2 =1

R 2 =0.4

-5 0 5 10 15 20 25 30 35 40 Observed temperature [ ° C] -5 0 5 10 15 20 25 30 35 40Modeled temperature [°C] March 11 - April 9, 2006

R 2 =1

R 2 =0.84

Figure C-3. Scatter plots of wind speed and temperat ure for all stations for the low pollution period of March 11 – April 9, 200669 November 26-29, 2005

0 5 10 15

0 12 24 36 48 60 72

time [h]wind speed [m/s]

CP - obs CP - mod GR - obs GR - mod WP - obs WP - mod DC - obs DC - mod SP -obs SP - mod WF - obs WF - mod

0 90 180 270 360

0 12 24 36 48 60 72

time [h]wind direction [deg]

CP - obs CP - mod GR -obs GR - mod WP - obs WP - mod DC - obs DC - mod SP - obs SP - mod WF - obs WF - mod Figure C-4. Comparison between predicted and obser ved wind speeds and directions at different locations – with strong synoptic winds70 November 21-24, 2005

0 2 4 6

0 12 24 36 48 60 72

time [h]wind speed [m/s]

CP - obs CP - mod GR - obs GR - mod WP - obs WP - mod DC - obs DC - mod SP -obs SP - mod WF - obs WF - mod

0 90 180 270 360

0 12 24 36 48 60 72

time [h]wind direction [deg]

CP - obs CP - mod GR -obs GR - mod WP - obs WP - mod DC - obs DC - mod SP - obs SP - mod WF - obs WF - mod Figure C-5. Comparison between predicted and obs erved wind speeds and di rections at different locations – with calm synoptic conditions71 24 November, 00h

Main flow 24 November, 12h

Main flow

Figure C-6. Simulated horizontal winds at 10m a nd surface temperature fields during nocturnal (left) and daytime conditions, November 24, 2005

The Root Mean Squared Error (RMSE) is given at different locations within the study domain in Table C-3. The model gives better predictions fo r wind speeds less than 6m/s but the errors double for the higher wind speeds. The model provides wind fl ow better in the southern part (DC, SP) than the northern and central parts (C P, WP, GR) of the study domain. The worst model performance is for WF, where the influence of lo cal effects from the Salt River channel and the South Mountains appears to be beyond the model’s grasp. The MM5 m odel gives better performance for the winter period (November 2005) than the spring period (M arch-April 2006), as the RMSE for the CP and DC sites in the former are twice that of the latter. Table C-3. RMSE for wind sp eed at different sites*

Range (m/s) CP WP GR WF DC SP November, 2005 All 2.56 3.25 3.1 4.0 1.95 3.32 <6 1.6 1.52 1.56 1.5 1.33 1.34 >6 6.17 4.52 3.87 6.37 4.96 5.13 March 11 – April 9, 2006 All 5.19 5.18 - 5.21 4.28 - <6 1.9 1.76 - 1.92 1.78 - >6 6.1 4.63 - 6.06 5.72 - * Abbreviations: CP – Central Phoenix; WP – West Phoenix; GR – Greenwood; WF – West 43 rd Avenue; DC – Durango complex; SP – South Phoenix.

Additional statistics are shown in Table C-4, giving the model's perf ormance at different locations in central Phoenix for both high and low pollution periods. The following indi cators were used for performance evaluation (see Appendix 1 for details).72 Table C-4. Summary of statistical measures, in comparing observations to model predictions of wind speed at different sites *

Site MAE RMSE RSQ IA MFB [%] November 2005 CP 1.02 2.56 0.36 0.55 -25.65 WP 0.12 3.25 0.36 0.67 7.74 GR 0.01 3.09 0.37 0.68 23.17 WF 0.35 3.98 0.38 0.60 -36.36 DC 0.20 1.95 0.34 0.70 -8.81 SP 0.25 3.32 0.33 0.70 5.03 March 11- April 9, 2006 CP 2.02 2.56 0.62 0.60 -52.82 WP 1.87 3.25 0.58 0.64 -47.39 WF 1.92 3.98 0.63 0.64 -53.63 DC 1.10 1.95 0.46 0.68 -2.42 * Abbreviations: CP – Central Phoenix; WP – West Phoenix; GR – Greenwood; WF – West 43 rd Avenue; DC – Durango complex; SP – South Phoenix; MAE - Mean Absolute Error, RMSE - Root Mean Square Error, RSQ - Coefficient of Determination (aka R2), IA – In dex of Agreement, MFB - Mean Fractional Bias73 Appendix D

Definitions of statistics:

The following indicators were used for performance evaluation. Here P is the predicted value, O the measured value, and P and O the mean values.

∑ O

- P N 1

= MAE N

1 = i i i (Mean Absolute Error)

∑ ) O

- P ( N 1

= MB N

1 = i i i (Mean Bias)

N ∑ ) O

- P (

= RMSE N

1 = i i i (Root Mean Square Error) 

# [ ]

∑ ∑ ) O

- O ( ) P

- P ( ∑ ) O

- O )( P

- P (

= R N

1 = i N

1 = i 2

i 2

i 2 N

1 = i i i 2 (RSQ Coefficient of determination)

∑ ) O

- O ( ∑ ) O

- P (

= R N

1 = i 2

i N

1 = i 2

i 2 (RSQ* Coefficient of determination)

∑ |) O

- O | + | O

- P (| ∑ ) O

- P ( 1 = IA i i i

2 2 (IA Index of Agreement)

% 100

. ∑ ) O + P ( ) O

- P (
N 2

= MFB N

1 = i i i i i (Mean Fractional Bias)

% 100

. ∑ O + P O P
N 2

= MFB N

1 = i i i i i (Mean Fractional Error)

100%

. ∑ O ∑ ) O

- P (

= NMB N

1 = i i N

1 = i i i (Normalized Mean Bias)

% 100

. ∑ O ∑ O

- P

= NME N

1 = i i N

1 = i i i (Normalized Mean Error)74 Appendix E

Mathematics of the case-crossover method This appendix begins with a gene ral discussion of the “odds ratio”, the typical statistic employed in reporting the results from the case-crossover method. The appendix continues with a discussion of conditional logistic regression, an d it concludes with the specific ma thematical formulation of the odds ratio used in the present work.

The odds ratio is a measure of effect size particularly important in logistic regression

.

It is defined as the ratio of the odds of an event occurri ng in one group to the odds of it occurring in another group, or to a sample-based estimate of that ratio. These groups might be men and women, an experimental group and a control group , or any other dichotomous classification. If the probabilities of the event in each of the groups are p (first group) and q (second group), then the odds ratio is:

An odds ratio of 1 indicates that the condition or event under study is equally likely in both groups. An odds ratio greater than 1 indicates that the cond ition or event is more likely in the first group. And an odds ratio less than 1 indicates that the condi tion or event is less likely in the first group. The odds ratio must be greater than or equal to zero. As the odds of the first group approaches zero, the odds ratio approaches zero. As the odds of the second group approaches zero, the odds ratio approaches positive infinity.

For example, suppose that in a sample of 100 men, 90 have drunk wine in the previous week, while in a sample of 100 women only 20 ha ve drunk wine in the same peri od. The odds of a man drinking wine are 90 to 10, or 9:1, while the odds of a woman drinking wine are only 20 to 80, or 1:4 = 0.25:1. The odds ratio is thus 9/0. 25, or 36, showing that men are much more likely to drink wine than women. Using the above formula for the calculation yields the same result:

The above example also shows how odds ratios are sometimes sensitive in stating relative positions: in this sample men are 90/20 = 4.5 times more li kely to have drunk wine than women, but have 36 times the odds. The logarithm of the odds ratio, the difference of the logits of the probabilities
, tempers this effect, and also makes the measure symmetric with respect to the ordering of groups. For example, using natural logarithms , an odds ratio of 36/1 maps to 3.584, and an odds ratio of 1/36 maps to

− 3.584.

The increased use of logistic regression in medical and social science research means that the odds ratio is commonly used as a means of expressing the results in some forms of clinical trials , in survey research , and in epidemiology , such as in case-control studies . It is often abbreviated "OR" in reports. When data from multiple surveys is combined, it will often be expressed as "pooled OR".75 Conditional logistic regression has been widely applied to model stratified data, such as casecrossover analysis and case-control analysis (Lin et al. 2005, Neas et al. 1999, Figuerias et al. 2005, and Redelmeier and Tibishirani 1997). Logistic re gression is a form of regression used when the response is binary. Let x be a vector of k predictor variables (or covariates). Let Y be a binary response ( y = 0 or y = 1). The conditional probability that Y = 1 is P ( Y = 1| x ) = p ( x ). The logit of the logistic regression model is

p p x x f β β β + + + =

... ) ( 1 1 0 x

so the logistic regression model is

) ( ) (
1 ) (

x f x f

e e p + = x

Let consider the system that has a binary response y and a two-level covariate x , x = 0 or x = 1. The logit transformation is defined in terms of p ( x ) as 

## ( ) 

## [ ] x x p x p x f 1 0 ) ( 1 ) ( ln ) ( β β + = − =

.
The conditional logistic regression works in sim ilar fashion as regular l ogistic regression. It considers the stratification structur e in the data specifying which in dividuals belong to which strata. Suppose there are K strata with k n subjects in the th k stratum, where k = 1, 2, 3,..., K There are

k n 1 case subjects, k n 0 control subjects, and k k k n n n 0 1 + = (Hosmer and Lemeshow 2000). The logistic regression model can be shown as

x β x β
x ′ + ′ +

+ = k k

e e p k 0 0 1 ) ( β β

where is a nuisance parameter with the cont ribution of all terms constant within the th k stratum, )
,..., , ( 2 1 k β β β = ′ β are the vector of coefficients of covariates, x

.
To analyze case-crossover design using conditional logistic regr ession, the PROC PHREG from SAS® statistical software package was used. C onditional logistic regression can be performed by the PHREG procedure by using the discrete logistic model and fo rming a stratum for each matched set. The dummy survival times are needed to be cr eated so that all the cases in a matched set have the same event time value, and the correspond ing controls are censo red at later times.

From the logistic regression mode l defined earlier. The odds of y = 1 when the covariates have values x 1 is defined as p (1)/(1- p (1)) and the odds of y = 0 when the covariates have values x 0 is defined as p (0)/(1- p (0)). The odds ratio (OR) is defined as the odds for x 1 to the odds for x 0 , odds ratios (OR) can be shown to be76 

## [ ] [ ]

( ) 1 0 1 0

0 1 0 0 0 0 1 0 1 0 1 0

1 1
1 1 1
1 ) 0 ( 1 ) 0 ( ) 1 ( 1 ) 1 ( OR

β β β β β β β β β β β β β β β β

e e e e e e e e e e p p p

p

= = = ⎟ ⎠ ⎞ ⎜ ⎝ ⎛

+ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛

+ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛

+ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛

+

= − −
=

− + + + + +

The odds ratio is an associate measure of how mu ch more likely (an odds ratio > 1), unlikely (an odds ratio < 1), or equally likel y (an odds ratio = 1) it is for the response to present under x 1 than under conditions x 0
.77 GLOSSARY

ADEQ Arizona Department of Environmental Quality

ADHS Arizona Department of Health Services

ADT Average Daily Traffic (vehicles per day)

AQS Air Quality System, EPA's archiv e of ambient air pollution data

AZHQ Arizona Health Query

CEFD Center for Environmental Fluid Dynamics at ASU

CHIR Center for Health Info rmation and Research at ASU

CMAQ Community Multi-scale Ai r Quality modeling system

COPD Chronic Obstructive Pulmonary Disease

EPA US Environmental Protection Agency

Eta The global weather model of the National Weather Service

GIS Geographic Information System

HC Hospital Counts

HIPAA Health Insurance Portab ility and Accountability Act

IDW Inverse Distance Weighting (a numerical interpolation technique)

IQR Interquartile Range

LST Local Standard Time

MAE Mean Absolute Error

MM5 Mesoscale Meteorological Model, version 5

MCDAQ Maricopa County Department of Air Quality

μ g/m 3 Micrograms per cubic meter, the units of concentration for PM

MLP Multi-Layer Perceptron

NCAR National Center for Atmospheric Research

NCEP National Center for Environmental Prediction78 NN Neural Network

PM Airborne Particulate Matter

PM 10 Particulate Matter 10 microns and smaller, aka “inhalable”

PM 2.5 Particulate Matter 2.5 microns a nd smaller, aka “fine particulates”

PM 2.5 - 10 Particulate Matter 2.5 to 10 mi cros, aka “coarse particulates”

[PM 10 ] Airborne concentration of PM 10
RMSE Root Mean Square Error

SCERP The Southwest Consortium for Environmental Research and Policy

SMOKE Sparse Matrix Operator Kernal Emissions, software to build modelready air pollutant emission fields