This evaluation comprises of two analytical steps. The first step is the simulation of PM2.5 concentrations originated from each developed shale gas well over the period of 2005 to 2017 in the Marcellus shale region of Pennsylvania. We identified the locations where annual mean concentration of PM2.5 changed due to shale gas development and estimated the number of people who could have experienced this exceedance. In the second step, we calculated the health risk associated to the simulated concentrations and determined the potential number of health problem cases, based on the results from epidemiological studies.
2.1. Model
The main model, implemented MATLAB (R2015a), is a health risk function which simulates the health risk due to exposure to PM2.5 emissions. For the purpose of this study, concentration-response function is a better method than dose-response function due to uncertainty in actual internal dose which serves as the input for the latter function. Concentration-response function (CRF) is as formulated below (Evans et al., 2013):
where is relative risk at exposure compared to the reference exposure , and represents the pollutant effect coefficient which indicates association between chronic emission exposure and cause of mortality or a disease. Coefficient is retrievable from epidemiological studies.
Results from CRF will be the input of the health impact function. In this study, we used a typical log-linear health impact function, introduced by Fann et al. (2012) as follows:
where indicates resulting change in the number of adverse health outcomes, is the baseline incidence rate associated to reference exposure , is the population affected by the change in air quality.
In this evaluation, we estimate exposure and the population affected by the change in air quality at each specific distance from every developed well between 2005–2017. For this purpose, we applied the Gaussian plume model developed by Banan and Gernand (2020) which simulates the concentration of non-reactive air pollutants such as primary PM2.5 within the vicinity of their source and sync the results with census data to provide an estimate of people experienced exposure to certain level of such pollutant concentrations. For the purpose of this analysis, we modified Banan and Gernand (2020) model to identify any location beyond PA’s setback from each developed well where the change in PM2.5 concentrations due to emissions from shale gas development was estimated to be positive. More details on the modified Gaussian plume model used in this evaluation is available in the Additional Text File.
2.2. Data Sources
This study focuses on estimation of changes in the incidences of health impacts associated with residents’ exposure to PM2.5 emissions from shale gas developments between 2005 and 2017 in Pennsylvania. Thus, we simulate such impacts in association with any developed well during this period of time in Marcellus shale region of Pennsylvania.
Well data comprising of location, number of wells per wellpad, and “SPUD Date” (the date that drilling operation started) was retrieved from the reports by Pennsylvania Department of Environmental Protection (Department of Oil and Gas Reporting website, 2019). We used the wind data provided by Iowa Environmental Mesonet (IEM, 2018), measured at seventeen weather monitoring stations in Pennsylvania between 2005 and 2017, comprising of hourly wind direction, wind speed, relative humidity, and cloud coverage.
Solid information is not available on the rate of pollutant emissions from shale gas development stages. We repeated the Monte Carlo simulation by Roy et al. (2014), accounting for drilling rate of 1000 ft/day and fracking rate of 1000 ft/day in three stages estimated in literature and technical reports (McKeon, 2011; Facts about Canada’s Oil and Natural Gas Industry, 2019; Coloradans for Responsible Energy Development, 2019). The final estimated mean and 95% CI for PM2.5 emission rates used in this evaluation were estimated to be 0.45 (0.1–1.3) and 1.06 (0.28–1.88) gr/hr from drilling and hydraulic fracturing, respectively (refer to Additional Text File for more details). In this manuscript, the PM2.5 concentration simulations and health risk evaluations are mainly associated to the high (97.5th percentile) emission rate value for either operation phase (unless otherwise noted). Detailed results for mean and low emission rates are available in Additional Excel File.
We used the latest U.S. Census block data (Census Data, 2010) to investigate potential number of people who might have experienced air quality changes due to emission from shale gas development. Corresponding data was population, block area, and block geographic location (latitude and longitude).
Pollutant coefficient (, Eq. 1) and baseline incidence rate (, Eq. 2) for different health outcomes are the two main inputs of our model. Epidemiologic studies have provided estimated pollutant coefficients for the health issues which evaluations show associations between exposure to PM2.5 emissions and increase in their incidence. Table-1 provides the list of pollutant coefficients and baseline incidence rates used in this evaluation. The U.S. EPA developed the environmental Benefits Mapping and Analysis Program (BenMAP) which has been used as an analytical tool by scientists, policy analysts, and decision makers for air quality management purposes and policy assessment and to estimate the human health impact associated with changes in air quality (RTI International, 2015). Therefore, we used the same reference studies included by BenMAP in order to keep the same ground for our analytical estimates as the ones used in policy and regulatory purposes before.
[Table-1: Pollutant coefficients and baseline incidence rates for health impacts, attributable to hazardous level of PM 2.5 emissions]
Currently, the health effects associated with exposure to PM2.5 emissions may not be attributed to specific source or individual components of PM2.5 emissions (EPA, 2010; Burnett et al., 2014; WHO, 2007). This study estimated the health impacts based on number of available pollutant effect coefficients for any source or component of PM2.5 emissions (listed in Table-1).
2.3. Analysis
This study provides an estimation of the health impacts imposed by shale gas development as a contribution to similar efforts by studies such as Litovitz et al. (2013) and Sovacool (2014). The output of our model is the estimated change in the number of health impact incidences due to exposure to PM2.5 emissions from shale gas development activities by county and year.
Assuming perfect compliance with the setback policy in Pennsylvania (152.4 m), our model indicates the points beyond this distance from wellsites with non-zero concentration. These points are specified by a grid size of 25-by-25 meters. Population density at each of the grids is estimated using Census data for the closest census block. The model assumes even distribution of population for each census block, similar to the approach previously considered by Czolowski et al. (2017). Thus, expected number of affected people is equal to the product of population density (persons per area of the block) at each grid and the grid size (625 m2).
For each wellpad, our model uses the “SPUD DATE” as well as the total drilled and fractured length of wells to specify the activities time window. It then uses the hourly wind data measurements at the closest weather monitoring station to that wellsite during the corresponding time window to simulate the PM2.5 concentrations at each hour of operation. The model assumes zero accumulation of emissions from hour to hour, but it accounts for accumulation of emissions dispersed from multiple sources at every hour to simulate PM2.5 concentrations at each grid.
The reference exposure ( in Eq. 1) at any location is assumed to be equal to the level of PM2.5 before the air quality was changed due to emissions from shale gas development. Therefore, the model estimates the increase in the incidence of any health impact from the expected number of incidences corresponding to no shale gas development in those areas. Incidence rate data was not available for all identified health impacts which were identified to be associated to chronic exposure to PM2.5 emissions. Thus, the model estimated the increase in the health impacts with available incidence rate.
For the purpose of this study, we used the unit of person-years in calculating the potential number of affected people by PM2.5 emissions from shale gas development over years. Person-years is an epidemiologic jargon which allows to accommodate for persons coming into and leaving the study when quantifying incidence rate of a health outcome. Therefore, it accounts for the possibility that same person be counted in multiple years.
Our model only accounts for emissions from diesel engines at the wellsite and does not consider emissions from trucks, fugitive dust, mineral dust from proppant handling, or road dusts. We also assumed a flat terrain (no buildings or similar structures) around the wellsite, as most of the shale gas wellsites in Pennsylvania are generally located in rural areas with few or no such a structure within their close vicinity. We assumed constant weather conditions for each hour of operation at each wellsite as well as the same hourly emission rates for drilling and hydraulic fracturing stages of all wells. The model assumes zero precipitation which is a typical assumption in Gaussian plume modeling. This study does not account for the role of delays due to weather, equipment failures, or other factors on the actual duration of development operations at the wellsites, as this information is not publicly available. Such delays would extend the period of emissions being released from the diesel engines in-use and could increase the annual average concentration within the vicinity of shale gas wellsite. This increase implies higher relative risk value (Eq. 1) and therefore, greater increase in the incidence of different health impacts (Eq. 2).