Domain Editor-in-Chief: Detlev Helmig; University of Colorado Boulder, US
Guest Editor: Stefan Schwietzke; University of Colorado Boulder, US

## 1. Introduction

Air quality across the United States, particularly in large cities, has improved substantially in recent times (e.g. Popp 2006; Frost et al., 2006; Russell et al., 2012; Duncan et al., 2013). Federal efforts and increasingly stringent environmental policies played a crucial role in paving the way for technological advancements. Modern vehicles are more fuel-efficient than ever before with improved combustion systems and a growing fleet of electric vehicles. Coal-fired power stations are being displaced with alternative energy sources including renewables and natural gas, that are certainly favourable from an environmental standpoint. But one such innovation—the combination of horizontal drilling with hydraulic fracturing, or fracking, to extract oil and natural gas (O&NG) from previously impenetrable shale formations—has emerged as a point of consternation, that could be offsetting these air quality benefits regionally. While the industry is clearly beneficial for the economy and energy security, it is known to emit greenhouse gases (GHGs) and other air pollutants from operations, which could be highly detrimental to the regional environment, climate, and public health (Howarth et al., 2011).

As a consequence of these concerns, numerous studies have sought to understand the air quality and climate impacts of the oil and natural gas industry in the U.S., the results of which so far have been compelling. Enhancements in ambient levels of methane (CH4), which is a highly potent GHG, have been discovered (e.g. Karion et al., 2013; Miller et al., 2013; Pétron et al., 2014; Caulton et al., 2014; Schneising et al., 2014). Emissions from O&NG activities have been found to have high concentrations of non-methane hydrocarbons (NMHCs) (e.g. Edwards et al., 2013; Oltmans et al., 2013; Thompson et al., 2014; Swarthout et al., 2015; Prenni et al., 2016; Gilman et al., 2013), which are comprised of a host of air pollutants that are regulated in some regions as a result of their damaging affects to the environment and public health. Sulphur dioxide (SO2), which when oxidised can contribute to new particle formation (Sipilä et al., 2010) and precipitation and particle acidity (Hegg and Hobbs, 1981; Weber et al., 2016), is also increasing (e.g. McLinden et al., 2012, 2014, 2016; Prenni et al., 2016). Numerous studies have also shown that regional ozone (O3), a health hazard and criteria pollutant regulated by the U.S. Environmental Protection Agency, is also influenced by O&NG emissions during both summer (e.g. Kemball-Cook et al., 2010; McDuffie et al., 2016; Pacsi et al., 2015; Rodriguez et al., 2009) and winter seasons (e.g. Ahmadov et al., 2015; Carter and Seinfeld, 2012; Edwards et al., 2013, 2014; Field et al., 2014; Helmig et al., 2014; Oltmans et al., 2013; Schnell et al., 2009).

Nitrogen oxides (NOx = NO2 + NO), which is comprised of a family of highly reactive and toxic gases, is another major air pollutant that is increasing in regions of shale O&NG activity (e.g. Duncan et al., 2016; McLinden et al., 2016; Prenni et al., 2016). NOx pollution has both atmospheric and public health implications as these compounds have long been reported to be a major precursor of tropospheric ozone (e.g. Liu et al., 1980; Jacob et al., 2016; Schultz et al., 1999) and cause respiratory problems (e.g. Detels et al., 1991; Kumar et al., 2007). As a result, the emissions of NOx are regulated in many countries including the United States (EPA, 1999). Given that there were an estimated 1.7 million conventional and unconventional O&NG wells in the U.S. in 2016 (The FracTracker Alliance, 2016), improving our understanding of NOx emissions from the industry is imperative in order to understand the associated environmental and health impacts. To date, however, only a few studies have focused on changes in NOx emissions in regions of shale oil and natural gas extraction (e.g. McLinden et al., 2012; Li et al., 2016; Duncan et al., 2016; Prenni et al., 2016; Abeleira & Farmer, 2017). In Canada, a series of studies have detected strong growth in NOx ambient mixing ratios over the oil sands, which are actively producing oil and natural gas by means of unconventional techniques (McLinden et al., 2012, 2014, 2016). In the United States, the Bakken Shale located in North Dakota is the primary study region. For example, Prenni et al. (2016) inferred significant increasing trends in NO2 mixing ratios—found to be consistent with the growth of regional oil and natural gas production—in the Bakken region based on ground based measurements. Similarly, Li et al. (2016) observed increases in the levels of tropospheric NO2 over Bakken using retrievals from the Ozone Monitoring Instrument (OMI) on board the NASA Aura satellite. Duncan et al. (2016) used the same satellite product to detect positive spatial changes (∼30%) in NOx over three different shale plays of the U.S. (Bakken, Permian, and Eagle Ford). Whilst the results from these works are significant, there are no such analyses as of yet that seek to quantify the long-term evolution of NOx over numerous O&NG regions across the U.S. Furthermore, improved attribution techniques to link positive NOx inferences to regional oil and natural gas activities are required.

