Abstract
Background
The United States Environmental Protection Agency's Toxic Release Inventory (TRI) data are frequently used to estimate a community's exposure to pollution. However, this estimation process often uses underdeveloped geographic theory. Spatial interaction modeling provides a more realistic approach to this estimation process. This paper uses four sets of data: lung cancer ageadjusted mortality rates from the years 1990 through 2006 inclusive from the National Cancer Institute's Surveillance Epidemiology and End Results (SEER) database, TRI releases of carcinogens from 1987 to 1996, covariates associated with lung cancer, and the EPA's RiskScreening Environmental Indicators (RSEI) model.
Results
The impact of the volume of carcinogenic TRI releases on each county's lung cancer mortality rates was calculated using six spatial interaction functions (containment, buffer, power decay, exponential decay, quadratic decay, and RSEI estimates) and evaluated with four multivariate regression methods (linear, generalized linear, spatial lag, and spatial error). Akaike Information Criterion values and P values of spatial interaction terms were computed. The impacts calculated from the interaction models were also mapped. Buffer and quadratic interaction functions had the lowest AIC values (22298 and 22525 respectively), although the gains from including the spatial interaction terms were diminished with spatial error and spatial lag regression.
Conclusions
The use of different methods for estimating the spatial risk posed by pollution from TRI sites can give different results about the impact of those sites on health outcomes. The most reliable estimates did not always come from the most complex methods.
Background
Environmental pollution data such as that collected by the United States Environmental Protection Agency's Toxic Release Inventory (TRI) have been used extensively for studies in environmental justice and medical geography [1]. These studies involved estimating an individual's or a community's exposure to pollution using the spatial information contained in the TRI database. Despite the use of this spatial information, the geographical theory used to guide the estimation of locationbased exposure to pollution has frequently been limited to basic containment and buffer analysis, especially at the national scale. The aim of this research is to improve the spatial analysis of TRI data by incorporating distance decay effects derived from spatial interaction modeling in order to provide a more realistic approach to the estimation of locationbased exposure to pollution, particularly airborne pollution. This is achieved by using several different functions for calculating this exposure and comparing the results when they are used in multivariate regression analyses with lung cancer mortality rates.
The different methods for estimating the risk at a location are evaluated because, while many studies have explored and demonstrated a link between environmental pollution and a variety of adverse societal and medical effects [1], understanding the nature of this relationship is equally important. As the variety of methods used to estimate these impacts attests, the nature of this relationship is not as well understood as the existence of the relationship. The form of this relationship greatly impacts the answers to questions that may arise from the discovery of a relationship, such as the extent to which rural counties experience adverse impacts from urban polluters. A visual cartographic comparison of some approaches has been explored by McMaster et al [2], although they do not make the statistical comparison carried out here.
Prior work
Spatial analyses of toxic pollution data, whether for environmental justice or for medical geography, have typically used a simple spatial estimate of exposure. The exposure has been recorded as a binary variable (exposed or not exposed) either through spatial containment such that a person is exposed if they live in the same census tract or county as an industrial site [1,38], or a spatial buffer such that a person is exposed if they live within a threshold distance (e.g., 1 mile) of an industrial site [5,912]. Variations on the latter use multiple buffers to approximate decreasing risk with increased distance, or select a small number of neighborhoods at increasing distances which can be treated as samples from multiple buffers. This enables the study to reflect decreased exposure as the distance increases [1,9,1315]. To provide a better measure of the impact of sites on a census tract, four studies [1619] use a raster grid that can account for whether a site is in the center of the tract, or near an edge, and whether any sites are just over the border in neighboring tracts. These raster grids reflect the density of TRI sites around each raster cell, although the density is calculated using a small buffer, such as the density of sites within a one mile radius of the cell. A gradual decay of impact as distance increases is still lacking.
Accounting for the volume of the release is another important factor missed by some TRI studies [6,10,16,18]. A binary approach that considers all TRI sites equally does not allow for gradients of risk, treats exposure to one site as equivalent to exposure to many sites, and does not account for the volume and toxicity of releases at each site. The release volumes vary by orders of magnitude (figure 1). Recognizing this, many researchers do account for varying release volumes from each site [7,9,12,2027]. They often use variations of the spatial containment and buffer models described above which can incorporate the release volume (equations 1 and 2 respectively). Here, k_{ij }is the impact of site i on county j, t_{i }is the volume of releases at site i, d_{ij }is the distance between site i and county j, and T is the threshold distance. As a result of this, most studies that use these techniques to account for the release volume still reflect a simple treatment of geography by not including distance decay effects. The toxic impacts of the different chemicals on human health vary as well [3,9,12,2225,27], although this variation is not addressed in the current study.
Figure 1. Histogram of the volume of TRI releases for 1987. Histogram of the TRI release volumes measured in pounds for 1987. Note the log scale on the horizontal axis.
To address these simplifying assumptions, Dent et al [26] have proposed using a GIS to combine atmospheric modeling with the release data and health outcomes and provide a detailed analysis of the potential effects and risks associated with TRI releases. MorelloFrosch et al [23,24] and Fisher et al [21] similarly incorporate atmospheric modeling in their analysis. These models are typically used for local, rather than nationalscale analysis. The United States Environmental Protection Agency's RiskScreening Environmental Indicators (RSEI) Model uses principles of atmospheric modeling to derive a level of risk across the entire United States [28]. This has been used by Abel [9] and Downey and Hawkins [27], and is used in this resesarch.
This background discussion is summarized by table 1, which shows that while distance decay approaches have been used, e.g. [20], containment and buffers are the most common with atmospheric modeling becoming more prevalent in local studies. Exponential and powerbased distance decay approaches as found in spatial interaction modeling have, to the author's knowledge, not been used at all.
Table 1. Prior work summary
Spatial interaction modeling
In this research, I use a spatial interaction modeling approach that is more flexible than the binary approaches commonly used in spatial analysis of TRI data, yet is fast enough and generic enough to apply to the thousands or millions of release sites involved in a national scale study. Spatial interaction modeling was developed in economic geography to estimate the level of economic interaction between two towns [2932]. The underlying assumptions are analogous to the physics theory of gravity. Just as two objects in space exert a stronger gravitational pull on each other as they increase in size and move closer to each other, two towns are expected to have a stronger level of economic interaction as the towns increase in size and as the distance between them decreases. These broad trends are applicable in many fields within geography, even though the specific functional form from physics (equation 3), may not be as useful as other functional forms.
These models are used to estimate the effect of each TRI site on each county. As the toxic release volume increases, the impact of that site on the county increases. Likewise, the impact of nearby sites is assumed to be greater than the impact of more distant sites, as from Tobler's First Law of Geography [33].
There are two common distance decay functions used to model spatial interaction, which control the rate at which the impact of a site decreases with distance. The first, taken from the physics model of gravity, is the power equation, in which the impact of a site is proportional to the size of the release and inversely proportional to the distance raised to a parameterized exponent (equation 3). Here, α and θ are positive constant parameters. The location of a county is given by its centroid. Because other functional forms may be more applicable than the gravitational form, exponential decay functions (equation 4) have also been developed and used.
The models in economic geography, such as those in Sen and Smith [31], give equations with a third term which in this work corresponds to the population of the county and a related positive constant β, such that the power model becomes and the exponential model becomes k_{ij }= t_{i}^{α }p^{β}_{i }exp(θ d_{ij}). Because I use ageadjusted rates of lung cancer rather than unadjusted counts for the dependent variable, these population terms are set to 1 and effectively removed from the equations.
The only application of a distance decay function to TRI data is a comparison of toxic releases and federally assisted housing which uses a quadratic distance decay function [20] (equation 5). This is referred to here as the Cutter function after the lead author of the publication in which it was first proposed. It uses a constant parameter, θ, controlling the rate of decrease, and a threshold distance beyond which the impact is zero. The equation given here modifies equation 1 from Cutter et al [18] to incorporate the volume of the release. As in the other equations, k_{ij }is the impact of site i on county j, t_{i }is the volume of releases at site i, d_{ij }is the distance between site i and county j, and T is the threshold distance.
Figure 2 shows the effect of increasing distance on all the models except the containment model. The parameters of the models shown are 1.0 for α, 2.0 for θ, and 100 for T, with a release volume of 10,000. More complex atmospheric models, which can incorporate distance decay concepts, have been used predominantly in studies at a local scale [21,23,24,26], with only the RSEI dataset used at the national scale [27].
Figure 2. Example graph of the distance decay functions. Example graph of the four distance decay functions examined in this study: a buffer, Cutter's quadratic decay function, a powerbased decay function, and an exponential decay function. All functions use 1.0 for α, 2.0 for θ, and 100 for T, with a release volume of 10,000.
Methods
Data used
Four sets of data are used in this paper. The first are lung cancer ageadjusted mortality rates from the National Cancer Institute's Surveillance Epidemiology and End Results (SEER) database [34]. These rates are from the years 1990 through 2006 inclusive. The second are TRI releases from 1987 to 1996. The years chosen for the TRI databases reflect a lag time between chronic exposure to toxic chemicals and the development of lung cancer. All data are temporally aggregated to the entire time series, rather than evaluating yearbyyear temporal lags. The third are risk estimates computed by the EPA's RSEI program to be used as a basis for comparison against the spatial interaction estimates The RSEI data are the riskrelated results calculated from airborne releases of chemicals that are flagged as carcinogenic and have a nonzero Inhalation Unit Risk.
The final dataset, the covariates, come from multiple sources. One source is the United States Census Area Resource File [35]. Thun et al [36] show variable risks for age, sex, and racial categories, so census data for the proportion of the population which is male and the proportion of the population which is nonwhite are included. Hendryx et al [37] note that lung cancer mortality is impacted by socioeconomic factors and access to health care, so additional covariates include the percent of the population with a less than high school education, the percent of the population with a college education, the percent of families below the poverty level, the unemployment rate, and the number of physicians per 1000 residents. Because smoking is the most significant risk factor for lung cancer [36], I also include the smoking rate of the county based on the BRFSS survey data from 2003 to 2006. Different regions of the United States have different rates of lung cancer [38], so the covariate data also includes spatial indicator variables recording whether a county is in the American south, northeast, Midwest, or western region, and whether a county is part of Appalachia, a regional designation from the Appalachian Regional commission which overlaps parts of the northeast, Midwest and southern regions. Due to a lack of data, information regarding personal movements is not included, although analysis comparing place of birth with place of death may partially account for this [39].
The Modifiable Areal Unit Problem [40] introduces difficulties into the interpretation at the county scale, especially in the larger counties in the western United States where the county centroid may be tens of miles away from the county's population center. Additionally, in these larger counties, the risk may vary within the county, and this variance is masked by calculating the risk at the county scale. However, some of the covariate data (e.g., the BRFSSderived smoking rate) is not available at a finer scale, necessitating a countylevel analysis. In the research presented here, the impacts of the Modifiable Areal Unit Problem and large county sizes are expected to have a similar impact across all models because all tests use the same spatial scale.
An examination using a synthetic dataset was considered, but the results of such a test would minimize the AIC in the situation reflecting the way the dataset is constructed (e.g., the impact falls off according to an exponential distance decay function), which may or may not correspond to a real world situation. Therefore, actual, rather than synthetic, data are used in this research.
Methods Applied
Three sets of releases from 1987 to 1996 in the TRI database are used. The first is all releases flagged as carcinogenic. The second is all releases of chemicals identified as inducing lung cancer. These chemicals are those from a parallel study [41] plus beryllium and lead, which were identified as related to lung cancer by the lead author of [41] in a private communication. The total list of chemicals is arsenic, beryllium, 1,3butadiene, cadmium, chromium, formaldehyde, lead and nickel. The third set of releases adds to the second set those releases identified as generic compound categories of elements in the first set. An example of this is a release of "arsenic compounds" in addition to releases of plain arsenic. The impacts of these three sets of releases on all counties in the contiguous United States were calculated using the containment, buffer, power, exponential and Cutter models given above. These release impacts are summed to create the cumulative impact on a county (equation 6), where k_{ij }is the impact of site i on county j and K_{j }is the cumulative impact on county j. Because the release amounts vary by several orders of magnitude and have an approximately lognormal distribution, as shown in figure 1, both the log_{10}transformed and the untransformed release volumes were tested.
The calculated impacts and covariates are then used in multivariate regression models calculated with the R software package [42]: ordinary least squares (OLS) linear regression, general linear model regression (GLM), spatial lag regression, and spatial error regression. The latter two incorporate spatial dependence in the regression model and are detailed below. There is spatial autocorrelation in the response variable (Moran's I = 0.69, p < 0.01), demonstrating the existence of spatial dependence and suggesting the applicability of the spatial regression techniques. Geographically Weighted Regression [40], which can vary the regression coefficients across the study area was not applied because it is unlikely that the nature of the relationship between toxic releases and mortality changes across the country.
This situation is not strictly one of evaluating a single function against a null hypothesis of zero impact from toxic emissions, but is rather evaluating many different functions against both each other and the null hypothesis. This makes the task more akin to model parameter optimization than traditional statistical hypothesis testing. It is considered here that, of the different functions and their parameterizations, the most appropriate representation is the one that minimizes the Akaike Information Criterion (AIC) of a regression test in which the modeled risk is one of the independent variables. In all the regression tests used in this comparison, the remainder of the independent variables are the demographic, behavioral and regional covariates, and the dependent variable is lung cancer mortality.
Experiments testing many parameterizations of the buffer, Cutter, power, and exponential functions were used to guide the results given here. These parameterizations are for the contiguous United States, and were evaluated on which parameterizations gave the lowest AIC values when combined with the covariate data in an OLS regression using the lung cancer mortality rate as the dependent variable. The same tests are conducted for the containment and RSEI approaches to risk estimation for comparison. AIC values were also computed for generalized linear model regressions, although none of the generalized linear model regressions produced lower AIC values than linear regression. As a result, the generalized linear model regressions are not further discussed. Preliminary experiments (not presented) demonstrated that the parameterizations that perform well for OLS regression also typically perform well for the spatial lag and spatial error regressions. The tests presented here evaluate α and θ values of 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, and 5.0, and distance thresholds of 5, 10, 15, ... 500 miles (8, 16, 24, ... 805 km) for the buffer and Cutter functions. After the parameterizations were found, spatial lag regression and spatial error regression models were computed. Also, linear regression models for each of the ruralurban continuum codes from the United States Department of Agriculture. Lastly, maps of proportion of the TRI impact for each county that originated in release sites located in urban areas were produced. This allows an examination of whether the impact of sites in urban areas is limited to those cities, or whether it extends far into the surrounding rural areas, and demonstrates that some environmental justice questions are not robust to the choice of risk model.
Spatial Regression Methods
The two spatial regression methods, spatial lag and spatial error regression, both account for the spatial autocorrelation that is almost always present in geographic data by adding a term to the regression equation. This spatial autocorrelation can be the result of diffusion effects of the dependent variable, which is unlikely in this situation, or the result of risk factors which have not been accounted for elsewhere in the regression model inducing spatial autocorrelation of the dependent variable [43]. The standard OLS regression (equation 7) estimates the dependent variable, which is the lung cancer mortality rate, as a linear combination of the independent variables, here the TRI interaction term and the demographic covariates. The county is j, the dependent variable at county j is Y_{j}, the independent variables at county j are X_{j,a}, ε_{j }is the error term, and β_{0 }... β_{n }are the regression coefficients.
Spatial lag regression incorporates the autocorrelation directly into the model by including a term where the dependent variable at county j is dependent not only on the independent variables at county j, but also on the dependent variable values of county j's neighbors [44]. The neighbors are defined by a weights matrix typically using one of the following three options: (1) all counties which share a border with county j as its neighbors, (2) all counties within a threshold distance of county j as its neighbors, or (3) the nearest neighbors of county j. Here, option (2) is used a distance threshold of 92 miles (148 km), which is the minimum distance that ensures all counties have at least one neighbor. Option (1) is also used considering shared corners (eg, the Four Corners meeting between Arizona, New Mexico, Colorado and Utah) as neighbors, which is called "queen contiguity" to assess the sensitivity of the results to the choice of weights matrix. Thus spatial lag regression uses the dependent variable values of the neighbor counties to calculate a new independent variable for county j, giving equation 8, where ρ is a coefficient describing the strength of the spatial autocorrelation, w_{j,k }is the spatial weight between counties j and k (typically 1 for neighbors and 0 for nonneighbors, but a distance decay form for the weights is possible), and N_{j }is the neighborhood of counties around county j. The ρ coefficient can be estimated in the same way that the β coefficients are estimated. A computationally efficient approach is given in Smirnov and Anselin [45].
Spatial error regression (equation 9) works similarly to spatial lag regression, except the autocorrelation term applies to the error terms of the neighboring counties rather than their dependent variable values [46]. Because of the circular dependence of the error terms (ie, if county j is in county k's neighborhood and vice versa, the value of ε_{j }is affected by the value of ε_{k }while the value of ε_{k }is affected by the value of ε_{j}), standard estimation techniques will not work. An estimation procedure for this is also given by Smirnov and Anselin [45].
Results
Table 2 presents the parameterizations that gave the lowest AIC values and thus are used for further analysis. The results shown use the lung carcinogens with compounds dataset. All three TRI release sets described above (all carcinogens, lung carcinogens, lung carcinogens with compounds) were evaluated as were the log_{10 }transformed release values, and the lung carcinogens and related compounds gave the lowest overall AIC values. While the containment approach was best fit with the logtransformed releases of lung carcinogens and related compounds, and the exponential and power functions had the best fits with releases of all carcinogens, the improvements were minimal; the differences in AIC are less than 2.0 for containment and the exponential function and approximately 10.0 for the power function. Therefore, to ensure consistency in the later tables, the untransformed releases of lung compounds with carcinogens are used.
Table 2. Model parameterizations
Table 3 shows the Rsquared values, the Akaike Information Criterion, and the probabilities that the spatial interaction terms are nonzero. For each regression model (OLS, spatial lag, or spatial error), the bestperforming distance decay function is highlighted in bold. In all cases, this was the buffer model. Table S1 in Additional file 1 shows the equivalent table for the queen contiguity weights matrix. The choice of weights matrix did not alter the results for most decay functions, only substantially increasing the AIC of the buffer model, but not enough to make another decay function better. As such, the change in weights matrix did not alter conclusions about which decay function performed best. Table 4 gives the full regression results for the overall best parameterization: the buffer model at 500 miles (804 km). This table also gives values of each independent variable's variance inflation factor (VIF). Since all the VIFs are less than 10, collinearity is not problematic in this model. As some of the covariates did not have significant coefficients in the bestperforming model, the least significant covariate was iteratively removed from the model until all independent variables were significant, producing the model in the right side of table 4. Similarly, nonlinear functions of each of the independent variables were also applied to each of the six parameterizations in Table 3, following [47]. While the fits are improved (minimized AIC = 22126.54 with the buffer function), the more complex regression models do not alter the conclusions about which spatial interaction models perform well and which perform poorly. Table S2 in Additional file 2 gives the best performing model results. Table 5 gives the Rsquared values of the OLS regressions by urbanrural code. The values of these codes are in table 6.
Table 3. Regression results
Additional file 1. Table S1  Results for the multivariate ordinary least squares, spatial lag, and spatial error regressions of ageadjusted lung cancer mortality versus covariates and risk estimates from releases of lung carcinogens and related compounds calculated with each spatial interaction model. Spatial regressions here use queen contiguity matrix to determine whether two counties are neighbors. Bold entries indicate which spatial interaction model performed best. Note that lower values for the Akaike Information Criterion are preferred.
Format: XLS Size: 20KB Download file
This file can be viewed with: Microsoft Excel Viewer
Table 4. OLS Regression results
Additional file 2. Table S2  Results for the multivariate ordinary least squares regression of ageadjusted lung cancer mortality versus nonlinear functions of both covariates and risk estimates from lung carcinogens and related compounds calculated with the buffer model. This is shown as it minimizes the AIC across all decay functions (Table 3). This is equivalent to Table 4 in the main document, but includes natural logarithms (e.g., log(pov)), squared values (e.g., pov2) and cubed values (e.g., pov3).
Format: XLS Size: 28KB Download file
This file can be viewed with: Microsoft Excel Viewer
Table 5. OLS regression Rsquared values by ruralurban code
Table 6. Definition of each ruralurban code
Maps of the percent of TRI impact in each county that is due to source locations in urban counties according to each of the distance decay functions are given in figure 3. The buffer, Cutter, exponential and power functions using the parameterizations in table 2 are shown. The darker counties have a greater percent of their impact from releases in urban counties, whether the total impact is high or low.
Figure 3. Percent of TRI impact from urban releases. Percent of TRI impact from urban releases using the different decay functions. Darker counties have a higher percentage of their impact coming from release sites in urban counties, while lighter counties have a higher percentage of their impact coming from release sites in suburban and rural counties.
Discussion
The buffer and Cutter functions outperform the containment, power and exponential functions (Table 3). This improvement is notable both for the Rsquared values of the OLS regressions and the Akaike Information Criterion (AIC) for all regressions. Because the AIC penalizes regression models with more parameters, lower values are preferred. These functions also outperformed the RSEI riskrelated results values, which was not the expected outcome. Other RSEI products, the hazard, modeled hazard, and modeled hazard*population product were also tested, but all performed worse than the buffer and Cutter functions. However, the power and exponential functions from the spatial interaction literature did not perform much better than leaving out the toxicity term, and occasionally even increased the AIC value, which may result from the coarse resolution of the countylevel dataset. These two functions may yet be useful at a finer scale.
The improvement of the buffer and Cutter functions over the RSEI data demonstrates that despite the difficulties posed by the Modifiable Areal Unit Problem and the size of large western counties obscuring variation of risk within the county, these spatial interaction approaches may still be an accurate reflection of the risks posed by TRI facilities. It should be cautioned that while this work demonstrates a relationship between TRI facilities and lung cancer, it does not yet indicate a causal link, nor does it indicate that the bestfitting risk estimation method, a large buffer around the TRI site, has the strongest causal relationship with lung cancer mortality.
Additionally, AIC values are better for the spatial regression techniques compared to the OLS regression values. However, including the TRI term in the spatial regression techniques does not lead to as much improvement over the base case of no interaction term. Even so, the improvement in the buffer model gives the spatial error regression of the buffer model the lowest AIC value. Both the Moran's I spatial autocorrelation statistic given above and the lower AIC values for the spatial regression techniques indicate that there is spatial dependence in lung cancer mortality. Moreover, this spatial dependence is not accounted for by the independent variables. This dependence is most likely the result of one or more additional spatial processes affecting lung cancer that are not accounted for in these data, rather than a simple diffusion or contagion process of lung cancer itself. The limited improvement from adding the TRI impacts strengthens the suggestion that there remain geographic processes affecting lung cancer that are not accounted for in these datasets. While it is not executed in this study, GWR may also reveal further evidence of confounding processes by revealing interactions with modeled covariates via nonstationary regression coefficients.
As with the different regression methods, the buffer and Cutter functions have the best Rsquared values across the entire range of ruralurban continuum codes (Table 5). Also, category 5, defined as counties containing a larger town (more than 20,000 residents) but which are not adjacent to a metropolitan area, has much higher Rsquared values than the other ruralurban codes across all models. It is not yet clear why this would be the case.
These results suggest that changing the method used to estimate risk will change the representation of the spatial impacts of the TRI sites on public health. As others have noted, the scope and scale of analysis can substantially impact the results [48], so researchers should be cautious when generalizing these findings at a county scale and national scope to more local scales and scopes. Nonetheless, researchers using the TRI dataset to estimate the health risks from pollution should carefully consider the method used to estimate the risk, as the most sophisticated model used here, the RSEI data, did not provide the lowest AIC values.
The maps in figure 3 display the percent of the TRI impact on each county from sources in urban areas calculated using the functions that performed best in the earlier results. As estimated by these models, the potential effects of pollution from urban TRI releases extend far beyond the limits of the urban areas. However, the extent varies depending on which function is used and how it is parameterized, highlighting the importance of using an appropriate function. In the power and exponential maps (figure 3a and 3b), the impacts from urban release sites are more limited to urban areas and the nearby rural communities. In both the buffer and Cutter maps (figure 3c and 3d), rural areas in the northeastern and southwestern United States, have between 75 and 100% of their estimated TRI impact from release sites in urban areas. These extended effects of urban areas are related to the large radii used in the distance decay functions. Additional work is needed to examine the environmental toxicology to determine whether the chemicals being released could travel such large distances or whether these models are simply capturing spatial dependence of the outcome that is induced by a confounding spatial process.
Future work will investigate the parameterization choices of the functions. This ad hoc approach to parameterizationexamining different possibilities of the α, θ and threshold parametersis not the ideal approach. Statistical approaches to finding the optimal α and θ parameters can be incorporated to improve the spatial interaction models that are generated [29,30,32]. A geostatistical approach can be applied to determine the decay function form and parameters. A correlogram plot comparing the distance between two counties and the difference between their mortality rates or their residuals from a regression function could be used to parameterize the function. Additionally, subsets of the correlogram could be examined separately to investigate anisotropy and nonstationarity. However, with both the ad hoc approach in this paper and a statistical modelfitting approach, using the data to optimize the parameters and then using those parameters to analyze the same data introduces circularity into the modelfitting process that would best be avoided.
A more theoretically sound approach would be to vary the α and θ parameters based on the properties of the toxic chemicals that are released. Varying α is similar to methods used somewhat frequently to account for the different toxicity of the chemicals released [2,8,11,2022], although the studies cited here use multiplicative rather than exponential modifiers (α t_{i }instead of t_{i}^{α}). In each case, higher values of α correspond to more toxic chemicals. Different studies have made this adjustment using different references, including American Conference of Governmental Industrial Hygienists Threshold Limit Values [3,22], a chronic toxicity index [12], an inhalation unit risk [23], a lifetime cancer risk [24], and the RSEI model [9,27]. Similarly, θ and T can be varied to reflect differences in airborne transport of the chemicals. If a chemical travels more easily and farther, lower values of θ and higher values of T can be used. These parameters can also be varied based on the direction from the release site to the affected community, thus incorporating anisotropy.
Ongoing work includes the refinement of atrisk population estimates using the LandScan USA population dataset [49] which can explore variation missed by countylevel populations unable to capture finescale risks. For example, if a chemical is only present in the atmosphere within a mile of the release site, any countybycounty analysis will be problematic because the spatial resolution of countylevel data are coarser than a square mile. The LandScan dataset provides population estimates at a 3 arcsecond resolution (roughly 90 meters). This can then provide improved estimates of the number of people within one mile of the release site instead of assigning the impact of a release site on the county as if everybody lived at the centroid of the county. This approach will have stronger effects on the power and exponential models because they have more rapid decreases in the impact as one travels farther from the release site (figure 2). This ongoing work also incorporates the adjustments given above varying the parameters to account for properties of the chemicals released and local climatic conditions to account for prevailing wind directions.
Conclusions
The research in this paper demonstrates that the use of simple containment techniques for estimating the spatial risk posed by pollution from TRI sites as well as the RSEI riskrelated results can give misleading results about the impact of those sites on health outcomes. This is done through a comparison of multivariate regression results using inputs of six different functions for estimating the impact of a release site on a county: containment, buffering, the quadratic distance decay function proposed by Cutter et al [20], an inverse power distance decay function, an exponential distance decay function, and the RSEI riskrelated results. The buffer and Cutter approaches consistently performed the best among these methods. The effects of this function choice are also demonstrated through mapping the percent of the overall impact that comes from urban TRI sites for all models except containment. As refinements to the parameterization process are made, the utility of more theoretically sound spatial interaction models will improve further.
Competing interests
The author declares that he has no competing interests.
Acknowledgements
Support for this report was provided by the Office of Rural Health Policy, Health Resources and Services Administration, PHS Grant No. 1 U1CRH106640100. The author also acknowledges Dr. Michael Hendryx for support and comments on an earlier draft of this manuscript. Lastly, the author thanks two anonymous referees for their valuable comments and suggestions.
References

Pastor M Jr, Sadd JL, MorelloFrosch R: Waiting to Inhale: The Demographics of Toxic Air Release Facilities in 21^{st}Century California.
Soc Sci Quart 2004, 85:420440. Publisher Full Text

McMaster RB, Leitner H, Sheppard E: GISbased Environmental Equity and Risk Assessment: Methodological Problems and Prospects.
Cart and Geog Info Sys 1997, 24:172189. Publisher Full Text

Bowen WM, Salling MJ, Haynes KE, Cyran EJ: Toward Environmental Justice: Spatial Equity in Ohio and Cleveland.
Annals of the Assn of Am Geog 1995, 85:641663. Publisher Full Text

Cohen MJ: The spatial distribution of toxic chemical emissions: Implications for nonmetropolitan areas.

Croen LA, Shaw GM, Sanbonmatsu L, Selvin S, Buffler PA: Maternal Residential Proximity to Hazardous Waste Sites and Risk for Selected Congenital Malformations.
Epidemiology 1997, 8:347354. PubMed Abstract  Publisher Full Text

Cutter SL, Solecki WD: Setting environmental justice in space and place: Acute and chronic airborne toxic releases in the southeastern United States.
Urb Geog 1996, 17:380399. Publisher Full Text

Daniels G, Friedman S: Spatial Inequality and the Distribution of Industrial Toxic Releases: Evidence from the 1990 TRI.

Shaw GM, Schulman J, Frisch JD, Cummins SK, Harris JA: Congenital Malformations and Birthweight in Areas with Potential Environmental Contamination.
Arch of Env Health 1992, 47:147154. Publisher Full Text

Abel TD: Skewed Riskscapes and Environmental Injustice: A Case Study of Metropolitan St. Louis.
Env Mgmt 2008, 42:232248. Publisher Full Text

Gragg RD III, Christaldi RA, Leong S, Cooper M: The location and community demographics of targeted environmental hazardous sites in Florida.

Kearney G, Kiros GE: A spatial evaluation of socio demographics surrounding National Priorities List sites in Florida using a distancebased approach.
Int J Health Geog 2009, 8:33. BioMed Central Full Text

Neumann CM, Forman DL, Rothlein JE: Hazard Screening of Chemical Releases and Environmental Equity Analysis of Populations Proximate to Toxic Release Inventory Facilities in Oregon.
Env Health Persp 1998, 106:217226. Publisher Full Text

Berry M, Bove F: Birth Weight Reduction Associated with Residence near a Hazardous Waste Landfill.
Env Health Persp 1997, 105:856861. Publisher Full Text

Knox EG, Gilman EA: Hazard proximities of childhood cancers in Great Britain from 195380.
J Epid Comm Health 1997, 51:151159. Publisher Full Text

Nordström S, Beckman L, Nordenson I: Occupational and environmental risks in and around a smelter in northern Sweden: I. Variations in birth weight.
Hereditas 1978, 88:4346. PubMed Abstract

Downey L: Spatial Measurement, Geography, and Urban Racial Inequality.
Social Forces 2003, 81:937952. Publisher Full Text

Mennis J: Using Geographic Information Systems to Create and Analyze Statistical Surfaces of Population and Risk for Environmental Justice Analysis.
Soc Sci Quart 2002, 83:281297. Publisher Full Text

Mohai P, Saha R: Reassessing racial and socioeconomic disparities in environmental justice research.
Demography 2006, 43:383399. PubMed Abstract  Publisher Full Text

Mennis JL, Jordan L: The Distribution of Environmental Equity: Exploring Spatial Nonstationarity in Multivariate Models of Air Toxic Releases.
Annals of the Assn of Am Geog 2005, 95:249268. Publisher Full Text

Cutter SL, Hodgson ME, Dow K: Subsidized Inequities: The Spatial Patterning of Environmental Risks and Federally Assisted Housing.
Urb Geog 2001, 22:2953. Publisher Full Text

Fisher JB, Kelly M, Romm J: Scales of environmental justice: Combining GIS and spatial analysis for air toxics in West Oakland, California.

Horvath A, Hendrickson CT, Lave LB, McMichael FC, Wu TS: Toxic Emissions Indices for Green Design and Inventory.

MorelloFrosch R, Pastor M Jr, Sadd J: Environmental justice and Southern California's 'riskscape'.
Urb Affairs Rev 2001, 36:551578. Publisher Full Text

MorelloFrosch R, Pastor M Jr, Porras C, Sadd J: Environmental Justice and Regional Inequality in Southern California: Implications for Future Research.

Pastor M Jr, MorelloFrosch R, Sadd JL: Breathless: Schools, Air Toxics, and Environmental Justice in California.
The Policy Studies J 2006, 34:337362. Publisher Full Text

Dent AL, Fowler DA, Kaplan BM, Zarus GM, Henriques WD: Using GIS to Study the Health Impact of Air Emissions.
Drug and Chem Toxic 2000, 23:161178. Publisher Full Text

Downey L, Hawkins B: SingleMother Families and Air Pollution: A National Study.
Soc Sci Quart 2008, 89:523536. Publisher Full Text

United States Environmental Protection Agency: RiskScreening Environmental Indicators (RSEI) Model. [http://www.epa.gov/opptintr/rsei/] webcite

Batty M, Mackie S: The calibration of gravity, entropy, and related models of spatial interaction.
Env and Planning 1972, 4:205233. Publisher Full Text

Baxter M: Estimation and inference in spatial interaction models.

Sen A, Smith TE: Gravity Models of Spatial Interaction Behavior. New York, NY: Springer; 1995.

Sheppard ES: The distancedecay gravity model debate. In Spatial Statistics and Models. Edited by Gaile GL, Willmott CJ. Boston, MA: D. Reidel Publishing Company; 1984:367388.

Tobler WR: A computer movie simulating urban growth in the Detroit region.
Econ Geog 1970, 46:234240. Publisher Full Text

SEER. Surveillance, Epidemiology, and End Results (SEER) Program: SEER*Stat Database: Incidence  SEER 13 Regs LimitedUse. Nov 2009 Sub (19922007) <Katrina/Rita Population Adjustment>. [http://www.seer.cancer.gov] webcite
National Cancer Institute, Cancer Statistics Branch, released April 2010, based on the November 2009 submission;

ARF: Area Resource File. Rockville, MD: 2006 U.S. Department of Health and Human Services, Health Resources and Services Administration, Bureau of Health Professions; 2005.

Thun MJ, Henley SJ, Burns D, Jemal A, Shanks TG, Calle EE: Lung Cancer Death Rates in Lifelong Nonsmokers.
J of Natl Cancer Inst 2006, 98:691699. Publisher Full Text

Blot WJ, Fraumeni JF Jr: Geographic Patterns of Lung Cancer: Industrial Correlations.

Hendryx M, O'Donnell K, Horn K: Lung Cancer Mortality Is Elevated in CoalMining Areas of Appalachia.
Lung Cancer 2008, 62:17. PubMed Abstract  Publisher Full Text

Sabel CE, Boyle PJ, Löytönen M, Gatrell AC, Jokelainen M, Flowerdew R, Maasilta P: Spatial Clustering of Amyotrophic Lateral Sclerosis in Finland at Place of Birth and Place of Death.
Am J Epidem 2003, 157:898905. Publisher Full Text

O'Sullivan D, Unwin DJ: Geographic Information Analysis. 2nd edition. Hoboken, NJ, USA: John Wiley & Sons, Inc; 2010.

Luo J, Hendryx M, Ducataman A: Association between Six Environmental Chemicals and Lung Cancer Incidence in the United States.
J Rural Health
in review

The R Project: The R Project for Statistical Computing. [http://www.rproject.org/] webcite

Ward MD, Gleditsch KS: Spatial Regression Models. Los Angeles, CA, USA: Sage Publications; 2008.

Robert Haining: Spatial Data Analysis: Theory and Practice. Cambridge, UK: Cambridge University Press; 2003.

Smirnov O, Anselin L: Fast maximum likelihood estimation of very large spatial autoregressive models: a characteristic polynomial approach.

Rogerson PA: Statistical Methods for Geography: A Student's Guide. 2nd edition. Los Angeles, CA, USA: Sage Publications; 2006.

Austin MP, Belbin L, Meyers JA, Doherty MD, Luoto M: Evaluation of statistical models used for predicting plant species distributions: Role of artificial data and theory.
Ecol Model 2006, 199:197216. Publisher Full Text

Baden BM, Noonan DS, Turaga RMR: Scales of Justice: Is there a Geographic Bias in Environmental Equity Analysis?
J Env Plan Mgmt 2007, 50:163185. Publisher Full Text

Bhaduri B, Bright E, Coleman P, Urban M: LandScan USA: A High Resolution Geospatial and Temporal Modeling Approach for Population Distribution and Dynamics.
GeoJournal 2007, 69:103117. Publisher Full Text