Spatio-temporal patterns in water surface temperature from Landsat time series data in the Chesapeake Bay, U.S.A.

Remote Sensing of Environment 168 (2015) 335–348

Contents lists available at ScienceDirect

Remote Sensing of Environment journal homepage: www.elsevier.com/locate/rse

Spatio-temporal patterns in water surface temperature from Landsat time series data in the Chesapeake Bay, U.S.A. Haiyong Ding a,b,⁎, Andrew J. Elmore b a b

School of Geography and Remote Sensing, Nanjing University of Information Science and Technology, Nanjing 210044, PR China Appalachian Laboratory, University of Maryland Center for Environmental Science, 301 Braddock Rd, Frostburg, MD, USA

a r t i c l e

i n f o

Article history: Received 25 September 2014 Received in revised form 25 June 2015 Accepted 9 July 2015 Available online xxxx Keywords: Water surface temperature Thermal remote sensing Landsat Chesapeake Bay Global change

a b s t r a c t Water temperature is a key factor used to assess biogeochemical cycles and aquatic habitat quality, and is typically monitored using in situ sensors deployed on distributed piers, buoys or mobile platforms. Due to multiple, spatially and temporally varying influences on surface water temperature, (e.g., solar radiation, water depth, tidal stage, turbidity, and industrial activities) multi-temporal remote sensing observations might effectively be used to link these drivers with observations. However, while low- to moderate-resolution sensors provide temporally continuous observations, large pixel sizes have proven problematic in coastal regions where shoreline influences on pixel radiance contaminate the radiation signal leaving the water surface. To alleviate this problem, we used a dense stack of Landsat TM/ETM+ thermal imagery organized by day of year of acquisition, thus producing a climatology of water temperature. We use these data to analyze the spatial patterns of water surface temperature climatology (e.g., average maximum and minimum temperatures) for the past thirty years. We also explore the impact of power plant thermal effluent on water surface temperature climatology of the Chesapeake Bay tributaries. Finally, we divide the Landsat record into 5-year intervals, and calculate the water surface temperature climatology for each period. Trends in water surface temperature over the Landsat record were then compared against air temperature records available from coastal NOAA weather stations. The resulting data exhibit broad scale patterns, such as water surface climatology differences between the main stem of the Bay and its tributaries. The results also include the influence of urbanization and industrialization such as increases in impervious surface area and thermal effluent from power plants. Trends of increasing water surface temperature and air temperature were found for more than 92% of the Bay. While water temperature was always ~2–3° cooler than the air temperature, water temperature has been increasing more rapidly than air temperature in some areas, particularly in the main stem of the Bay and in the Potomac estuary. Therefore, there is a detectable impact of global change on the Chesapeake Bay in the form of an increase in water temperature, which can only partially be explained by increasing air temperatures. © 2015 Elsevier Inc. All rights reserved.

1. Introduction In estuarine systems, water temperature has a strong control on biogeochemical and ecological processes, ultimately regulating the metabolism of these highly productive waters and determining habitat area for aquatic life (Demars & Manson, 2013; Gergory, Beesley, & Kirk, 2000; Handcock et al., 2006; Hedger, Malthus, Folkard, & Atkinson, 2007; Jacobs, Heusinkveld, Kraai, & Paaijmans, 2008; Kaushal et al., 2010; Lisi, Schindler, Bentley, & Pess, 2013; Paaijmans, Takken, Githeko, & Jacobs, 2008; Scheitlin, 2013). For example, through water temperature the timing of key lifecycle events (e.g., spawning) can be altered, emphasizing the importance of water temperature data for understanding fisheries (Lisi et al., 2013). Biogeochemical processes ⁎ Corresponding author at: School of Geography and Remote Sensing, Nanjing University of Information Science and Technology, Nanjing 210044, PR China. E-mail address: [email protected] (H. Ding).

http://dx.doi.org/10.1016/j.rse.2015.07.009 0034-4257/© 2015 Elsevier Inc. All rights reserved.

are influenced by water temperature. For example, eutrophication of coastal waters is an increasing problem in many estuaries of the world and is often viewed as being driven solely by increases in nutrient loading, but in fact spatial patterns in eutrophic conditions can be influenced by water temperature because the rate of chemical processes that lead to entropic conditions scale with water temperature (Lomas & Glibert, 1999; Miller & Harding, 2007). Primary production of organic matter and respiration rates increase with higher temperatures, which leads to higher oxygen demand in the water column under sufficient light and nutrient conditions (Alcântara et al., 2010; Demars & Manson, 2013; Hette-Tronquart, Roussel, Dumont, Archaimbault, & Pont, 2013; Lomas & Glibert, 1999). In contrast, the solubility of oxygen in water decreases as temperature increases. Increasing water surface temperatures can be driven by climate change (Preston, 2004), runoff from impervious surfaces (Tan & Cherkauer, 2013), and thermal effluent from industrial processes (Fisher & Mustard, 2004). Particularly under conditions of climate

336

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