We present here a quantitative assessment of the evolution of nitrogen dioxide levels over multiple regions of shale oil and natural gas operations in the United States and study the relationship between observations and regional oil and natural gas production. For this purpose, we use satellite-based measurements of tropospheric NO2 from OMI. Using remote sensing data for this study is highly advantageous given that industrial operations tend to be located in remote regions, where ground-level monitoring stations are rare, even in developed nations. Further, the OMI-NO2 dataset has been shown by many studies to serve as an effective proxy for NOx as it correlates well with surface NO2 and it can be used to identify and monitor NOx changes from stationary point sources (e.g. Kim et al., 2006; Lin et al., 2010; Russell et al., 2012; McLinden et al., 2012; Duncan et al., 2013; McLinden et al., 2014; Vinken et al., 2014; Duncan et al., 2016; Li et al., 2016).

In this paper, we firstly describe the OMI-NO2 observations, as well as auxiliary datasets used in our analysis (Section 2). In Section 3, we describe our methods in processing OMI-NO2 and extracting trends through an additive decomposition approach, as well as statistically testing the relationship between observed trends and regional oil or gas production. We also describe our methodology to maximise signals from O&NG activities that incorporates the location of wells, large cities, and power stations. Finally, we present our results that show: (1) strong spatial and temporal increases in tropospheric NO2 over three intense oil-producing shale plays, which correlate with the expansion in regional O&NG activity, and (2) a weakening of decreasing NO2 trends that is likely associated with regional O&NG activities (Section 4).

## 2. Data

### 2.1. OMI-NO2

Our analysis uses measurements of NO2 from OMI on-board the NASA EOS-Aura spacecraft. OMI is a Dutch-Finnish spectrometer that was launched in July 2004, which provides hyperspectral measurements in 740 wavelength bands within an ultraviolet and visible (UV-Vis) spectral range (270–500 nm) (Levelt et al., 2006). The instrument provides near-daily global coverage through 14 orbits per day.