change and coastal urbanization, the importance of monitoring Water Surface Temperature (WST) cannot be understated, as many aspects of estuarine management and restoration are dependent on good, spatially referenced, WST data (Preston, 2004). The temperature of shallow estuaries fluctuate widely between summer and winter, but the spatial pattern of this oscillation is complex (Cowan & Boynton, 1996; Kaushal et al., 2010; Paaijmans et al., 2008) due to tidal cycles, freshwater discharge patterns, and industrial use of water during power generation (Demars & Manson, 2013; Takeshige, Takahashi, Nakata, & Kimura, 2013; Tan & Cherkauer, 2013). These multiple influences on water temperature operate across a range of temporal and spatial scales, suggesting that an ideal monitoring network must be both spatially dense and extensive, and operate continuously over time. Clearly these are characteristics that are difficult to accomplish. In situ water temperature data taken from piers, buoys, and mobile platforms have been used extensively to understand variation in water temperature over a range of temporal scales (Demars & Manson, 2013; Kaushal et al., 2010). However, in situ measurement is expensive and time consuming to deploy over large water bodies. Thermal infrared water leaving radiance acquired from remote sensing sensors has been successfully employed to retrieve WST (Alcântara et al., 2010; Handcock et al., 2006; Hedger et al., 2007; Politi, Cutler, & Rowan, 2012; Takeshige et al., 2013; Tan & Cherkauer, 2013), however, in coastal systems the large pixel sizes of most sensors with appropriate temporal resolution is too large to effectively study near-shore spatial and temporal patterns. In recent years, there have been several applications of Landsat-class data to coastal WST, which significantly reduces the impact of shoreline contamination on the analysis of WST (Fisher & Mustard, 2004; Hook et al., 2004; Lamaro, Mariñelarena, Torrusio, & Sala, 2013; Sobrino, Jimenez-Munoz, & Paolini, 2004; Wloczyk, Richter, Borg, & Neubert, 2006). The objectives of this paper were to (1) describe the annual climatology of WST in the Chesapeake Bay and provide appropriate validation and uncertainty analysis, and (2) describe spatio-temporal trends in WST (1984–2010), including how these trends might relate to urbanization, thermal effluent from power generation facilities, and trends in air temperature due to global climate change. 2. Research sites and data 2.1. Research sites The Chesapeake Bay (Fig. 1) is one of the largest and most productive estuaries in the world (Najjar et al., 2010), supporting an important fishery and opportunities for recreation. As in all estuaries, Chesapeake Bay water is a mixture of fresh river water and coastal oceanic water. This combination of water sources leads to a multitude of habitats and is a key factor determining the high biodiversity and productivity of the Bay (Cowan & Boynton, 1996). The Bay is productive habitat for many aquatic organisms such as blue crab, oyster, and commercially and recreationally important fishes (Cowan & Boynton, 1996; Miller & Harding, 2007; Najjar et al., 2010; Preston, 2004; Scheitlin, 2013). Since 1975, urban land cover, characterized by high impervious surface area, has increased by more than 100% (from ~5% to N10%) in portions of the coastal plain adjacent to the Chesapeake Bay, which has influenced physical, chemical, and biological aspects of water quality (Elmore & Guinn, 2010). 2.2. Data 2.2.1. Satellite data Landsat Thematic Mapper (TM) and Enhanced Thematic Mapper Plus (ETM+) scenes (Row/Path: 15/33) were acquired from the United States Geological Survey (USGS) that covered the northern two-thirds of the Bay and all of the Bay in Maryland, USA. A total of 196 images available from the Landsat 4, 5, and 7 archives were used

in this research (Fig. 2). The Landsat ETM+ sensor has a low gain and high gain band for the thermal channel (10.4 to 12.5 μm) (Lamaro et al., 2013). For this work focused on water temperature the high gain band was used. All images were provided by the Landsat Ecosystem Disturbance Adaptive Processing System (LEDAPS) (Masek et al., 2006; Masek et al., 2008), in which all thermal data were calibrated to top of atmosphere brightness temperature. Thermal bands acquired using TM sensor and ETM+ sensors have 120-m and 60-m spatial resolution, respectively, however, as part of the LEDAPS processing all thermal data were re-sampled to a 30-m grid. As part of LEDAPS, a cloud quality band was generated, which was used here to omit the cloud-covered and cloud shadow pixels. LEDAPS also produces a land-water band that was used to mask out the land area. Atmospheric correction was made using the single-channel algorithm (SC) (Jimenez-Munoz et al., 2009). Land surface temperature Ts was retrieved using Ts ¼ γ


 1 ðψ1 Lsen þ ψ2 Þ þ ψ3 þ δ 3

ð1Þ

where Lsen is at-sensor registered radiance, and ε is the surface emissivity; γ and δ are the parameters dependent on the Planck's function (Jimenez-Munoz et al., 2009). SC algorithm retrieves water surface temperature using the historical thermal-band data and fitting the atmospheric functions which were assumed to be dependent only on the atmospheric water content, and the parameters, ψ1, can be approximated using the a polynomial of the water vapor content according to the following fitting equation 2

3 2 c11 ψ1 4 ψ2 5 ¼ 4 c21 ψ3 c31

c12 c12 c32

32 2 3 c13 w c13 54 w 5 c33 1

ð2Þ

where w is the water vapor content with unit g·cm−2, and coefficients cij are estimated by simulation (Jimenez-Munoz et al., 2009). Atmospheric water vapor content was acquired from the University of Wyoming online database of atmospheric sounding data (Table 1). The SC algorithm (Jimenez-Munoz et al., 2009) has the potential to lead to inaccuracies in water surface temperature at very lower values compared with accuracy at more moderate values. In this study, there were 17 scenes with water content lower 0.5 (more than 82% of these scenes were acquired in winter). 2.2.2. In situ water temperature and weather data To validate Landsat-based WST observations we compiled a concurrently-collected suite of in situ water temperature observations from the Patuxent river near Solomans, MD (data collection methods are described in (Kaushal et al., 2010)). Water temperature observations were made daily at the end of a 200-m long pier at 0.5 m depth. The available in situ data was used in a statistical comparison with concurrently-collected Landsat observations (linear regression) to constrain our understanding of Landsat measurement uncertainty and to document the improvement in uncertainty achieved through the atmospheric correction procedures employed (Fig. 3). In subsequent spatio-temporal analyses, we assume that errors are normally distributed through space and time and therefore do not impact spatial or temporal patterns when averaged over large areas or time. We also sought to compare water temperature with air temperature, both directly against each other and individually against time, to determine if trends in water temperature corroborated with trends in air temperature. Weather observation data reporting annual mean air temperature were obtained from several national weather service stations surrounding the Bay. Data from eight weather stations were used: Annapolis (MD), Baltimore (MD), Cambridge (MD), Chestertown (MD), Fredericksburg (VA), Millington (MD), Princess (MD), and Royal Oak (MD).

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Fig. 1. The study area and points of interest.