The dataset used for this study is a high-resolution (0.1° latitude × 0.1° longitude) Level-3 OMI-NO2 product containing daily geolocated measurements. The measurements are obtained in the broad visible spectral window between 405–465 nm and processed using the version 2.1 retrieval algorithm. Our data are tropospheric vertical column densities (VCDs), which is the sum of NO2 molecules between the Earth’s surface and the tropopause per unit area (i.e. molec. cm–2). This is calculated by subtracting the stratospheric content from the total NO2 column (Duncan et al., 2016). The uncertainty in the individual OMI-NO2 slant column is ∼0.75 × 1015 molec. cm–2 (Lamsal et al., 2014, references therein). The entire archive of OMI data is publicly available and can be retrieved from the website of NASA Goddard Earth Sciences Data Active Archive Center (http://disc.sci.gsfc.nasa.gov/). OMI began its earth observations in late 2004, and so we select and analyse data between the timeframe of January 2005 to December 2015, which provides a complete archive for 11 years. OMI developed a partial field of view blockage in 2007 that affects a portion of its 60 cross-track scan positions—known as the row anomaly (Li et al., 2016). To mitigate against this, we use data from positions 5–23, which are unaffected by the row anomaly. Coverage from OMI during winter months can be sparse or entirely absent as a result of snow, ice, or polar night (Li et al., 2016). This may be particularly problematic when sampling from high northern latitudes, and we take this into account in our analysis as explained in Section 3. We also acknowledge that there is now a newer version (3.0) of the OMI-NO2 data product. However, works such as those presented in this paper, which study long-term trends, are not affected by the use of the earlier version of the data (Krotkov et al., 2017).

### 2.2. Auxiliary datasets

Additional auxiliary datasets are also applied in this study for the purpose of exploring the relationship between the OMI-NO2 VCD and regional oil or gas production. Data on national production statistics and drilling rig counts between 2005 and 2015 are taken from the U.S. Energy Information Administration’s database (EIA, 2017a). These are filtered and sorted according to the extraction basin. Spatial data on the extents of shale basins are from the EIA (2016a). Locations of both conventional and unconventional oil and natural gas wells across the U.S. are downloaded from the The FracTracker Alliance (2016).

Emissions of nitrogen oxides are well-understood to be closely linked with (a) population density (Schneider et al., 2015, references therein) and (b) coal-fired power plants (de Foy et al., 2015, references therein). Moreover, emissions of NOx from large cities and coal-fired power plants are known to be decreasing across the United States (Duncan et al., 2016, references therein), which could dampen the signal from O&NG activities. Therefore, in order to minimise the NO2 contribution from additional sources such as large urban areas and power plants, we incorporate locations of coal-fired power stations, as well as cities with a population greater than 400,000 inhabitants, into our analysis. Population data are acquired from the U.S. Census Bureau (2016) and locations of coal power stations from EIA (2016b). In addition, prescribed and forest fires can emit significant amount of NO2 over the U.S. (e.g. van der Werf et al., 2010). Therefore, we analyse fire activity data from the National Interagency Fire Center (NIFC, 2016) between 2005 and 2015 to identify major fire events and/or long-term trends of fires that may interfere with our results.

## 3. Methods

The methods applied in this work are based-upon those that have already been applied extensively in other studies (e.g. Lamsal et al., 2015; Schneider et al., 2015). We provide a summary of the analysis techniques here and describe in more detail the data filtering method that we use to maximise the signal from O&NG activities.

### 3.1. Site Selection

According to The FracTracker Alliance (2016), there were ∼1.7 million active oil and natural gas wells across 35 of the 50 U.S. states in 2016. In this study, we focus on those areas where total production of oil or gas was higher than 10 million barrels per day (mb d–1) and 95 million cubic meters per day (mcm d–1), respectively, at the close of 2015. This includes production from both conventional and unconventional methods. Only 8 shale plays fit this criteria: Bakken, Eagle Ford, Permian, Niobrara-Codell, Marcellus, Haynesville, Barnett, and Utica (EIA, 2017a). Figure 1 shows the distribution of O&NG wells across the United States in 2015 and highlights the location of these main 8 O&NG regions, which form the basis of our study. Marcellus and Utica are located within the same proximity over eastern U.S., and so we consider them as one study area (henceforth referred to as Marcellus-Utica), thus making a total of 7 study regions for this work. At the end of 2015, the sum of production from these shale plays accounted for over 95% and 75% of the total oil and natural gas production in the U.S., respectively, making them the most prolific production regions in the country. Therefore, studying these regions provide us a robust representation of the NOx emissions in areas of shale O&NG extraction in the U.S.

Figure 1

Location of wells across the continental United States in 2015. Oil (red) and gas (pink) wells are shown. Seven box regions denoted in are: Bakken [102–104° W, 47–49° N], Barnett [97–99° W, 32–34° N], Eagle Ford [97–100° W, 28–30° N], Haynesville [93–95° W, 31.3–32.8° N], Marcellus-Utica [75–83° W, 37–43° N], Niobrara-Codell [104–105° W, 40–41° N], and Permian [101–104° W, 31–33° N]. DOI: https://doi.org/10.1525/elementa.259.f1

Figure 2 shows a breakdown of oil and natural gas production for December 2015 from each of the 7 study regions. Clearly, there is a primary commodity in each region of either oil or gas with the exception of Niobrara-Codell where production is relatively balanced. The engineering infrastructure required in regions where oil is the primary commodity (e.g. Bakken, Eagle Ford, Permian) is different to that in primarily gas producing regions (e.g. Barnett, Haynesville, Marcellus-Utica). For example, levels of associated gas flaring and venting in the Bakken shale play are known to be large due to a lack of sufficient infrastructure needed to capture, transport and utilise the natural gas (EIA, 2016c). Additional operational differences across each site could include varying methods for flowback fluid storage and removal, as well as the presence of regional processing facilities. Further, oil production sites typically comprise more compressors and heavy-duty vehicles (e.g. trucks to transport oil) (Lattanzio, 2013). Therefore, we anticipate that there may be differences in NOx pollution depending on the main product of the shale play.

Figure 2

Total oil and gas production from each study region at the close of our study period (December 2015). DOI: https://doi.org/10.1525/elementa.259.f2

### 3.2. Data Analysis

We processed the Level-3 files containing raw daily OMI-NO2 data by aggregating the values to monthly averages per grid cell (0.1° × 0.1°; ∼10 km ×10 km). At this stage, we also filtered for outliers by removing those data that fell outside of the 5th and 95th percentiles since we are interested in studying long-term trends, and so this process of filtering makes the results less sensitive to the effects of extreme values, likely associated with some transient synoptic features, local sources or retrieval artifacts (Lamsal et al., 2015).

Next, we extracted data from each of the 7 study regions using the geographic domains noted in Figure 1, and calculated monthly averages across the studied regions. During the data extraction, we applied a data filtering criteria to maximise the signal from O&NG operations and minimise the impact of other local sources. We included data from only those grids that contained at least five wells per pixel. Further, we excluded grids that were within a 10-km perimeter (roughly the width of 1-pixel) from large cities (population of > 400,000 inhabitants) and coal-fired power stations. This removed around 58% of the data on average within the study domains. To include statistically robust monthly averages in the analysis, we only considered those months that contained at least five daily OMI-NO2 VCD observations.

To calculate trends from a seasonal dataset, we deseasonalised the monthly averages time series using the seasonal trend (STL) additive decomposition procedure developed by Cleveland et al. (1990). STL uses a series of applications of the locally weighted regression smoother (LOESS); it is a powerful and versatile filtering procedure that allows for the inference of trends by removing the inherent variability in the data, which arises mainly as a result of seasonal effects. This technique has been widely used to extract the interannual and trend signals from non-stationary and noisy time series in climate sciences (Dufresne et al., 2013), including those products generated from OMI (e.g. Lamsal et al., 2015; Shukla et al., 2017).

The STL method separates the time series of monthly NO2 values into three independent components: a seasonal component S(t), a deseasonalised component D(t) and a residual or noise component R(t), where t denotes time (month). Hence, we can express the time series as:

(1)

The term D(t) quantifies the deseasonalised time series, as well as the low-frequency variations in the time series. S(t) describes the variation in the data at or near the seasonal frequency, which in our case is one cycle per year. Lastly, R(t) is the noise beyond the trend and seasonal components.

As mentioned above, OMI data over high northern latitudes during winter months can be limited as a result of snow and ice. As a result, this led to a few missing months in the generated time series Y(t), particularly that which was generated for the Bakken region, where 5 of the 132 months were missing. The few missing months were filled in with the mean of the entire time series to enable the STL procedure to be performed, which requires a complete dataset. This filling approach, when compared with alternative methods such as taking winter averages or interpolation, did not impact the final result of the time series significantly given that the number of missing months are small relative to the size of the data array.

Finally, we determined the absolute (molec. cm–2 yr–1) and relative (% yr–1) changes from the decomposed time series. Usually, most oil and natural gas operations begin with an exploratory phase, where production levels are low, after which extraction is ramped up if deemed economically and technically feasible. As such, the magnitude of NOx emissions would vary between the two phases. Initially, relatively low levels of NOx would be generated from heavy-duty road vehicles (e.g. transport trucks) and equipment such as off-road diesel engines that power pumps, compressors, and drilling rigs, as well as gas flares. In the next phase, more of these equipment are installed and are operated at higher intensities, which would in turn emits more NOx (Lattanzio, 2013). These operational phases can also be witnessed in the modern oil and natural gas industry of the U.S. Though the industry has existed for over a decade, it has been the period after 2010 that can be considered to be most significant across almost all shale plays, where an exponential increase in unconventional O&NG production can be witnessed. Therefore, we focused on two distinct periods for computing the rate of change in NOx, which are: (1) before (2005–2009) and (2) after (2010–2015) the exponential increases. Absolute changes were computed by applying two least-squares regression models onto the deseasonalised time series D(t), one for each period before and after the exponential increases, and subsequently taking the derivatives of the model function. The corresponding uncertainty was taken as the standard error in the gradient estimator of the regression model. Relative changes were computed following Schneider et al. (2015) as:

$\Delta =\frac{D\left(t\right)}{{\overline{D}}_{t}}×100$
(2)

where ${\overline{D}}_{t}$ is a long-term mean of Dt. By using a long-term mean as opposed to a single reference year at the start of the time series, the results are less sensitive to potential outliers (Schneider et al., 2015). The derived changes were deemed to be statistically significant at the 95% confidence level if the rates of change were larger than twice the corresponding standard error.

## 4. Results and Discussion

### 4.1. Spatial analysis of tropospheric NO2 observations

Figure 3 shows absolute changes (×1015 molec. cm–2) in the OMI-NO2 VCD between 2005 and 2015, which reveals the major spatial changes across the United States. It should be noted that the spatial analysis incorporates all grid cells (i.e. grid cells co-located with urban areas and power plants have not been removed at this stage). Anthropogenic emissions of NOx have been known to be decreasing across the U.S. for the past decade (Lamsal et al., 2015). Results derived from our processing of NO2 measurements from OMI are consistent in that they also show large-scale decadal decreases of NOx over the U.S. of similar magnitude. In almost all cases, wherever there is a large city and/or power station, there is a decrease in the OMI-NO2 VCD of typically greater than 2 × 1015 molec. cm–2. This is far more obvious over eastern U.S. (70–90° W), where there is a higher population and number of power stations density per unit area.

Figure 3

Absolute changes (×1015 molec cm–2) in the OMI-NO2 VCD between 2005 and 2015. Large cities (red triangles) and power stations (grey circles) are overlaid. DOI: https://doi.org/10.1525/elementa.259.f3

Figure 4 shows relative changes (%) in the OMI-NO2 VCD between 2005 and 2015 for our study regions: Bakken, Permian, Eagle Ford, Barnett, and Haynesville, Marcellus-Utica, and Niobrara-Codell. The domains used for the analysis in this study are marked by grey boxes. From our spatial analysis, we observe widespread increases in the OMI-NO2 VCD in three particular regions: Bakken, Eagle Ford, and Permian. These regions are the most prolific zones of oil production in the U.S. (Figure 2) and as such have the highest well density per unit area. Furthermore, the three regions are relatively more remote and distant from urbanised regions, in which O&NG pollution may be the dominant source. Across these shale plays, some areas exhibit increases of more than 1 × 1015 molec. cm–2 (up to 30%) in 2015 as compared with background levels from 2005. These results are consistent with the observations of Duncan et al. (2016) and Li et al. (2016). All other study regions (Barnett, Haynesville, Marcellus-Utica, and Niobrara-Codell) do not display such stark increases. Rather, we infer decreases across these four shale plays from our spatial analysis, in accordance with national trends in NOx emissions. These regions are closely located with urban environments. For example, across the Marcellus-Utica shale—primarily located in the state of Pennsylvania (PA)—there are a large number of sparsely distributed coal-fired power plants, and the state as a whole has a relatively large population density. This explains the evidently dramatic declines of more than 30% over some regions. Despite our efforts to maximise the signal from O&NG activities, then, it is likely that emissions from surrounding urbanised regions are dominating signals.

Figure 4

Relative changes (%) in OMI-NO2 VCD between 2005 and 2015 across all of our study regions. Wells (black circles), urban boundaries from large cities (dotted lines), and power stations (pink circles) are overlaid. DOI: https://doi.org/10.1525/elementa.259.f4

### 4.2. Relationship between OMI-NO2 and shale activity

To study the relationship between OMI-NO2 observations and O&NG activities, we correlate the annual deseasonalised NO2 values against regional oil or gas production. In this analysis and onwards, we remove the effect of urban areas and power plants from the analysis (i.e. we present filtered data). Figure 5 shows the results of the correlation, in which statistically significant correlations can be observed across all regions. It is important to note that the OMI-NO2 VCD in highly populated regions (e.g. Marcellus-Utica, Niobrara-Codell, and Barnett) are substantially larger than in more remote areas (e.g. Bakken). In other words, relative background levels are much higher. To further determine the influence of drilling operations, we also correlate annual NO2 against number of drilling rigs (Figure S1), and obtain similar relationships across all regions.

Figure 5

Relationship between deseasonalised annual OMI NO2 VCD and regional annual oil or gas production. Deseasonalised annual OMI NO2 VCD (solid lines), with 95th percentile confidence bounds (dashed lines), and regional annual oil or gas production (bars) are shown. Note that annual OMI-NO2 is plotted with different axis ranges to make the relationship more apparent. DOI: https://doi.org/10.1525/elementa.259.f5

For the primarily oil producing zones of Bakken, Eagle Ford, and Permian, we find a strong correlation between OMI-NO2 VCD and annual oil production, with a discernible increase in OMI-NO2 after 2010 that coincides with the timeframe in which regional production surged. These results are statistically significant (r2 = 0.6–0.9; p ≪ 0.01). A similar relationship, but slightly weaker correlation, is found between OMI-NO2 VCD and number of drilling rigs (r2 = 0.4–0.6), which also includes the oil producing basin of Niobrara-Codell. Across these three regions, gas flaring—the process by which natural gas associated with oil reserves is combusted into the atmosphere—was observed by NASA Suomi NPP VIIRS in 2010 (Hillger et al., 2013). Furthermore, Duncan et al. (2016) noted elevated OMI-NO2 levels in these three regions, which spatially matched the VIIRS images of gas flaring activities. There is a strong indication that O&NG activities (e.g. gas flaring, compressors and other drilling equipment, and heavy duty vehicles) are giving rise to high, regional concentrations of nitrogen oxides. Though we cannot conclusively discount other contributing factors such as urban changes and power plant activities.

Over the Barnett, Haynesville, Niobrara-Codell, and Marcellus-Utica shale plays, we do not find such a strong relationship and our correlation tests indicate that there is no statistically significant relationship between NO2 changes and regional gas production (p > 0.01). However, we find a clear change in the decreasing tendency of OMI-NO2 that coincides with the increase in gas production for Barnett, Haynesville and Marcellus-Utica, and in oil production for Niobrara-Codell. One such possible explanation for the lack of a relationship with production rates may be the relatively low-levels of production in these regions (e.g. Niobrara-Codell, Barnett, and Haynesville). Any changes in local ambient NOx concentrations might be too small to be detected from OMI. In the Marcellus-Utica region, where production levels are high, it may be that increasing NO2 signals from oil and natural gas have been dampened by strong regional NOx decreases as a result of the more efficient power stations and high population density regions. Furthermore, natural gas is a commercial product across all four of these basins (Figure 2), and so it is logical to predict that gas flaring levels in these regions would be relatively lower since it would not be economically favourable. We also note that the regulations on gas flaring in North Dakota and Texas have been relatively lax over the past decade or so, whereas flaring is tightly regulated in Colorado. This may also be another reason for the lack of a statistically significant relationship over the Niobrara-Codell basin.

### 4.3. Influence of O&NG operations on decadal trends

Figure 6 shows the OMI-NO2 time series Y(t) along with the derived trend component D(t) from the STL procedure. Computed trends from D(t) before and after the intensification of oil and natural gas production are also shown in relative (% yr–1) units. We consider the relatively low oil and natural gas production period as 2005–2009 and high production period as 2010–2015 when oil and natural gas production across the nation soared exponentially or had peaked, as shown in Figure 5. In addition to all of our study regions, we also present the national trend including all basins for the lower 48 states inferred from OMI-NO2 for a comparison. All trends presented are statistically significant at the 95% confidence level (Section 3.2), with the exception of Bakken during the low O&NG production phase.

Figure 6

Decadal OMI NO2 observations over the study shale regions. Monthly OMI NO2 (grey circles), deseasonalised monthly OMI-NO2 (D(t)) (black dashed line), and trends for relatively low- (green line) and high- (red line) production periods are shown. Relative changes for the periods are shown in the inset. DOI: https://doi.org/10.1525/elementa.259.f6

Our analysis shows declining trends (between –8 and –2% yr–1) across all regions in the low production phase in line with the national trend (–4.6 ± 0.1% yr–1), with the exception of Bakken, where a small positive change of 0.3 ± 0.2% yr–1 is inferred. Contrariwise, in the high production period, we observe distinctive changes. In the Bakken, Eagle Ford, Permian, and Niobrara-Codell shale plays, the deseasonalised OMI-NO2 time series begins to rise dramatically relative to the low-production phase at rates of 4.5 ± 0.3, 3.0 ± 0.3, 3.1 ± 0.2 and 1.3 ± 0.2% yr–1, respectively. Li et al. (2016) studied long-term OMI-NO2 changes over the Bakken region and derived an absolute magnitude of 0.017 ± 0.005 molec. cm–2 yr–1 between 2005 and 2015. In comparison, our analysis shows a slightly larger increase in the same period (0.030 ± 0.002 molec. cm–2 yr–1). This discrepancy is likely associated with the difference in methodologies applied in the two studies. For example, Li et al. (2016) focus on summertime trends only and so do not incorporate winter months into their analyses. Also, Li et al. (2016) maximise signals attributable to O&NG activities by setting a brightness threshold to their study area using night time images from satellites. Further, Li et al. (2016) adopted a different technique to deseasonalise the OMI-NO2 time series.

The trends in the high-production phase observed for other three regions—Marcellus-Utica, Barnett, and Haynesville—do not exhibit a similarly dramatic growth. Instead, the trend begins to flatten or resumes declining but at a shallower gradient (–2 to 0.4% yr–1), which follows our observations of the national trend (–0.3 ± 0.1% yr–1). Results for the Marcellus-Utica region are noteworthy given that a slope of magnitude 8 ± 0.3% yr–1 from the low-production phase flattens out to 0.4 ± 0.2% yr–1 as production is ramped up between 2010 and 2015. Our analysis suggests that declining trends in tropospheric NO2 VCD have been halted over this basin since the instigation of regional oil and gas expansion (Figure 5). Therefore, there may be an implication here that O&NG activities have caused a slow down in the observed declining trend of surface NOx levels.

Our analysis has showed that the high-production phase trends in regions where oil is the dominant product are different than those in primarily gas producing regions, particularly when compared to the national trends. Further, there is a notable change in the slope between the low- and high-production phase of the national trend, too. The disparity in the observed trends could be due to the inherent differences in operational equipment (e.g. gas flares, drilling compressors etc.) and production site characteristics (e.g. more large vehicles on oil-production sites for transportation of oil). Other variables contributing to the differences in the observed results could include regional changes in population and the closing of local power plants. The flattening of NO2 decreases after 2010 across the U.S. is consistent to previous observations (e.g. Russell et al., 2012; Tong et al., 2015), and is mainly attributed to a recovery after the 2008 Global Economic Recession.

Regional fires, population changes, and power plant pollution can also affect tropospheric NO2 trends as discussed previously. In our analysis, the influence of individual fires is minimised since we consider only OMI-NO2 observations within the 95th percentiles, which removes the influence of extreme values. To further determine whether our trends are driven by any interannual variability in fire NOx emissions, we analyse long-term trends (2005–15) of fire activity in all U.S. states corresponding to our geographical domains (Figure S2). This analysis does not show any consistent relationship between our derived NO2 trends and the decadal changes in fire activity. Populations across each state studied in this work are reported to have been increasing steadily within the time frame of this work (U.S. Census Bureau, 2016). At the same time, production from coal-fired power stations has declined sharply in the U.S. since 2008, due to competition with natural gas production and insufficient growth in electricity demand (EIA, 2017b). In this work, we applied a data filtering methodology to minimise the influence of population changes and power plant pollution on the decadal trends; in spite of which, however, some residual influence may have remained as mentioned above.

Figure 7 shows the changes in the deseasonalised OMI-NO2 for all 7 study regions in the low and high oil and natural gas production phases, and summarises results shown in Figure 6. We calculate the relative (%) changes by computing the difference between 2005–2009 and 2010–2015, i.e. by subtracting the annual OMI NO2 VCD in 2005 from that of 2009 and in 2010 from 2015, respectively. We observe that the percentage difference in the OMI-NO2 VCD increases significantly (6–18%) across all four oil producing shale plays—Bakken, Eagle Ford, Permian, and Niobrara-Codell. Changes in the Niobrara-Codell shale play are relatively lower in comparison to the other three oil producer regions, which is likely due to the combination of lower levels of oil production and being co-located to an urbanised region, which may weaken the NO2 signal from O&NG activities. Additional factors that can influence the results across these regions include: rates of population change, differences in the age and type of equipment in each basin, as well as the O&NG regulations in each state. It is interesting to note the dramatic NO2 change between the low (–32%) and high (–2%) production phases over the Marcellus-Utica region. Over this region, the percentage NO2 change in the high production phase is 16 times smaller than in the low production phase. This change is significantly larger than that observed over the lower 48 states (–14% and –3% in the low and high production phases, respectively; not shown). This result underlines the potentially negative effect that emissions of NOx can have in offsetting the benefits of anthropogenic NOx reductions over these regions.

Figure 7

Relative changes (%) in OMI NO2 VCD between the relatively low-production (2005–2009) and high-production (2010–2015) periods. DOI: https://doi.org/10.1525/elementa.259.f7

### 4.4. Quantification of tropospheric NO2 changes associated with shale oil operations

To further quantify the influence of oil extraction on tropospheric NO2 and study further implications, we compare the enhancement of the deseasonalised OMI-NO2 VCD to the regional annual oil production during the high-production phase across the major oil production regions: Bakken, Eagle Ford, and Permian. We focus our analysis on these three regions since our previous results show a significant correlation between OMI-NO2 VCD and annual oil production (r2 = 0.6–0.9; Figure 5).

Figure 8 presents an initial relationship between annual oil production rates and NO2 enhancements (ΔNO2). We consider the enhancement as the change in the deseasonalised OMI-NO2 VCD relative to background levels (i.e. the first three years of the low-production phase: 2005–2007). We find a significant linear relationship between regional oil operations and tropospheric NO2 enhancements (r2 = 0.50; p ≪ 0.01). Our analysis indicates that there is an increase of approximately 0.15 × 1015 molec. cm–2 of NO2 for an increase in the daily oil production of a million barrels.

Figure 8

Relationship between annual oil production rates and OMI-NO2 enhancements (ΔNO2) during the high oil production phase (2010–2015). DOI: https://doi.org/10.1525/elementa.259.f8

Emissions of NOx, and associated NMHCs, can drive elevated surface ozone, in particular downwind from the sources (Rodriguez et al., 2009). Helmig et al. (2016) estimate an increase of 2.5 ppbv in summertime surface ozone over the northern hemisphere, as a result of increases in NMHCs from O&NG operations from 2004–2014. In their modelling analysis, emissions of nitrogen oxides were assumed to be constant. However, our analysis suggests that NOx have been increasing over the past decade as a result of O&NG operations in some regions of the U.S. To understand if NOx emitted from oil and natural gas operations is a concern for attaining ozone air quality standards, it is important to accurately characterise NOx emissions from these operations, and consider them in future modelling analyses. It should be noted, however, that although we observe levels of NO2 to be increasing across some of our study regions, the actual emissions of NOx from O&NG operations are still relatively smaller than those from urban sources and the ambient NOx mixing ratios are still relatively low compared to cities.

## 5. Conclusions

In summary, we analysed the changes in tropospheric NO2 over regions of shale oil and natural gas activity in the U.S. in the period between 2005 and 2015 using measurements of tropospheric NO2 from the Ozone Monitoring Instrument aboard the EOS-Aura aircraft. A total of 7 shale plays were studied in detailed in this work, and the national trend averaged over the lower 48 states was inferred for a comparison. Our results show that tropospheric NO2 pollution has risen dramatically since 2010 in the Bakken, Eagle Ford, Niobrara-Codell, and Permian shale plays, which correlates well with regional oil production and/or number of drilling rigs. Over Marcellus-Utica, we find that declining trends in tropospheric NO2 VCD have been halted. It is likely that the practice of gas flaring and oil operations in these regions is increasing ambient NOx, although urban populations and distribution, as well as the location and fuel types of operating power plants may be contributing to this trend. No such increase was found over the other two study regions—Barnett, and Haynesville. Hydraulic fracturing operations in the United States are a fairly recent venture. In most cases, the industry in each state has seen a rapid increase in only the last five years or so. Therefore, NOx patterns in proximate areas to oil and natural gas activity must be monitored over the next few years.

This paper has demonstrated yet another example of the advantages of using satellite observations to study NOx pollution from point sources located in remote regions. Data acquisition, global coverage, and high spatial resolution are all attractive features that make this sort of work feasible, however it is important to validate these studies using bottom-up methodologies. Other local and regional sources of NOx can affect tropospheric NO2. In this work, we have developed a method to maximise the signal from O&NG activities by incorporating locations of wells, power stations, and large cities, and an analysis of fire activity. Yet despite our efforts, drawing out signals from some oil and natural gas fields can be incredibly challenging, particularly if regional production levels are modest, which can dampen the signal to the point they are below the instrumental detection limits. Further complications can arise if the O&NG operations are closely located to other sources, where the signal can also be swamped by regional trends.

Despite the recent volatility in oil and natural gas markets, shale operations in the U.S., and indeed other regions of the world, are only expected to grow given that operational efficiencies continue to improve, which are slashing production costs. As shale oil and natural gas production continues to expand, our results show that tropospheric NO2 levels may continue to rise, unless more efficient control technologies are applied to these operations. Monitoring NOx emissions from shale oil and natural gas operations must continue to further evaluate the influence of such as activities on future regional air quality and climate.

## Data Accessibility Statement

OMI-NO2 data used in this study is publicly available and can be retrieved from the website of NASA Goddard Earth Sciences Data Active Archive Center (http://disc.sci.gsfc.nasa.gov/).