Fig. 2. Distribution of study images across year and time of year.

337

338

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Table 1 Atmospheric water content for example image dates (http://wather.uwyo.edu/upperair/ sounding.html, Station: 72402 WAL Wallops Island Observations). Landsat scene data

Atmospheric water content (mm)

Atmospheric water content (g/cm2)

Aug. 27, 1984 May 24, 1990 Jul. 9, 1995 Mar. 8, 2000 Sept. 22, 2005

31.79 17.17 22.87 21.04 19.85

3.179 1.717 2.287 2.104 1.985

3. Data analysis 3.1. Curve fitting Seasonal variation in earth surface temperature exhibits a sinusoid curve over the course of a single year with minimum temperatures generally occurring shortly after January 1 and maximum temperatures occurring in mid to late July. When any reasonably large number of Landsat acquisitions are organized by day of year (DOY), this seasonal cycle is readily apparent in the resulting WST time series (Fig. 4). Therefore, an annual temperature cycle model (ATC) (Weng & Fu, 2014) was applied to fit the WST data, at any day of year, d, using the following fitting model: WST ðdÞ ¼ a  sinð2dπ=365 þ bÞ þ c

ð3Þ

where a, b, c are parameters to be estimated (Table 2). Parameter a is the amplitude of the sine function, which represents the degree of fluctuation from the mean temperature. Parameter c is the mean temperature in the data record being fit, while parameter b represents the phase of the sine function. Two approaches for organizing the Landsat input data were employed. In the first, all of the 196 Landsat images available were used in the curve-fitting model to arrive at an overall average WST climatology for the entire 28-year Landsat record (Fig. 5). In the second approach, the fitting model was applied to sequential five-year subsets of the 28-year record so that trends in WST could be evaluated. The sequential five-year subsets began with the period 1984–1988 and advanced in one-year increments (1985–1989, 1986–1990, etc.) so that a timeseries of WST climatology could be derived. From this time series, the annual mean WST (corresponding to each individual 5-year subset) was extracted to compare with the time series of annual mean air temperature recorded by the nearest weather stations. The number of Landsat scenes in each 5-year window varied through time and ranged from 27 to 66.

Fig. 4. Example annual temperature cycle (ATC) model fitting of the temperature data.

In both types of analysis, pixels covered by cloud or in cloud shadow were masked and removed from the analysis. Thus, for each grid cell a different number of pixels was used to fit Eq. (3). The advantage of using these long time series data was realized in that a sufficient number of clear-sky pixels were found for each grid cell despite cloud contamination in specific Landsat scenes. Although the cloud mask successfully removed most undesirable pixels from the dataset, there were still occasional pixels with very high or very low WST, likely due to water disturbances or view obstructions such as boats or floating debris. To eliminate the influence of such WST observations on the fit of the model, we used a weighted curve fitting method involving two steps. First, the ATC model was applied to get the fitting parameters and the residuals. Residuals were calculated, representing the distance between the measured WST and the fitted curve, and their inverse was used as a first estimate of the uncertainty for each WST measurement. Then, in the second step, we used a weighted ATC model where the weight was given as pi ¼

di N X

ð4Þ dj

j

where di was the inverse of the absolute value of the i-th residual given as:

di ¼

1 : jresiduali j

Fig. 3. Plot of remotely-estimated water surface temperature and brightness temperature versus in situ measured water temperature.

ð5Þ

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348 Table 2 Interpretation of the model parameters. Parameter Representation a b c

339

Table 3 Statistics of the fitting parameters. Range

Units

Amplitude of the sine function, representing the (−20, 20) °C maximum deviation from the mean temperature Phase of the sine function, representing the DOY N0 Radian of min and maximum water temperatures Mean value of the sine function, representing (0, 30) °C the average water temperature

3.2. Estimates of model parameter uncertainty The Monte Carlo bootstrapping method was used to estimate the fitting parameters uncertainties for a randomly selected location in the middle of the Bay. Using standard bootstrapping methods, the ATC was first fit to the entire WST timeseries for the test location. Secondly, random samples were identified and removed (with replacement by a different random sample) to construct a new time series. Following this process we generated 1000 bootstrapped time series for the test location and fit the WST model to each, generating a list of 1000 fitting parameters. We used the bootstrapping results to calculate uncertainty estimates for each model parameter and also to evaluate colinearity between model parameters that might impact the independence of model results.

Parameters

a (°C)

b (Radian)

c (°C)

Minimum value Maximum value Mean value Standard deviation

10.63 18.96 12.48 0.5178

7.21 7.7 7.38 0.0749

12.23 23.18 13.93 0.5719

Due to strong auto-correlation in WST across the Bay we did not feel it was necessary or appropriate to calculate trends in WST for every grid cell. Therefore, for a standard, uniform grid with cells spaced every 40 × 40 pixels covering the entire Bay (i.e., each grid cell contained 1600 Landsat pixels), we performed a linear regression between the mean annual WST derived from the ATC model (calculated via successive 5-year intervals) and year. The slope of each of these relationships was mapped to reveal spatial patterns of the change trend in WST. At select locations, we identified areas (labeled ‘water samples’ in Fig. 1) that had corresponding air temperature records in coastal locations suitable for comparison with WST trends. At these locations we plotted both air and WST against year, and air temperature against WST. In each case, regression statistics were generated and reported to gain insight into the coupling of air and water temperature increases over the past 28 years.

4. Results 3.3. Analysis of resulting WST geo-spatial data 4.1. Water surface temperature comparison with in situ water temperature Following insights gained from the literature (Fisher & Mustard, 2004) we derived four key climatological factors from the ATC (Eq. (3)) at each grid cell:(1) maximum annual temperature, (2) minimum annual temperature, (3) the DOY of the maximum annual temperature and (4) the DOY of the minimum annual temperature. Each of these values were determined from Eq. (3) (with fitted parameters) and represent a mean value over the period spanned by data used to fit the parameters. In other words, when all Landsat data were used to fit the parameters, the average maximum annual temperature is for the period 1984–2011. When data from each of the 5-year windows was used, the mean annual temperature is for that corresponding period.

Water temperature was measured in situ at one location for each remote sensing image. A direct comparison of in situ and remote observations yielded a correlation of 0.9574 (R2 = 0.8867). We recognize that water temperature at 0.5 m depth is not always comparable to surface water temperatures and that the difference is dependent on the sea state (Fisher & Mustard, 2004). Under windy conditions the top 1 m of water can be considered well mixed and in this case WST and in situ water temperatures are similar. However, under calm conditions the surface temperature can be elevated above in situ water temperatures at 0.5 m. Previous research has applied an empirically-derived offset to each remote sensing image to account for this effect (Fisher &

Fig. 5. The fitting parameters of ATC model for WST.

340

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Table 4 Comparison of parameters and their 95% lower upper bounds estimated from fitting model and bootstrapping method. Parameter

Fitting model estimate

95% lower and upper bound from the fitting model

Mean bootstrapping estimate

Bootstrap 95% lower and upper bound

a (°C) b (Radian) c (°C)

13.2892 7.4733 14.3325

(13.2330, 13.3453) (7.47, 7.4766) (14.2993, 14.3656)

13.3376 7.4773 14.3474

(13.1785, 13.4001) (7.4669, 7.4802) (14.2669, 14.3982)

Mustard, 2004), but for our work there were insufficient in situ observations to build a correction factor for each image. The relationship between WST and in situ water temperature was slightly steeper than

that between brightness temperature and in situ water temperature, which demonstrates that the atmospheric correction methods employed had the desired effect of accounting for the imperfect atmospheric transmission of water leaving radiance. Further, this correction factor was greater at higher temperatures expected during the humid summer months. An inspection of Fig. 3 shows that 5–10 images reported anomalously elevated WST above in situ observations, which might be due to sea state conditions favoring elevated WST during the warmest days of the study period. Across all image dates (with the exception of the 5–10 warmest observations just mentioned), we did not see patterns in accuracy associated with season or temperature, further supporting the use of these data for building WST climatologies.

Fig. 6. Maximum and minimum water surface temperature and DOY of maximum and minimum water surface temperature.

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

4.2. Results of the curve fitting

341

Table 5 Water surface temperature of samples around the power plants (Unit: °C).

We successfully fit the ATC model to every pixel in the analysis area and calculated extreme temperatures for each pixel. Visually, the ATC model was an appropriate approximation to the trend of WST, and the analysis area mean Root Mean Square Error (RMSE), which indicates how well data points fit the ATC. Each of the parameters (a, b, and c) provide information regarding the pattern of WST in Chesapeake Bay (Fig. 5). The range of each parameter differed, depending on its sensitivity to geographic variation in temperature (Table 3). For example, the range of parameters a and c was larger than b, presumably owing to larger geographic variation in mean, min, and max temperature than in the timing of minimum and maximum temperatures each year (i.e. the phase of the sine curve.)

Power stations

Dominion Virginia Power Morgantown Generating Station Calvert Cliffs Nuclear Power Plant Chalk Point Generating Station Brandon Shores Generating Station/Herbert A. Wagner Generating Station Baltimore Gas & Electric Co Constellation Energy Commodities Group

Near the station

Medium to station

Far to station

Max T

Min T

Max T

Min T

Max T

Min T

28.73 29.25 27.40 29.74 30.09

2.08 3.45 2.41 3.10 4.35

27.90 28.06 26.76 28.85 28.91

1.39 2.63 2.40 2.59 3.31

27.90 27.58 26.44 28.66 27.38

1.29 2.22 2.39 2.04 1.30

32.42 2.89 29.62 2.71 28.27 1.67 32.34 3.04 29.17 2.59 28.33 0.64

4.3. Uncertainty analysis of the fitting parameters 4.4. Maximum, minimum temperature and dates Bootstrapped estimates of uncertainty (95% confidence bounds) for the fitting parameters were consistent with the parameters derived from curve fitting model (Table 4). This indicated that there was low uncertainty in the estimated parameters, and the fitting results induced from a robust curve fitting computation could be used for further analysis of the entire research site. The fitting parameters estimated from bootstrapped re-sampling were normally distributed, suggesting the standard deviation of these distributions can be used as an uncertainty estimate for each parameter. Linear correlation analysis among the parameters indicated only week colinearity between the parameters; the correlation coefficients between a and b, a and c, b and c were 0.402, 0.534 and 0.559, respectively. The bootstrapped parameters provided only small deviations in the shape and amplitude of the ATC model compared with parameters derived from the original WST data. This further verified that uncertainty of the fitting parameters was quantifiable and small relative to the geographic variation.

Maximum average temperatures (1984–2011) were greatest in the shallow tributaries where fresh water inputs from rivers dominate the WST signal (Fig. 6). From these narrow and shallow, yet still tidally influenced, water bodies, maximum water temperatures decreased with distance along the tributaries and into the main stem of the Bay. This gradient ranged from above 30 °C to 25 °C, with the coldest maximum temperature identified in the Bay main stem, along the eastern side of the bay. Likewise, the date of onset of maximum temperature correlated with position along the tributaries, with the very top of each tributary warming earlier than deeper downstream water. Maximum temperatures in the upper tributaries were reached by mid-July (DOY = 193), whereas maximum temperatures in the main stem were not reached until nearly mid-August (DOY = 220). Minimum average temperatures were more spatially variable, especially in the main stem of the bay, and appeared to correlate more strongly with latitude than water depth. The warmest minimum

Fig. 7. Fitting parameters for WST samples in the vicinity of a power station.

342

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Fig. 8. Air and water surface temperature trends and linear regression analysis air temperature and water surface temperature at eight representative locations.

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Fig. 8 (continued).

343

344

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Fig. 8 (continued).

temperatures were found in the main stem of the bay, near the confluence with the Potomac River. Here, minimum water temperatures in winter were between 2 °C and 9 °C. The coldest minimum temperatures were at the top of the bay near the confluence with the Susquehanna River. As with maximum temperatures, the date of onset of minimum temperatures was generally earlier at the top of the tributaries than in the main stem of the Bay, but this effect was much stronger in tributaries with large urban centers (e.g., the Potomac and Patapsco rivers). 4.5. Thermal influence of power plants We studied the impact on WST of eight power plants along the Chesapeake Bay shoreline that use Bay water for cooling. We extracted samples at increasing distance around each power plant to give a visual and quantitative demonstration of the impact of thermal effluent on WST. A clear gradient with distance from the power generation facility was found at most plants, with water close to the facility exhibiting an elevated water temperature relative to water more distant from the facility (Fig. 7). For example, at Brandon Shores Generating Station and Herbert A. Wagner Generating Station, the mean temperature (parameter c) spanned 5 °C; WST along a transect from the power station to the

middle of the Bay were 19.48 °C, 17.15 °C, 15.45 °C, 14.52 °C and 13.89 °C. Although a similar pattern was found at other power stations, there was variability, and some power stations did not show this spatial pattern (Table 5), presumably due to other larger influences on WST in the area. Maximum and minimum WST around each power plant was extracted and tabulated (Table 5). Locations proximate to each power plant typically exhibited higher minimum temperatures compared with locations further from the power plants. This was most noticeable along the Potomac River estuary, in the inner harbor of Baltimore, and adjacent to Calvert Cliffs nuclear power plant on the western shore of the Bay main stem (Fig. 6). The effect of power plant thermal effluent on maximum water temperatures was not as strong as it was on minimum water temperatures, but it was still measurable for several of the power plants. 4.6. Linear regression between WST, air temperature and year Eight weather stations where air temperature was recorded over the past 28 years were used to compare water and air temperatures (Figs. 1&9). At each station, the mean air temperature was calculated

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

at five-year intervals, which was consistent with the interval of annual mean WST calculations. In general, all sites exhibited increasing air and water temperatures over the analysis period (Fig. 8). However, the slope and significance of the linear regression between WST and air temperature varied at different locations (Table 6; Fig. 8), indicating spatial variability in the coupling of air and water temperatures. The largest and the smallest effect of air temperature on water temperature were found at Fredericksburg station (slope: 2.01) and Annapolis station (slope: 0.43), respectively. Because air and water temperature were both increasing through the analysis period, larger slopes (effects) indicate that water temperature increased more rapidly than air temperature. Increases in WST over the past 28 years were found at all stations, but a more complete census of the Bay showed that the increase in water temperature was weakest at the top of the bay near the confluence with the Susquehanna river. The trend of increasing water temperature over time strengthened towards the main stem of the bay, and therefore strengthened with greater influence of oceanic water (Fig. 9). However, the spatial patterns also show locally more rapid water temperature warming in the Patapsco, Patuxent, and Potomac river estuaries, suggesting that these western tributaries are generally warming faster than the eastern tributaries and potions of the Bay main stem. 5. Discussion 5.1. Curve fitting model performance The extreme temperatures (maximum, minimum), the DOY of their occurrence, and trends over time all required the calculation of the ATC parameters with sufficiently low uncertainty to reveal spatio-temporal trends. Our bootstrapped uncertainty range was + 0.0543 °C/− 0.0073 °C for maximum and minimum temperatures, respectively, which was 1 to 2 orders of magnitude lower than the spatio-temporal variation seen in the data. The variation in parameters a and c (amplitude and mean temperature, respectively) was largest of the three parameters, and thus controls the uncertainty in derived values (for example, max and min temperature). When compared with in situ observations of water temperature, our WST observations from Landsat explained 89% of the variance and

345

resulted in a slope close to 1. Sea state at the time of in situ observation is likely to influence WST comparison with in situ observations, yet we had limited a information to use to assess the strength of this effect or correct for its influence. However, the logical approach would be a deeper analysis of the many types of existing in situ water temperature observations. Despite wide variation in the time and method of acquisition it might be possible to bring some light to this issue. 5.2. Spatial pattern of WST climatology There was small variation in the maximum temperature across the main stem of the Chesapeake Bay, ranging from 24 °C to 26 °C. However, there was larger a difference in maximum WST between the tributaries and the main stem of the Bay. This was true for most tributaries, however, the Patapsco river was a moderate exception. The Patapsco river receives the majority of runoff from the highly urbanized area of Baltimore, MD, and maximum water temperature here was slightly cooler than tributaries located just a short distance to the north or south. Therefore, while we must conclude that water depth and water contributions from rivers are the primary influences on the maximum water temperature, the cooler summer water temperatures of the Patapsco river hints that increased runoff generated by impervious surfaces might actually lead to cooler WST. Certainly greater in-depth analysis is required to clarify this response. A clear gradient in minimum temperature was also found (Fig. 6) that ran from north (low minimum temperature) to south (higher minimum temperature). This gradient was not strongly correlated with distance up each tributary (as was found for maximum temperature) and therefore is not as strongly associated with water depth. Instead, distance to the open ocean is the primary control on minimum water temperature. Deviations from this gradient were seen in the tributaries with a major urban center, which likely exhibited higher than average minimum WST presumably due to warmer runoff from impervious surfaces in winter (Carlson & Arthur, 2000). In the most heavily urbanized tributary (the Patapsco river) it indeed seems to be the case that summer maximum temperatures are lower and winter minimum temperatures are higher than near-by less-urbanized tributaries. Dates of maximum temperature demonstrated that most of the water in the Bay attained a maximum temperature between DOY 213 and 219

Table 6 Summary of the results from linear regression analysis of different locations. Stations

Linear regression

Slope

95% confidence interval

Intercept

95% confidence interval

Mean square error of residual

Princess

WST ~ Air Air ~ year WST ~ year WST ~ Air Air ~ year WST ~ year WST ~ Air Air ~ year WST ~ year WST ~ Air Air ~ year WST ~ year WST ~ Air Air ~ year WST ~ year WST ~ Air Air ~ year WST ~ year WST ~ Air Air ~ year WST ~ year WST ~ Air Air ~ year WST ~ year

1.56 ± 0.25 0.024 ± 0.005 0.059 ± 0.007 1.94 ± 0.29 0.031 ± 0.005 0.072 ± 0.012 2. ± 0.324 0.028 ± 0.005 0.076 ± 0.012 1.74 ± 0.41 0.02 ± 0.005 0.075 ± 0.009 1.43 ± 0.26 0.04 ± 0.006 0.08 ± 0.009 0.43 ± 0.21 0.08 ± 0.013 0.054 ± 0.018 1.54 ± 0.26 0.03 ± 0.005 0.063 ± 0.010 1.02 ± 0.30 0.02 ± 0.005 0.037 ± 0.01

(1.03,2.11) (0.014,0.035) (0.045,0.073) (1.33,2.55) (0.022,0.042) (0.047,0.098) (1.334,2.68) (0.017,0.04) (0.051,0.102) (0.89,2.59) (0.011,0.033) (0.057,0.093) (0.92,2.04) (0.025,0.05) (0.061,0.101) (−0.015,0.88) (0.05,0.11) (0.015,0.094) (1.0,2.08) (0.02, 0.04) (0.043,0.084) (0.39,1.65) (0.01,0.3) (0.017,0.058)

−5.98 ± 3.28 −36.3 ± 9.79 −103.95 ± 13.2 −13.63 ± 4.14 −49.17 ± 9.68 −131.3 ± 24.7 −12.60 ± 4.37 −43.54 ± 10.71 −138. ± 24.66 −10.50 ± 5.72 −30.04 ± 10.41 −135.8 ± 17.5 −9.18 ± 4.23 −59.05 ± 12.35 −147.3 ± 18.89 7.20 ± 2.96 −154.6 ± 25.19 −95.61 ± 36.82 −7.34 ± 3.35 −50.39 ± 10.38 −114.8 ± 19.6 0.28 ± 3.81 −34.01 ± 10.63 −61.95 ± 19.85

(−12.80,0.82) (−56.6, −-16.0) (−131., −76.4) (−22.2, −5.04) (−69.2, −29.1) (−182., −80.1) (−21.6, −3.54) (−65.8, −21.3) (−189., −86.8) (−22.3, 1.36) (−51.6, −8.43) (−172., −99.5) (−17.96, −0.4) (−84.6, −33.4) (−186., −108.) (0.89,13.51) (−208, −100.9) (−174., −17.1) (−14.3, −0.38) (−71.9, −28.8) (−155., −74.1) (−7.62,8.20) (−56.0, −11.96) (−103., −20.7)

0.08 0.027 0.050 0.15 0.027 0.176 0.17 0.033 0.175 0.21 0.031 0.08 0.187 0.044 0.103 0.1720 0.065 0.139 0.125 0.03 0.111 0.124 0.03 0.113

Royal

Fredericksburg

Cambridge

Baltimore

Annapolis

Chestertown

Millington

346

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Fig. 9. Spatial patterns in the slope of a regression line between WST and Year for each pixel region.

(i.e., the begin of August). Dates of minimum temperature were attained at the end of January to February. Unlike temperature magnitude, the dates WST reached minimum and maximum temperatures showed very similar spatial patterns. The main stem always lagged the tributaries by one or two weeks. This effect was controlled by the phase parameter, which influences a shift of the sine function in the ATC. Pixels located near the shore exhibited a greater phase parameter compared with pixels further off shore (Fig. 10). For deeper water off the shore, more time with an air-water temperature differential is required to attain the same water temperature compared with more shallow near-shore water. Similarly, WST began to decrease in autumn near shore first, followed by deeper offshore water. While these patterns do not necessarily reveal new understanding of the physical oceanography in the Chesapeake Bay, the fact they are consistent with what we would expect lays credence to the methods used here and the usefulness of the data for further research (such as analysis of trends over time discussed in the next section). The influence of cities was most apparent in the record of minimum water temperature in the shallow tributaries adjacent to the largest urban areas in the region. We also saw spatially elevated WST increases

in mean and minimum water temperatures (primarily minimum temperatures) in the vicinity of some power plants. While this effect was not as strong as has been observed in other locations along the U.S. east coast (Fisher & Mustard, 2004), it was a clear pattern in some areas (Figs. 6 & 7). This suggests that future work looking to quantify the impact of thermal effluent pollution on the Bay could benefit from using a Landsat time series approach. 5.3. WST and air temperature trends Water surface and air temperatures increased from 1984 to 2007 at all temperature stations analyzed, but the strength of this increase varied by location (Fig. 8). Further, there was variation in the relationship between mean air temperature and mean WST making it difficult to assign causality (Fig. 8, Table 4). The slope between water temperature (y-axis) and air temperature (x-axis) ranged from 2.01 to 0.43, with values greater than 1 indicating a more rapid increase in water temperature than air temperature. For example, the Cambridge and Fredericksburg sites had greater slopes (1.74 and 2.01, respectively),

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

347

such as urbanization and power generation that has produced water with an elevated temperature (particularly in winter) in the regions of cities located at the top of tributaries to the Chesapeake Bay; and (2) the influence of ocean water on the water temperature, particularly in the main stem of the Bay. Future research in this area would benefit from integrating relatively high-resolution observations of WST with models to attain the goal of mechanistically linking the various drivers of WST change with observations.

Acknowledgement This project was supported by Chinese Scholarship Council (File No. 2011832037). We thank Steven Guinn for providing Python code used in some of the remote sensing image processing, and Hongsheng Bi and Sean McGuire for providing in situ temperature data. The paper also benefited from conversations with Laura Harris.

References Fig. 10. Annual temperature cycle differences for pixels at different distances from shore.

demonstrating that WST has increased faster than air temperatures at these locations. For Princess location, which is located in the middle of the main stem of the Bay, the temperature of ocean water has a larger impact on WST than it would at any of the other sites studied. The influence of ocean water at this site has a moderating effect on the WST climatology (i.e., this site exhibits the coolest maximum temperature and warmest minimum temperature), but is also associated with a more rapid warming since 1984. The Fredericksburg site, which is located near the transition from fresh to saline estuarine water in the Potomac River estuary, exhibits a similarly faster rate of increase in WST relative to air temperature. While this might also be due to increasing influence of warmer ocean water over time, it is also possible that impervious surfaces and industrial activities such as power generation in the large Washington DC metropolitan area have amplified the impact of increasing air temperatures on WST. The air and water temperature records from the other tributary with a large city (Baltimore) also exhibited a rapid increase in WST since 1984, but here air temperature has also increased rapidly. More data on the history of urbanization, distance to the ocean or tidal variation, water quality, water turbidity, and water depth are needed to assign causality to water temperature variation at different regions. However, relationships between air temperature, water temperature and year demonstrated that the eight stations studied here exhibited an obvious increase in WST over the past 28 years, motivating future research into drivers. 6. Conclusions Water temperature is one of the most important controls on the functioning of aquatic ecosystems, motivating the study of spatiotemporal patterns in water surface temperature. Based on the Landsat thermal images and applying a curve fitting model to all available WST data organized by DOY, our work enables the study of broad spatial patterns in WST climatologies as well as trends over time. The curve fitting approach resulted in model output that was free of disturbances caused by clouds and other effects that influence the measurement of WST in any given year or location, thus making general inference more attainable. Increasing trends in air and water temperature were found at all stations studied here. Because WST data were derived from spatially continuous remote sensing imagery, a map of the rate in WST increase was produced (Fig. 9), revealing interesting spatial patterns in this rate. While it is difficult to assign causality to increases in WST, the spatial patterns supported a mix of drivers including (1) human activity

Alcântara, E.H., Stech, J.L., Lorenzzetti, J.A., Bonnet, M.P., Casamitjana, X., & Assireu, A.T. (2010). Remote sensing of water surface temperature and heat flux over a tropical hydroelectric reservoir. Remote Sensing of Environment, 114, 2651–2665. Carlson, T.N., & Arthur, S.T. (2000). The impact of land use land cover changes due to urbanization on surface microclimate and hydrology: A satellite perspective. Global and Planetary Change, 24, 49–65. Cowan, J.L.W., & Boynton, W.R. (1996). Sediment-water oxygen and nutrient exchanges along the longitudinal axis of Chesapeake Bay: Seasonal patterns, controlling factors and ecological significance. Estuaries, 19, 562–580. Demars, B.O., & Manson, J.R. (2013). Temperature dependence of stream aeration coefficients and the effect of water turbulence: A critical review. Water Research, 47, 1–15. Elmore, A.J., & Guinn, S.M. (2010). Synergistic use of Landsat Multispectral Scanner with GIRAS land-cover data to retrieve impervious surface area for the Potomac River Basin in 1975. Remote Sensing of Environment, 114, 2384–2391. Fisher, J.I., & Mustard, J.F. (2004). High spatial resolution sea surface climatology from Landsat thermal infrared data. Remote Sensing of Environment, 90, 293–307. Gergory, J.S., Beesley, S.S., & Kirk, R.W.V. (2000). Effect of springtime water temperature on the time of emergence and size of Pteronarcys californica in the Henry's Fork catchment, Idaho, U.S.A. Freshwater Biology, 45, 75–83. Handcock, R.N., Gillespie, A.R., Cherkauer, K.A., Kay, J.E., Burges, S.J., & Kampf, S.K. (2006). Accuracy and uncertainty of thermal-infrared remote sensing of stream temperatures at multiple spatial scales. Remote Sensing of Environment, 100, 427–440. Hedger, R.H., Malthus, T.J., Folkard, A.M., & Atkinson, P.M. (2007). Spatial dynamics of estuarine water surface temperature from airborne remote sensing. Estuarine, Coastal and Shelf Science, 71, 608–615. Hette-Tronquart, N., Roussel, J.M., Dumont, B., Archaimbault, V., & Pont, D. (2013). Variability of water temperature may influence food-chain length in temperate streams. Hydrobiologia, 718, 159–172. Hook, S.J., Chander, G., Barsi, J.A., Alley, R.E., Abtahi, A., & Palluconi, F.D. (2004). In-flight validation and recovery of water surface temperature with Landsat-5 thermal infrared data using an automated high-altitude lake validation site at Lake Tahoe. IEEE Transactions on Geoscience and Remote Sensing, 42, 2767–2776. Jacobs, A.F., Heusinkveld, B.G., Kraai, A., & Paaijmans, K.P. (2008). Diurnal temperature fluctuations in an artificial small shallow water body. International Journal of Biometeorology, 52, 271–280. Jimenez-Munoz, J.C., Cristobal, J., Sobrino, J.A., Soria, G., Ninyerola, M., & Pons, X. (2009). Revision of the single-channel algorithm for land surface temperature retrieval from Landsat Thermal-Infrared data. IEEE Transactions on Geoscience and Remote Sensing, 47, 339–349. Kaushal, S.S., Likens, G.E., Jaworski, N.A., Pace, M.L., Sides, A.M., & Seekell, D. (2010). Rising stream and river temperatures in the United States. Frontiers in Ecology and the Environment, 8(9), 461–466. Lamaro, A.A., Mariñelarena, A., Torrusio, A.E., & Sala, S.E. (2013). Water surface temperature estimation from Landsat 7 ETM+ thermal infrared data using the generalized single-channel method: Case study of Embalse del Río Tercero (Córdoba, Argentina). Advances in Space Research, 51, 492–500. Lisi, P.J., Schindler, D.E., Bentley, K.T., & Pess, G.R. (2013). Association between geomorphic attributes of watersheds, water temperature, and salmon spawn timing in Alaskan streams. Geomorphology, 185, 78–86. Lomas, M.W., & Glibert, P.M. (1999). Temperature regulation of nitrate uptake: A novel hypothesis about nitrate uptake and reduction in cool-water diatoms. Limnology and Oceanography, 44, 556–572. Masek, J.G., Huang, C.Q., Wolfe, R., Cohen, W., Hall, F., Kutler, J., et al. (2008). North American forest disturbance mapped from a decadal Landsat record. Remote Sensing of Environment, 112, 2914–2926. Masek, J.G., Vermote, E.F., Saleous, N., Wolfe, R., Hall, F.G., Huemmrich, F., et al. (2006). A Landsat surface reflectance data set for North America, 1990–2000. Geoscience and Remote Sensing Letters, 3, 68–72. Miller, W.D., & Harding, L.W. (2007). Climate forcing of the spring bloom in Chesapeake Bay. Marine Ecology-Progress Series, 331, 11–22.

348

H. Ding, A.J. Elmore / Remote Sensing of Environment 168 (2015) 335–348

Najjar, R.G., Pyke, C.R., Adams, M.B., Breitburg, D., Hershner, C., & Kemp, M. (2010). Potential climate-change impacts on the Chesapeake Bay. Estuarine, Coastal and Shelf Science, 86, 1–20. Paaijmans, K.P., Takken, W., Githeko, A.K., & Jacobs, A.F. (2008). The effect of water turbidity on the near-surface water temperature of larval habitats of the malaria mosquito Anopheles gambiae. International Journal of Biometeorology, 52, 747–753. Politi, E., Cutler, M.E.J., & Rowan, J.S. (2012). Using the NOAA Advanced Very High Resolution Radiometer to characterise temporal and spatial trends in water temperature of large European lakes. Remote Sensing of Environment, 126, 1–11. Preston, B.L. (2004). Observed winter warming of the Chesapeake Bay Estuary (1949–2002): Implications for ecosystem management. Environmental Assessment, 34, 125–139. Scheitlin, K. (2013). The Maritime influence on diurnal temperature range in the Chesapeake Bay Area. Earth Interactions, 17, 1–14.

Sobrino, J.A., Jimenez-Munoz, J.C., & Paolini, L. (2004). Land surface temperature retrieval from Landsat TM 5. Remote Sensing of Environment, 90, 434–440. Takeshige, A., Takahashi, T., Nakata, H., & Kimura, S. (2013). Long-term trends in sea surface temperature in coastal water in relation to large-scale climate change: A case study in Omura Bay, Japan. Continental Shelf Research, 66, 73–82. Tan, J., & Cherkauer, K.A. (2013). Assessing stream temperature variation in the Pacific Northwest using airborne thermal infrared remote sensing. Journal of Environmental Management, 115, 206–216. Weng, Q., & Fu, P. (2014). Modeling annual parameters of clear-sky land surface temperature variations and evaluating the impact of cloud cover using time series of Landsat TIR data. Remote Sensing of Environment, 140, 267–278. Wloczyk, C., Richter, R., Borg, E., & Neubert, W. (2006). Sea and lake surface temperature retrieval from Landsat thermal data in Northern Germany. International Journal of Remote Sensing, 27, 2489–2502.