Actual evapotranspiration (ET) is perhaps the most difficult quantity to directly measure among the major water balance components. Because of the high cost and labor constraints associated with the direct measurement of ET, empirical data-driven modeling has frequently been used to estimate ET. Beyond the widely used traditional type regression that has the effect of producing ‘global’ parameter estimates, assumed to be uniform throughout an area, we utilized a more localized spatially non-stationary technique – the geographically weighted regression (GWR) – to estimate mean monthly ET in the Passaic River Basin (PRB). We identified the key environmental controls of ET and developed new sets of spatially varying empirical ET models based on variable combinations that produced the best-fit model. The analysis showed that temporal and spatial variabilities in ET over the PRB are driven by climatic and biophysical factors. We found that the key controlling factors were different from month to month, with wind speed being dominant throughout the year in the study basin. A monthly mean ET index map was further generated from the model to illustrate areas where ET exceeds precipitation. This will among others enable water loss due to evapotranspiration to be accounted for in future water supply plans for the basin.

  • New set of empirical monthly actual ET models was developed.

  • Use of spatially varying regression technique removes frequently overlooked stationarity assumption of ordinary least square models (OLS).

  • The GWR technique reveals any rapid changes in environmental variables over an area, necessitating further investigation.

  • Readily available data can be used to reasonably quantify monthly ET under mean climatic conditions.

Graphical Abstract

Graphical Abstract
Graphical Abstract

Actual evapotranspiration (hereafter referred to as ET), involving soil–water evaporation (E) and vegetation transpiration (T), plays an integral role in the transfer of water and energy within the hydrologic cycle (Gowda et al. 2007). It accounts for more than 60% of precipitation input on a global scale (Ma & Szilagyi 2019), and about 50–70% in the United States (Sun et al. 2002; Brooks et al. 2012). Apart from being used as an index for climate change, location-specific data on ET has practical significance in water resources management planning and monitoring of hydrologic systems (e.g. storage changes) in river basins (Verstraeten et al. 2008; Senay et al. 2016). However, accurate quantification of this important variable at the watershed scale is often costly and time-consuming (Gasca-Tucker et al. 2007), with underlying uncertainties due to the complex array of ET processes (Li et al. 2014), and data limitation. Because it occurs in the gaseous state, unlike precipitation and streamflow, ET is also the most difficult quantity to directly measure among the major water balance components.

The literature is replete with different techniques for quantifying ET at various spatial and temporal scales. Generally, direct measurement of ET includes the lysimeter method, water balance approach, sapflow, eddy covariance (EC), Bowen Ratio methods, and stable isotope techniques (Lu et al. 2003; Williams et al. 2004; Liu et al. 2013; Sanford & Selnick 2013; Gebler et al. 2015). While it is acknowledged that each method exhibits some inherent limitations, the eddy covariance technique is by far the most accurate (Fang et al. 2016). It continuously measures site-level fluxes to generate data series at high temporal resolution. On the other hand, recent advances in optical remote sensing (RS) have also achieved sufficient accuracy in estimating ET and are widely used at large spatial scales (e.g. Zhang et al. 2016; Reitz et al. 2017; Chen & Liu 2020; Ma & Zhang 2022). In estimating ET, remote sensing models employ surface energy balance formulations to partition incident solar radiation into soil heat flux, sensible heat flux, and latent heat flux. At a relatively high spatial resolution (∼4 km) covering broader temporal record, ET has been estimated with the TerraClimate dataset (Abatzoglou et al. 2018). The water balance-based ET (WBET) product from TerraClimate has received wide acceptance in the remote sensing community due to its demonstrated skill in capturing hydroclimatic variables across different regions (e.g. Hu & Hu 2019; Salhi et al. 2019; Xu et al. 2019; Zhao & Gao 2019). It was developed by combining high spatial resolution climatological normals of temperature and precipitation from the WorldClim dataset with time-varying datasets from the Climate Research Unit data (CRU Ts4.0) and Japanese 55-year Reanalysis (JRA-55) (Abatzoglou et al. 2018). The performance of resulting RS – WBET model products are often evaluated using EC flux observations (e.g. Ruhoff et al. 2013; Fang et al. 2016). However, Kalma et al. (2008) after reviewing 30 published validation studies underscored that more sophisticated physical and analytical methods do not necessarily perform better in estimating ET than empirical and statistical techniques.

Because of the high cost and labor constraints associated with the direct measurement of ET (e.g. eddy covariance method), coupled with its limited coverage, empirical data-driven modeling has largely been used to estimate ET (e.g. Lu et al. 2003; Sanford & Selnick 2013; Valipour 2015; Fang et al. 2016). These studies cover multiple scales, use different data types, and address various research questions. Data from such studies are typically drawn from specific geographical units where a single regression is estimated for ET based on a relationship between one or more controlling factors and a dependent variable. For example, Valipour (2015) used a linear regression technique to estimate evapotranspiration from 11 temperature-based models in 31 provinces of Iran and revealed that the best model for estimating evapotranspiration performed well (R2 values > 0.99) in only 11 out of the 31 provinces. While this traditional type of regression has seen wide utility, it has the effect of producing ‘average’ or ‘global’ parameter estimates which are always assumed to be uniform over the study area – an implicit assumption that has frequently been overlooked. This assumption is inherently deficient particularly when applied to environmental variables such as ET which are spatially non-stationary over a large area coverage. Fotheringham et al. (2003) noted that relationships which are not stationary, when applied in a conventional regression model, create problems for the interpretation of estimated parameters. As such, a more localized approach (i.e. geographically weighted regression (GWR)), may be an ideal alternative. GWR belongs to the family of local statistics, comprised of multi-valued estimates as opposed to global statistics. As location changes, local statistics can take on different values. The strength of GWR lies in its ability to explain spatially varying relationships by essentially allowing model parameters to vary over space. Within a highly heterogeneous space, this type of regression will provide an opportunity to make significant progress toward understanding and predicting patterns in response variables on the basis of influential environmental factors.

One terrain inviting, and perhaps requiring a more localized regression (i.e. GWR) approach is the Passaic River Basin (PRB) (Figure 1). The PRB is a highly diverse and complex basin with substantial heterogeneity in land use, soils, geology, reservoirs, vegetation, slope, and topography. It has been distinctly separated into three sections: The Highlands Area, the Central Basin, and the Lower Valley. Trending northeast–southwest in the basin is the mountainous heavily forested Highlands, which is by far the largest and well noted for its pristine and environmental integrity. The Central Basin is described as the hydrologic centerpiece of the PRB (USACE NY District 1987), with its many wetlands serving as a buffer to flood waters generated from the Highlands region. The Lower Valley, as opposed to the Highlands, is a densely populated, highly industrialized urban belt forming the eastern flank of the basin and features the tidal, saline lower reach of the Passaic River. In the New Jersey Water Supply Plan, 2017–2022 (NJDEP 2017), there is the recognition that the biophysical arrangement within New Jersey's five water regions, including the PRB, represents a complex array of competing interests and issues in respect of water use and demand. Consequently, the amount of water loss due to evapotranspiration from reservoirs, vegetation, and soil, is unaccounted for in the water supply plan due to the unavailability of reliable and complete data. Yet ET is considered as the primary determinant of available water resources. It is the water that would otherwise become streamflow if not released into the atmosphere; and this water loss constitutes on average 50% of the approximately 49-inches precipitation that occurs over the PRB (Newcomb et al. 2000).
Figure 1

Physiographic map of the Passaic River Basin.

Figure 1

Physiographic map of the Passaic River Basin.

Close modal

It is therefore crucial that in order to accurately assess water resources availability in the PRB, the factors that control the amount and timing of ET in the basin be fully explored, based on which ET can be modeled and reliably quantified. This has become even more relevant at this time in light of the pressing need to predict future water stress and risks in the context of climate change (e.g. NJDEP 2017). Given the widely diverse and heterogeneous character of the basin, it will be necessary to identify and analyze the role that internal (e.g. biophysical) characteristics and external (e.g. climatic) conditions play in the spatial distribution of actual evapotranspiration within a spatiotemporal framework.

In that regard, the objectives of this study are to (1) use the classical ordinary least square (OLS) method to identify and determine the major internal (i.e. leaf area index (LAI), elevation) and external (i.e. mean temperature, precipitation, dew point, mean vapor pressure deficit (VPD), solar radiation, wind speed) controls on ET at monthly time scale; and (2) employ GWR to develop spatially varying mathematical models from readily available datasets that can be used to estimate monthly ET within the PRB.

Study area

The non-tidal PRB is an oval-shaped area of about 2,135 km2 (824 square miles), of which about 84% is located in New Jersey and the rest in New York State (Figure 1). The surface elevation in the basin ranges from below sea level at 0.2 m (0.66 ft) to 454 m (1,490 ft). As previously mentioned, physiographically, the basin is divided into three main regions: the series of parallel ridges that trend northeast/southwest forming the Highlands; the Central Basin, comprised of large areas of swamps and meadows; and the roughly flat Lower Valley.

According to Paulson et al. (1991), New Jersey is located in a modified continental climate zone (i.e. with hot summer and cold winter). Five different climate zones have been identified: North, Central, Southwest, Pine Barrens, and Coastal zones, with PRB located within the North and Central zones. Temperatures in the North are about 10 °F colder than the coastal zone with above 90 °F commonly observed within the Central zone during warm seasons. Annual precipitation typically ranges from 1,016 to 1,321 mm (40–52 inches), with peak values observed along the coast during summer and winter (Paulson et al. 1991). On average, approximately 1,245 mm (49 inches) of precipitation occurs in the PRB annually (Newcomb et al. 2000). Of this, about 50% is lost to the atmosphere through evapotranspiration, another 5% becomes surface runoff, and the remaining 45% becomes available as recharge to groundwater aquifers (Mitchell 1992 as cited in NJ Watershed Basins, n.d.).

Land use/land cover patterns in the basin are dominated by forest type vegetation. According to the 2011 National Land Cover Dataset (NLCD; https://www.mrlc.gov/viewer/), approximately 43% of the area is forested; 36% developed; 14% Woody Wetlands; 3% open water; and the remaining 4% comprising other land cover types. Land use decisions in the basin continue to encourage flood events, with a direct and indirect consequence on water resources (i.e. in terms of flow, water quality, and water quantity). Driven by socioeconomic and biophysical factors, many communities in the lower urbanized section of the basin are close to being built out. Meanwhile, the Highlands region is currently undergoing suburban and rural development (NY Rapid Assessment Profile 2011). Coupled with historically unprecedented warming projected over New Jersey (State Climate Summaries, https://statesummaries.ncics.org/chapter/nj/), it is likely that climate and land use change will pose a serious threat to the water supply systems in the basin. As climate change continues and is expected to accelerate the hydrological cycle, water resource availability and ecosystem services will be directly affected through alteration in evapotranspiration (ET) processes, and indirectly through vegetation water use.

Environmental factors and data sources

Gridded datasets from Parameter-elevation Regressions on Independent Slopes Model (PRISM, Oregon State University, http://prism.oregonstate.edu), TerraClimate (Abatzoglou et al. 2018) (http://www.climatologylab.org/terraclimate.html), and Moderate Resolution Imaging Spectroradiometer (MODIS) LAI product (LAADS/DAAC, https://ladsweb.modaps.eosdis.nasa.gov/search/) provided the environmental variables necessary to develop the ET models. PRISM has received wide acceptance because of its ability to reasonably reproduce weather patterns over areas with complex topography such as the study region and serves as the official spatial climate datasets of the United States Department of Agriculture (USDA; Daly et al. 2008). As actual observations of ET are not available for the basin, the water balance-based TerraClimate ET was used as a proxy. According to Abatzoglou et al. (2018), the accuracy of TerraClimate has been proven by its strong validation with station-based observations from meteorological networks including the Global Historical Climate Network (GHCN), Snow Telemetry (SNOTEL), and Remote Automatic Weather Stations (RAWS). TerraClimate fields of annual reference evapotranspiration also tracked well with reference evapotranspiration from FLUXNET stations (Abatzoglou et al. 2018). In the assessment of water balance through different sources of precipitation and actual ET datasets, Neto et al. (2022) found that against the other three sources, TerraClimate emerged as the most sensitive to variations in the spatial distribution of precipitation and actual ET variables. This demonstrated skill of the TerraClimate product in capturing both point and spatial observations is also corroborated by Soleimani-Motlagh et al. (2022), Wiwoho & Astuti (2022), and Lemenkova (2022), and therefore gives credence to its use as a proxy in this study.

On the other hand, biophysical variables suggested to sufficiently explain variability in ET included precipitation (ppt, mm), air temperature (temp, °C), dew point (dewpt, °C), solar radiation (srad, MJ m−2), vapor pressure deficit (vpd, 100 Pa), wind speed (ws, m/s), leaf area index (LAI), and surface elevation (elev, m) (Lu et al. 2003; Sanford & Selnick 2013; Reitz et al. 2017). Table 1 provides details of the data source, resolution, and time periods of each environmental variable. Further processing of the datasets was carried out using the raster calculator in ArcGIS 10.2.2 to compute mean monthly values. After converting each variable file from raster to feature class, a spatial join was performed to obtain attributes of all the environmental variables in a single layer. Rather than resampling which tends to compromise data values, a one-to-one join operation based on intersect matching was used in the spatial join tool. This ensured that spatial relationship was still established between records of coarse versus finer resolution datasets. In the end, a single feature layer for each month containing attributes of all dependent and independent variables at 800 m grid was prepared for further analysis in R (R 2020).

Table 1

Environmental variables used in the study

DataAbbreviationSpatial resolutionTime frameSource
Dependent variable     
 Actual evapotranspiration ET 4 km 1981–2010 TerraClimate 
Independent variables     
 Elevation ELEV 800 m 1981–2010 PRISM 
 Leaf area index LAI 1 km 2000–2010 MODIS 
 Precipitation PPT 800 m 1981–2010 PRISM 
 Minimum temperature TEMPmin 800 m 1981–2010 PRISM 
 Maximum temperature TEMPmax 800 m 1981–2010 PRISM 
 Dew point DEWPT 800 m 1981–2010 PRISM 
 Maximum vapor pressure deficit VPDmax 800 m 1981–2010 PRISM 
 Minimum vapor pressure deficit VPDmin 800 m 1981–2010 PRISM 
 Solar radiation SRAD 4 km 1981–2010 TerraClimate 
 Wind speed WS 4 km 1981–2010 TerraClimate 
DataAbbreviationSpatial resolutionTime frameSource
Dependent variable     
 Actual evapotranspiration ET 4 km 1981–2010 TerraClimate 
Independent variables     
 Elevation ELEV 800 m 1981–2010 PRISM 
 Leaf area index LAI 1 km 2000–2010 MODIS 
 Precipitation PPT 800 m 1981–2010 PRISM 
 Minimum temperature TEMPmin 800 m 1981–2010 PRISM 
 Maximum temperature TEMPmax 800 m 1981–2010 PRISM 
 Dew point DEWPT 800 m 1981–2010 PRISM 
 Maximum vapor pressure deficit VPDmax 800 m 1981–2010 PRISM 
 Minimum vapor pressure deficit VPDmin 800 m 1981–2010 PRISM 
 Solar radiation SRAD 4 km 1981–2010 TerraClimate 
 Wind speed WS 4 km 1981–2010 TerraClimate 

Regression diagnostics and evaluation

Diagnosing a regression model involves using a class of techniques for detecting problems raised either by the model or the datasets that may compromise the predictive ability of the model. In arriving at an optimal model for predicting ET, environmental variables were assessed using: (1) Pearson correlation metrics (r), (2) statistical significance, (3) Bayesian Information Criterion (BIC), (4) Mallow's Cp, and (5) adjusted R-square (R2). The presence of multi-collinearity that reveals near-redundancy among explanatory variables was also explored. Although multi-collinearity was detected, the exclusion or otherwise of explanatory variables from the model was guided by our understanding of the key processes that influenced the response variable (i.e. ET) as well as their statistical significance. Thus, because all variables were significant, and in order not to compromise the predictive ability of the model, we refrained from removing any variables just on the bases of multi-collinearity.

Further diagnostic analysis such as the Best Subsets regression technique was performed to provide a reasonable subset of variables to include in the model. Unlike other selection techniques such as the Stepwise, the Best Subsets approach models every possible combination of variables that exist. For this study, the 10 independent variables underwent this technique, and the best significant model was selected according to a reasonable adjusted r-square, a lower BIC, and a minimum Cp. Note that underpinning our objectives was to identify the model with the minimum number of variables (i.e. not more than 5) that potentially provided an optimal model accuracy and whose variables are readily available from standard meteorological monitoring stations and regional remote sensing products.

Regression modeling approach

Classical regression modeling approaches such as OLS ride on the assumption that the statistical properties of variables (e.g. mean, standard deviation, and covariance) are constant over space and time for a given area of interest. As such, parameter estimates of an OLS model are considered global. The downside to this is that important variations in the spatial pattern and relationship between independent and dependent variables are lost. The standard expression for spatially invariant OLS model is given as:
(1)
where represents a proxy ET variable, is the intercept and is the slope coefficients for the independent variables of interest , respectively, that control ET, and is the error term.
In respect of environmental variables that are spatially dynamic, a better approach such as GWR was developed by Fotheringham et al. (2003). Their concept of GWR is based on Tobler's first law in Geography (Tobler 1970) which states that ‘everything is related to everything else, but near things are more related than distant things’. Unlike the traditional OLS regression framework, GWR allows local rather than global parameters to be estimated for each of the predictor variables. So that Equation (1) is rewritten as:
(2)
where are coordinate locations of the point in space, is the local regression coefficient for the independent variable of the point, is the independent variable of the point, and is the error term of the point.
Parameter estimates in a GWR model depend on a spatial weighting function and the selected bandwidth for model calibration. The weighting function determines the weight to assign to a local observation based on its closeness to a sample point. As formulated by Fotheringham et al. (2003), the local regression coefficients for the data point are estimated as:
(3)
where X is the matrices of independent variables, represents an estimate of b in Equation (2) and is an matrix (where n is the sample size) whose off-diagonal elements are zero and whose diagonal elements represent the geographic weighting for each n observed data around the data point, and Y is the vector of the dependent variable.
Here, the weighting scheme, as expressed in Equation (4) below, is applied to all observations in a certain sliding neighborhood around each data point. The size of this sliding neighborhood is the bandwidth and it determines which nearby observation is to be considered when calibrating coefficient for a data point. Bandwidths can either be constant (fixed kernel) or variable (adaptive kernel) depending on the density of data points at a location. Bandwidths are smaller where sample points are dense and larger at sparse locations (Fotheringham et al. 2003; Tu & Xia 2008). A fixed bandwidth was chosen in this study because of the gridded structure of the datasets. Two schemes of weighting function can be used in calculation: Gaussian and the bi-square. In this study, the Gaussian scheme, suited for modeling numerical rather than binary response variables, was chosen. It is expressed as:
(4)
where is the weight observation j exerts at data point i, is the distance between observation i and j, and b is the kernel bandwidth. Weight rapidly approaches zero when the distant observations are greater than the kernel bandwidth (Tu & Xia 2008). The optimal bandwidth was selected automatically with the spgwr R package based on the corrected Akaike information criterion (AICc).

Evaluation of model performance

Performance of the modeled actual ET was evaluated based on two approaches: (1) an external basin-scale validation and (2) an internal metric-based validation. For the external validation, precipitation, discharge, and storage data were utilized to estimate basin-level water balance ET for comparison with geographically weighted ET (GW ET). Total annual precipitation was based on the gridded 4 km PRISM dataset (Daly et al. 1994). Discharge data at the outlet of the study basin (i.e. Passaic River at Little Falls NJ, site 01389500) was obtained from the USGS water data website (USGS, https://waterdata.usgs.gov/nwis/), and change in storage from Gravity Recovery and Climate Experiment (GRACE) data (Watkins et al. 2015; Wiese et al. 2019). Specifically, we used version 2 RL06 grids processed by the National Aeronautics and Space Administration (NASA) Jet Propulsion Laboratory (JPL, https://grace.jpl.nasa.gov/data/get-data/jpl_global_mascons/). Constrained by the available time range of GRACE data as well as the end of our analysis period, we carried out our comparison for 2003–2010. The data were processed and aggregated into annual precipitation, annual discharge, and annual change in storage over the PRB to compute ET using the water balance equation (e.g. Ma & Szilagyi 2019) below:
(5)
where P, Q, and are annual basin precipitation, basin discharge, and change in storage, respectively. Because of the tendency of water balance not closing at monthly scale, we chose to do the validation over the annual scale although our model was developed on a monthly basis (e.g. Ma et al. 2021). Indeed, WBET validation is not new. It has been used extensively to verify remote sensing-based ET at regional and watershed scales (e.g. Zhang et al. 2010; Senay et al. 2011; Velpuri et al. 2013; Senay et al. 2016; Ma et al. 2021). In their study, Senay et al. (2016) used WBET to validate Landsat8-based ET estimates.

For the internal variable evaluation, three metrics were used to assess the performance of the geographically weighted ET model. They are:

  • The local coefficient of determination (local R2), which denotes how much the model explained the variability in the independent variables.

  • The prediction error (PE), which expresses the difference between Terraclimate ET and modeled geographically weighted ET.

  • The corrected Akaike Information Criterion (AICc), which estimates the quality of models relative to each other for a given set of data.

Model comparison analysis

Comparison of annual WBET with aggregated annual GWET is depicted in Figure 2 and Table 2. GWET compared reasonably well with the WBET over the basin for the 2003–2010 validation period. The calculated water balance ET was, however, slightly pronounced as against the modeled ET with a mean bias of approximately 11%. This could be attributed to the numerous open water bodies within the basin. The runoff–rainfall coefficients (Q/P) of less than 0.55 in Table 2 are indicative of water balance closure. Velpuri et al. (2013) found that coefficients more than 0.55 suggest dominant regional groundwater flows, and such conditions could introduce errors in the water balance computation, thereby affecting closure. Similarly, computations should not result in negative water balance (where combined (Q + S) > P) (Table 2).
Table 2

Comparison of water balance-based ET with modeled geographically weighted ET

YearPrecipDischargeGrace storageQ/P(Q + S) > PWB ETGW ET
2003 1,645.21 697.50 7.07 0.42 704.57 940.64 757.66 
2004 1,366.83 455.57 −1.73 0.33 453.84 913.00 796.89 
2005 1,343.87 531.78 0.79 0.40 532.58 811.30 780.92 
2006 1,392.75 531.44 0.91 0.38 532.36 860.39 815.67 
2007 1,383.73 432.44 −2.56 0.31 429.88 953.85 798.63 
2008 1,344.85 477.88 2.31 0.36 480.19 864.66 768.27 
2009 1,281.54 403.08 −0.28 0.31 402.80 878.74 809.98 
2010 1,327.27 491.78 −2.14 0.37 489.64 837.63 770.69 
YearPrecipDischargeGrace storageQ/P(Q + S) > PWB ETGW ET
2003 1,645.21 697.50 7.07 0.42 704.57 940.64 757.66 
2004 1,366.83 455.57 −1.73 0.33 453.84 913.00 796.89 
2005 1,343.87 531.78 0.79 0.40 532.58 811.30 780.92 
2006 1,392.75 531.44 0.91 0.38 532.36 860.39 815.67 
2007 1,383.73 432.44 −2.56 0.31 429.88 953.85 798.63 
2008 1,344.85 477.88 2.31 0.36 480.19 864.66 768.27 
2009 1,281.54 403.08 −0.28 0.31 402.80 878.74 809.98 
2010 1,327.27 491.78 −2.14 0.37 489.64 837.63 770.69 
Figure 2

Comparison of the annual (2003–2010) geographically weighted modeled ET with ET water balance-based ET in the PRB.

Figure 2

Comparison of the annual (2003–2010) geographically weighted modeled ET with ET water balance-based ET in the PRB.

Close modal
Over the entire PRB, the GWR model successfully predicted the magnitude and seasonal ET patterns derived from the observed environmental variables (i.e. PPT, TEMP, DEWPT, SRAD, VPD, WS, LAI, ELEV). Figures 3 and 4 respectively illustrate ET from TerraClimate data and that predicted by the GW regression model over the PRB. Their close semblance demonstrates the clear spatial and seasonal patterns featured in the area. For approximately 95% of the basin's area, the GWR model explained between 52—100% of the variations in monthly ET (Figure 5; Table 3). The local coefficient of determination (local R2) was all statistically significant at p-value < 0.001. These results showed the robustness and reliability of the geographically weighted ET model. There were, however, few locations in the area that showed somewhat large differences between the predicted and TerraClimate ET values. These were especially observed in the months of December, February, and March where predicted ET was either lower or higher than the TerraClimate ET. Overall, the geographically weighted ET model was a good predictor when absolute PE values were less than 0.5 (Figure 6; Table 3).
Table 3

Summary of monthly local coefficient and prediction error from the GWR ET model

JanFebMarAprMayJunJulAugSepOctNovDec
     Local coefficient of determination     
Minimum 0.374 −0.156 0.215 0.074 −0.122 −0.339 0.089 −0.169 0.275 −2.372 −1.590 −0.351 
1st quantile 0.809 0.729 0.717 0.739 0.715 0.662 0.669 0.709 0.767 0.694 0.759 0.770 
Median 0.883 0.858 0.815 0.815 0.846 0.824 0.791 0.864 0.822 0.879 0.911 0.880 
Mean 0.867 0.820 0.793 0.801 0.808 0.778 0.776 0.813 0.810 0.805 0.851 0.846 
Standard Deviation 0.098 0.148 0.124 0.113 0.157 0.180 0.141 0.171 0.100 0.225 0.162 0.143 
3rd quantile 0.947 0.941 0.883 0.883 0.936 0.926 0.900 0.954 0.879 0.966 0.973 0.960 
Maximum 1.000 1.000 0.996 0.998 0.999 1.000 0.999 1.000 0.983 1.000 1.000 1.000 
     Prediction error     
Minimum −0.980 −1.443 −4.550 −1.021 −1.622 −1.087 −1.654 −0.651 −1.959 −1.069 −1.161 −1.444 
1st quantile −0.031 −0.050 −0.235 −0.090 −0.043 −0.028 −0.058 −0.023 −0.094 −0.024 −0.031 −0.042 
Median 0.000 0.000 0.001 −0.001 −0.001 0.000 −0.001 0.000 −0.002 0.000 0.000 0.000 
Mean −0.001 0.000 0.008 0.004 0.000 0.000 0.002 0.000 0.000 −0.001 0.001 −0.002 
Standard Deviation 0.112 0.148 0.803 0.220 0.205 0.121 0.203 0.087 0.238 0.124 0.154 0.224 
3rd quantile 0.034 0.044 0.223 0.093 0.041 0.027 0.056 0.022 0.092 0.020 0.032 0.044 
Maximum 0.903 1.715 6.410 1.376 2.329 1.134 2.117 0.614 1.646 1.367 1.495 1.768 
JanFebMarAprMayJunJulAugSepOctNovDec
     Local coefficient of determination     
Minimum 0.374 −0.156 0.215 0.074 −0.122 −0.339 0.089 −0.169 0.275 −2.372 −1.590 −0.351 
1st quantile 0.809 0.729 0.717 0.739 0.715 0.662 0.669 0.709 0.767 0.694 0.759 0.770 
Median 0.883 0.858 0.815 0.815 0.846 0.824 0.791 0.864 0.822 0.879 0.911 0.880 
Mean 0.867 0.820 0.793 0.801 0.808 0.778 0.776 0.813 0.810 0.805 0.851 0.846 
Standard Deviation 0.098 0.148 0.124 0.113 0.157 0.180 0.141 0.171 0.100 0.225 0.162 0.143 
3rd quantile 0.947 0.941 0.883 0.883 0.936 0.926 0.900 0.954 0.879 0.966 0.973 0.960 
Maximum 1.000 1.000 0.996 0.998 0.999 1.000 0.999 1.000 0.983 1.000 1.000 1.000 
     Prediction error     
Minimum −0.980 −1.443 −4.550 −1.021 −1.622 −1.087 −1.654 −0.651 −1.959 −1.069 −1.161 −1.444 
1st quantile −0.031 −0.050 −0.235 −0.090 −0.043 −0.028 −0.058 −0.023 −0.094 −0.024 −0.031 −0.042 
Median 0.000 0.000 0.001 −0.001 −0.001 0.000 −0.001 0.000 −0.002 0.000 0.000 0.000 
Mean −0.001 0.000 0.008 0.004 0.000 0.000 0.002 0.000 0.000 −0.001 0.001 −0.002 
Standard Deviation 0.112 0.148 0.803 0.220 0.205 0.121 0.203 0.087 0.238 0.124 0.154 0.224 
3rd quantile 0.034 0.044 0.223 0.093 0.041 0.027 0.056 0.022 0.092 0.020 0.032 0.044 
Maximum 0.903 1.715 6.410 1.376 2.329 1.134 2.117 0.614 1.646 1.367 1.495 1.768 
Figure 3

Monthly actual evapotranspiration maps of the Passaic River Basin from the TerraClimate dataset.

Figure 3

Monthly actual evapotranspiration maps of the Passaic River Basin from the TerraClimate dataset.

Close modal
Figure 4

Modeled monthly actual evapotranspiration maps of the Passaic River Basin from geographically weighted regression.

Figure 4

Modeled monthly actual evapotranspiration maps of the Passaic River Basin from geographically weighted regression.

Close modal
Figure 5

Map of local coefficient of determination for GWR over the PRB.

Figure 5

Map of local coefficient of determination for GWR over the PRB.

Close modal
Figure 6

Map of residual error for GWR over the PRB.

Figure 6

Map of residual error for GWR over the PRB.

Close modal

Table 4 shows how the geographically weighted ET model compares with the traditional OLS model in terms of their Akaike Information Criterion (AIC) (Akaike 1974). The AIC determines which model performs best in associating explanatory variables to the dependent variable. A model with smaller AIC value is considered better as it is more likely to minimize information loss in contrast to the ‘true’ model that generates the observation data (Burnham & Anderson 2003). The smaller AIC values obtained by the geographically weighted ET model in this study further corroborate the point that spatially non-stationary local models are an ideal alternative to the global OLS models when explaining spatially varying relationships over the large scale.

Table 4

Monthly comparison of Akaike Information Criterion of ordinary least square and geographically weighted regression models

OLS AICcGW AICc
Jan 4,507 − 898 
Feb 8,642 2,490 
Mar 14,829 10,490 
Apr 6,920 1,545 
May 7,982 1,495 
Jun 5,849 − 1,217 
Jul 8,172 1,307 
Aug 3,504 − 4,059 
Sep 7,629 2,082 
Oct 6,269 − 1,212 
Nov 5,991 − 38 
Dec 9,088 2,776 
OLS AICcGW AICc
Jan 4,507 − 898 
Feb 8,642 2,490 
Mar 14,829 10,490 
Apr 6,920 1,545 
May 7,982 1,495 
Jun 5,849 − 1,217 
Jul 8,172 1,307 
Aug 3,504 − 4,059 
Sep 7,629 2,082 
Oct 6,269 − 1,212 
Nov 5,991 − 38 
Dec 9,088 2,776 

Potential sources of uncertainty associated with the modeled ET can be linked to the corresponding uncertainties inherent in the TerraClimate data, remote sensing LAI data from MODIS, and other meteorological data from PRISM. For instance, any existing noise or errors in the MODIS LAI product will be propagated in the ET mathematical model. Additionally, relatively coarse grid cells from the MODIS LAI (1 km) as well as TerraClimate data (4 km) were applied to an 800 m grid for the study basin to derive the ET model. At such smaller resolution, the MODIS LAI, for instance, may not sufficiently capture sub-grid scale vegetation signals in the basin, especially where the topography and land cover are highly complex and spatially heterogeneous, and thus introduce some errors in the model. Likewise, the coarse grid size of the radiative and aerodynamic variables obtained from the TerraClimate products may present some degree of errors in the derived ET model. In the course of our analysis, it was observed that the accuracy of ET estimates was largely dependent on the meteorological datasets. For the most part, these datasets, obtained from PRISM at a resolution of 800 m, guaranteed a reasonably accurate ET model.

Spatiotemporal patterns of ET in the PRB

The geographically weighted ET model was carried out over the diverse physiographical and biophysical terrain of the PRB from 1981 to 2010 at 800 × 800 m grid resolution. Based on the ET maps (Figures 2 and 3), there is evidence of sharp basin-wide variation and longitudinal gradients of ET, with lower values trending from northwest to higher values in the southeast. Generally, the ET patterns reflect the climatic conditions at the time and varies depending upon land use/land cover type, topography, and availability of water. Lower values are typically observed in the north climate zone where the temperatures are relatively low, with high topographic relief and largely forested vegetation cover. Mean annual ET value of about 748 mm/year is estimated in the upper mountainous region. In contrast, the high ET values are seen in the central climate zone where temperatures are intermediate to high, and in highly developed urban communities. Estimated mean annual ET is approximately 793 mm/year. Although some uncertainties are present and methodologies may differ, the magnitude and spatial patterns of the estimated ET rates over the PRB are consistent with other investigations in the vicinity (e.g. Sumner et al. 2012).

Mean monthly ET rates over the basin from 1981 to 2010, based on the TerraClimate as well as the GW-modeled ET, show distinct seasonal fluctuations (Figures 2 and 3). The temporal patterns of ET and their spatial variability reflect the controlling effects of prevailing climatic conditions as well as the vegetation distribution. In order of magnitude, ET rates in the PRB are generally greater in the summer months (June–July–August), followed by spring (March–April–May), fall (September–October–November), and winter months (December–January–February). In the summer months, particularly in July, the lower southeast half of the basin shows significantly high ET values because of relatively warm temperature and high radiation (Figure 7).
Figure 7

Mean monthly spatial ET rates in the lower and upper halves of the PRB.

Figure 7

Mean monthly spatial ET rates in the lower and upper halves of the PRB.

Close modal

In contrast, low ET rates are observed in the largely forested upper northwest half because of the precipitation deficit. Although vegetation and water bodies typically dominate the upper half, the winter months show significantly low ET values as a result of relatively low mountainous temperatures and sparse deciduous vegetation. Following the trajectory of ET rates in Figure 7 (monthly time series), it can be seen that ET values in all the months are higher in the lower SE half than in the upper NW portion of the basin. The lowest mean ET values are observed in December (Lower SE: 5.65 mm; Upper NW: 2.53 mm) and the highest values in July (Lower SE: 118 mm; Upper NW: 115 mm).

Dominant controls on monthly ET

Pearson correlation coefficient values show that elevation correlates strongly with ET across all seasons (Table 5). This is followed by energy inputs (DEWPT > TMAX > TMIN > SRAD), aerodynamic input (i.e. WS), atmospheric demand (VPD_MAX > VPD_MIN), vegetation biomass (i.e. LAI), and water input (i.e. PPT) in that order. After elevation, the two most important drivers on monthly ET were DEWPT and TMIN in the winter and spring, DEWPT and TMAX in the summer, and SRAD and TMAX in the fall. Interestingly, the strength of the correlation of explanatory variables with ET is relatively weak in the summer (i.e. < 0.5) across the study basin, although these variables are, to a large extent, strongly correlated among themselves. It was also noted that water availability associated with PPT correlated negatively with ET in all seasons except summer.

Table 5

Pearson correlation coefficient between ET and independent environmental variables

AETELEVDEWPTLAIPPTSRADTMAXTMINVPD_MINVPD_MAXWSAETELEVDEWPTLAIPPTSRADTMAXTMINVPD_MINVPD_MAXWS
AET      Winter          Spring     
ELEV − 0.89          − 0.88          
DEWPT 0.89 − 0.93         0.8 − 0.87         
LAI − 0.49 0.43 − 0.42        − 0.5 0.45 − 0.41        
PPT − 0.37 0.47 − 0.4 0.26       − 0.24 0.33 − 0.37 0.1       
SRAD 0.86 − 0.77 0.8 − 0.39 − 0.17      0.44 − 0.48 0.55 − 0.13 − 0.15      
TMAX 0.86 − 0.91 0.88 − 0.38 − 0.3 0.84     0.77 − 0.86 0.84 − 0.37 − 0.17 0.64     
TMIN 0.88 − 0.91 0.96 − 0.45 − 0.48 0.7 0.79    0.82 − 0.82 0.88 − 0.45 − 0.37 0.28 0.7    
VPD_MIN 0.62 − 0.56 0.6 − 0.37 − 0.48 0.3 0.38 0.79   0.39 − 0.31 0.27 − 0.32 − 0.18 − 0.15 0.13 0.66   
VPD_MAX 0.77 − 0.81 0.73 − 0.33 − 0.23 0.78 0.97 0.64 0.23  0.63 − 0.7 0.61 − 0.29 − 0.01 0.65 0.93 0.45 0.03  
WS 0.7 − 0.51 0.51 − 0.4 − 0.35 0.43 0.41 0.65 0.76 0.32 0.69 − 0.54 0.44 − 0.42 − 0.07 0.16 0.42 0.63 0.58 0.32 
AET      Summer          Fall     
ELEV − 0.42          − 0.83          
DEWPT 0.4 − 0.85         0.74 − 0.86         
LAI − 0.32 0.5 − 0.46        − 0.55 0.5 − 0.5        
PPT 0.34 0.22 − 0.2 0.11       − 0.28 0.42 − 0.51 0.28       
SRAD 0.03 − 0.27 0.31 0.02 − 0.22      0.75 − 0.74 0.69 − 0.39 − 0.36      
TMAX 0.39 − 0.91 0.85 − 0.44 − 0.24 0.32     0.8 − 0.86 0.75 − 0.42 − 0.31 0.83     
TMIN 0.36 − 0.8 0.86 − 0.49 − 0.12 0.01 0.78    0.66 − 0.74 0.88 − 0.46 − 0.38 0.4 0.53    
VPD_MIN 0.13 − 0.36 0.29 − 0.33 0.02 − 0.4 0.33 0.73   0.18 − 0.13 0.24 − 0.18 − 0.24 − 0.09 0.66   
VPD_MAX 0.34 − 0.8 0.62 − 0.37 − 0.23 0.33 0.94 0.57 0.23  0.59 − 0.58 0.33 − 0.22 − 0.07 0.69 0.87 0.1 − 0.33  
WS 0.37 − 0.66 0.54 − 0.49 − 0.05 − 0.05 0.64 0.75 0.66 0.56 0.78 − 0.59 0.57 − 0.46 − 0.08 0.43 0.48 0.71 0.54 0.24 
AETELEVDEWPTLAIPPTSRADTMAXTMINVPD_MINVPD_MAXWSAETELEVDEWPTLAIPPTSRADTMAXTMINVPD_MINVPD_MAXWS
AET      Winter          Spring     
ELEV − 0.89          − 0.88          
DEWPT 0.89 − 0.93         0.8 − 0.87         
LAI − 0.49 0.43 − 0.42        − 0.5 0.45 − 0.41        
PPT − 0.37 0.47 − 0.4 0.26       − 0.24 0.33 − 0.37 0.1       
SRAD 0.86 − 0.77 0.8 − 0.39 − 0.17      0.44 − 0.48 0.55 − 0.13 − 0.15      
TMAX 0.86 − 0.91 0.88 − 0.38 − 0.3 0.84     0.77 − 0.86 0.84 − 0.37 − 0.17 0.64     
TMIN 0.88 − 0.91 0.96 − 0.45 − 0.48 0.7 0.79    0.82 − 0.82 0.88 − 0.45 − 0.37 0.28 0.7    
VPD_MIN 0.62 − 0.56 0.6 − 0.37 − 0.48 0.3 0.38 0.79   0.39 − 0.31 0.27 − 0.32 − 0.18 − 0.15 0.13 0.66   
VPD_MAX 0.77 − 0.81 0.73 − 0.33 − 0.23 0.78 0.97 0.64 0.23  0.63 − 0.7 0.61 − 0.29 − 0.01 0.65 0.93 0.45 0.03  
WS 0.7 − 0.51 0.51 − 0.4 − 0.35 0.43 0.41 0.65 0.76 0.32 0.69 − 0.54 0.44 − 0.42 − 0.07 0.16 0.42 0.63 0.58 0.32 
AET      Summer          Fall     
ELEV − 0.42          − 0.83          
DEWPT 0.4 − 0.85         0.74 − 0.86         
LAI − 0.32 0.5 − 0.46        − 0.55 0.5 − 0.5        
PPT 0.34 0.22 − 0.2 0.11       − 0.28 0.42 − 0.51 0.28       
SRAD 0.03 − 0.27 0.31 0.02 − 0.22      0.75 − 0.74 0.69 − 0.39 − 0.36      
TMAX 0.39 − 0.91 0.85 − 0.44 − 0.24 0.32     0.8 − 0.86 0.75 − 0.42 − 0.31 0.83     
TMIN 0.36 − 0.8 0.86 − 0.49 − 0.12 0.01 0.78    0.66 − 0.74 0.88 − 0.46 − 0.38 0.4 0.53    
VPD_MIN 0.13 − 0.36 0.29 − 0.33 0.02 − 0.4 0.33 0.73   0.18 − 0.13 0.24 − 0.18 − 0.24 − 0.09 0.66   
VPD_MAX 0.34 − 0.8 0.62 − 0.37 − 0.23 0.33 0.94 0.57 0.23  0.59 − 0.58 0.33 − 0.22 − 0.07 0.69 0.87 0.1 − 0.33  
WS 0.37 − 0.66 0.54 − 0.49 − 0.05 − 0.05 0.64 0.75 0.66 0.56 0.78 − 0.59 0.57 − 0.46 − 0.08 0.43 0.48 0.71 0.54 0.24 

The Best Subset regression technique afforded the opportunity to select at most five out of the 10 variable combinations that produced the best-fit ET model for each month. As represented in Table 4, the estimated coefficients from OLS regression revealed the most influential variables that control monthly ET. For example, four variables (WS, TMAX, LAI, and ELEV) explained 93% of the variability in ET in November. Wind speed (WS) shows a strong positive influence on ET followed by TMAX, LAI, and ELEV in that order. Indeed, WS appears to be the most dominant ET controlling factor throughout the year, appearing in 10 out of 12 months. Quite notably, details of the relatively weak correlation observed between environmental variables and ET in the summer are revealed clearly in the month of June. While it may appear that only 39% of the variability in ET has been explained by the five variables, it is worth pointing out that all the 10 variables included in the model could only explain approximately 41% of variability in ET for June. This counters what may be construed as the model's inability to predict ET in the summer. It clearly indicates that the model is quite robust and largely represents the key processes that strongly influence monthly ET in the study basin. More so, corresponding Mallow's Cp and BIC were relatively low compared with other model sets. Lower Cp and BIC are indicative of relatively precise and best models, respectively. Thus, the surprisingly large unexplained variance in the relationship for June may be attributed to some inherent data error propagated from the source data.

Geographically weighted monthly ET models

Based on the intercept and coefficient values, a mathematical expression of the different sets of key independent variables controlling monthly ET could be formulated from Table 6. However, such an expression will only result in a basin-wide average value of ET, concealing the complex nuances of the spatially varying environmental factors that influence monthly ET in the basin.

Table 6

Best-fit ET models according to key monthly environmental variables controlling ET

AprMayJun
VariableJanFebMarEstimated CoefficientJulAugSepOctNovDec
(Intercept) −6.0217 *** −14.1613 *** 20.7849 *** 59.2585 *** 85.7725 *** 92.6907 *** 291.3125 *** 49.7324 *** 25.6541 *** 43.5208 *** 9.0569 *** −10.7300 *** 
TMIN 0.3401 *** – – – 0.8166 *** – 7.4594 *** – – – – −1.6440 *** 
TMAX −7.6322 *** −9.9541 *** – 0.5113 *** – – −15.0457 *** 4.4419 *** 1.9057 *** – 1.0366 *** – 
WS 3.4407 *** 5.7222 *** 8.9909 *** 5.5493 *** 4.3073 *** – – 3.9956 *** 7.8044 *** 5.7030 *** 6.7321 *** 8.0780 *** 
LAI – −0.9655 *** – – – −0.0591 *** – −0.0297 *** −0.1012 *** −0.0996 *** −0.2760 *** – 
ELEV – – – −0.0120 *** −0.0079 *** −0.0058 *** – – – −0.0057 *** −0.0054 *** −0.0096 *** 
VPD_MAX 16.4220 *** 20.4178 *** 6.2951 *** – – 0.0872 *** 6.3268 *** −2.1246 *** – 0.3785 *** – 2.6220 *** 
VPD_MIN – – – – −1.1314 *** −0.6552 *** −9.3819 *** – −2.2326 *** −1.8201 *** – – 
DEWPT 5.0595 *** 6.5406 *** 4.0319 *** – – – – −2.5923 *** – – – 4.1380 *** 
PPT – – – – – 0.1202 *** 0.1485 *** – −0.0816 *** – – – 
Cp 78 157 365 204 258 206 650 71 382 325 248 255 
BIC −8,885 −8,254 −6,747 −7,315 −5,465 −1,598 −7,343 −3,581 −7,153 −5,966 −8,687 −8,464 
R-squared 0.93 0.92 0.87 0.89 0.81 0.39 0.89 0.67 0.89 0.84 0.93 0.93 
AprMayJun
VariableJanFebMarEstimated CoefficientJulAugSepOctNovDec
(Intercept) −6.0217 *** −14.1613 *** 20.7849 *** 59.2585 *** 85.7725 *** 92.6907 *** 291.3125 *** 49.7324 *** 25.6541 *** 43.5208 *** 9.0569 *** −10.7300 *** 
TMIN 0.3401 *** – – – 0.8166 *** – 7.4594 *** – – – – −1.6440 *** 
TMAX −7.6322 *** −9.9541 *** – 0.5113 *** – – −15.0457 *** 4.4419 *** 1.9057 *** – 1.0366 *** – 
WS 3.4407 *** 5.7222 *** 8.9909 *** 5.5493 *** 4.3073 *** – – 3.9956 *** 7.8044 *** 5.7030 *** 6.7321 *** 8.0780 *** 
LAI – −0.9655 *** – – – −0.0591 *** – −0.0297 *** −0.1012 *** −0.0996 *** −0.2760 *** – 
ELEV – – – −0.0120 *** −0.0079 *** −0.0058 *** – – – −0.0057 *** −0.0054 *** −0.0096 *** 
VPD_MAX 16.4220 *** 20.4178 *** 6.2951 *** – – 0.0872 *** 6.3268 *** −2.1246 *** – 0.3785 *** – 2.6220 *** 
VPD_MIN – – – – −1.1314 *** −0.6552 *** −9.3819 *** – −2.2326 *** −1.8201 *** – – 
DEWPT 5.0595 *** 6.5406 *** 4.0319 *** – – – – −2.5923 *** – – – 4.1380 *** 
PPT – – – – – 0.1202 *** 0.1485 *** – −0.0816 *** – – – 
Cp 78 157 365 204 258 206 650 71 382 325 248 255 
BIC −8,885 −8,254 −6,747 −7,315 −5,465 −1,598 −7,343 −3,581 −7,153 −5,966 −8,687 −8,464 
R-squared 0.93 0.92 0.87 0.89 0.81 0.39 0.89 0.67 0.89 0.84 0.93 0.93 

Significant codes (p-values): 0, ‘***’ 0.001, ‘**’ 0.01, ‘*’ 0.05, ‘.’ 0.1, ‘ ’ 1.

A new set of mathematical models derived for monthly ET according to the GW regression analysis is presented in Table 7. Rather than the single, fixed value applied uniformly over space by the OLS technique, the GW regression model provides a range of values of intercepts and coefficients that reflects the local variations in environmental variables as they relate to ET in space (Table 7). By this, ET has been spatially mapped to key environmental variables across the PRB for each month, depicting the range of values in both intercept and coefficient. This will provide a reasonable first approximation of evapotranspiration rates and their spatial distribution to aid in the quantification of ET for water resource planning and decision making to serve communities in the largely diverse PRB.

Table 7

Monthly GWR ET models for the PRB using best-fit environmental variables

MonthGeographically weighted model
January  
February  
March  
April  
May  
June  
July  
August  
September  
October  
November  
December  
MonthGeographically weighted model
January  
February  
March  
April  
May  
June  
July  
August  
September  
October  
November  
December  

To avoid being redundant, we only show, in Figure 8, the spatial map of results relating to April as revealed in the geographically weighted ET model. Tabular data for each month that includes coordinates locations can be provided upon request. A visual distribution of the coefficients of respective variables for selected months is also presented in Figure 9. Unlike the OLS model, the results of the GWR technique indicate a major degree of spatial and temporal variation in the relationship between ET and predictor variables in the study basin. As shown in Figure 8, the influence of wind speed on ET is both positive and negative. For the most part, the lower section of the basin showed a more dominant positive effect on ET whereas portions of the upper part showed a negative effect. Positive coefficient at the lower section is intuitive, given its location near the Atlantic Ocean where wind speed is high. The rate of evapotranspiration would typically be greater in very windy areas. The mostly negative coefficient observed in the upper Highlands region is a reflection of the largely forested vegetation in the area and suggests that at a reduced wind speed, ET decreases. Although elevation and maximum temperature feature as key variables in the area, their relatively low coefficients suggest that, unlike wind speed, they are less influential in driving mean long-term ET in the month of April. Thus, unlike the OLS approach, the changing pattern in coefficients for each controlling factor is clearly depicted by the GWR analysis (Figures 8 and 9).
Figure 8

Spatial distribution of intercept and coefficients of key environmental variables predicting ET in April over the PRB.

Figure 8

Spatial distribution of intercept and coefficients of key environmental variables predicting ET in April over the PRB.

Close modal
Figure 9

Violin plot distribution of the coefficients of key environmental variables predicting ET for selected months in the PRB.

Figure 9

Violin plot distribution of the coefficients of key environmental variables predicting ET for selected months in the PRB.

Close modal
Our results generally support the well-known notion that long-term mean ET for a region is controlled by precipitation and potential evapotranspiration (Budyko 1947). However, it appears that ET is insensitive to precipitation for the most part of the year in the study basin. While this could indicate that the PRB is rarely under water stress, the summer months may be an exception (see Table 6). Although our analyses reveal that wind speed is the main driving force behind long-term mean monthly ET, precipitation appears to be the limiting factor in the summer, particularly the month of June. Overall, a combination of biophysical and climatic factors contribute to long-term ET on a monthly scale. More importantly, the spatial heterogeneity that characterizes the PRB brings to bear the complex challenge in appropriately quantifying ET. However, this has been made possible because of the GWR technique adopted in this study. Flowing from our results, a monthly ET index (i.e. ET/PPT) was reproduced to show how much of precipitation is lost to evapotranspiration across the PRB (Figure 10).
Figure 10

Monthly ET index (%) map over the Passaic River Basin based on GWR.

Figure 10

Monthly ET index (%) map over the Passaic River Basin based on GWR.

Close modal

While the GWR technique can be more appropriate than the global regression approaches, it is not without concerns. There are concerns raised over issues such as kernel and bandwidth selection including those shared with conventional regression techniques. However, the key issue here is that the local variation in the relationship between significant predictor variables and ET would have gone unnoticed in the commonly used global, OLS method. The monthly ET model outputs generated by the GW regression analysis offer new insights into the key external and internal controls on ET. Within the PRB, these models can assist in effectively quantifying ET under typical climatic conditions and to a large extent, accurately estimate mean monthly ET at similar locations. Practical models such as this have been widely used to estimate ET given the huge cost involved in measuring ET on large temporal and spatial scales (Sun et al. 2011a, 2011b). Here, we have established a reasonably accurate non-stationary monthly ET model that uses readily available input data over the PRB.

The 30-year continuous ET and other gridded hydroclimatic and biophysical observations from PRISM/TerraClimate/MODIS provided the means to carefully examine and identify the key environmental variables controlling ET at a fine temporal (i.e. monthly) and spatial (i.e. 800 m grid resolution) scale. From this, the best subset variables for each month were used to develop 12 empirical geographically weighted ET models that will potentially estimate monthly ET over the PRB with reasonable accuracy under mean climate conditions. These variables are readily available from any standard meteorological monitoring stations and regional remote sensing products. Key conclusions from our analyses may be summarized as follows:

  • (1)

    Temporal and spatial variabilities in mean monthly ET over the PRB are significantly controlled by climatic (i.e. TEMP, WS, DEWPT, VPD, PPT) and biophysical (i.e. LAI, ELEV) drivers. The analysis revealed that key controlling factors may be different from month to month, with wind speed taking dominance throughout the year in the study basin. Precipitation, while appearing insignificant in the course of the year, appears to be a limiting factor in the summer months.

  • (2)

    The ability to successfully map ET to key environmental variables in such a complex terrain in this study demonstrates the superiority of the GWR over the OLS approach in modeling spatially varying environmental relationships. Given the usually unquestioned limitation of the OLS technique, one further refinement that could be included in studies concerned with spatial analyses should be the use of a more appropriate localized (i.e. geographically weighted) method when attempting to explain spatial relationships.

  • (3)

    Modeled spatially varying monthly ET developed from this study offer convenient and cost effective means to empirically estimate monthly water loss from similar ecosystems. The ET index map generated for the PRB illustrates areas where ET exceeds precipitation, especially in the summer months, and hence is useful for water resource planning and decision making by water managers in the basin. Moreover, reliable quantification of ET has been made possible in the study basin. As such, the amount of water loss due to evapotranspiration can be accounted for in future water supply plans for the basin.

It is hoped that this work, being the first of its kind in the study basin to the best of our knowledge, will form the foundation for future climate impact studies in the basin and the region at large.

All relevant data are included in the paper or its Supplementary Information.

The authors declare there is no conflict.

Abatzoglou
J. T.
,
Dobrowski
S. Z.
,
Parks
S. A.
&
Hegewisch
K. C.
2018
TerraClimate, a high-resolution global dataset of monthly climate and climatic water balance from 1958–2015
.
Scientific Data
5
,
170191
.
Akaike
H.
1974
A new look at the statistical model identification
.
IEEE Transactions on Automatic Control
19
(
6
),
716
723
.
Brooks
K. N.
,
Ffolliott
P. F.
&
Magner
J. A.
2012
Hydrology and the Management of Watersheds
.
John Wiley & Sons, Ames, IA
.
Budyko
M. I.
1947
On the water and heat balances of the land surface
.
Meteorol. Gidrol
.
5
.
Burnham
K. P.
&
Anderson
D. R.
2003
Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach
.
Springer Science & Business Media, Fort Collins, CO
.
Daly
C.
,
Neilson
R. P.
&
Phillips
D. L.
1994
A statistical-topographic model for mapping climatological precipitation over mountainous terrain
.
Journal of Applied Meteorology and Climatology
33
(
2
),
140
158
.
Daly
C.
,
Halbleib
M.
,
Smith
J. I.
,
Gibson
W. P.
,
Doggett
M. K.
,
Taylor
G. H.
,
Curtis
J.
&
Pasteris
P. P.
2008
Physiographically sensitive mapping of climatological temperature and precipitation across the conterminous United States
.
International Journal of Climatology: A Journal of the Royal Meteorological Society
28
(
15
),
2031
2064
.
Fang
Y.
,
Sun
G.
,
Caldwell
P.
,
McNulty
S. G.
,
Noormets
A.
,
Domec
J. C.
,
King
J.
,
Zhang
Z.
,
Zhang
X.
,
Lin
G.
&
Zhou
G.
2016
Monthly land cover-specific evapotranspiration models derived from global eddy flux measurements and remote sensing data
.
Ecohydrology
9
(
2
),
248
266
.
Fotheringham
A. S.
,
Brunsdon
C.
&
Charlton
M.
2003
Geographically Weighted Regression: The Analysis of Spatially Varying Relationships
.
John Wiley & Sons, New York, NY
.
Gasca-Tucker
D. L.
,
Acreman
M. C.
,
Agnew
C. T.
&
Thompson
J. R.
2007
Estimating Evaporation From a Wet Grassland
.
Gebler
S.
,
Franssen
H. H.
,
Pütz
T.
,
Post
H.
,
Schmidt
M.
&
Vereecken
H.
2015
Actual evapotranspiration and precipitation measured by lysimeters: a comparison with eddy covariance and tipping bucket
.
Hydrology and Earth System Sciences
19
(
5
),
2145
.
Gowda
P. H.
,
Chavez
J. L.
,
Colaizzi
P. D.
,
Evett
S. R.
,
Howell
T. A.
&
Tolk
J. A.
2007
Remote sensing based energy balance algorithms for mapping ET: current status and future challenges
.
Transactions of the ASABE
50
(
5
),
1639
1644
.
Kalma
J. D.
,
McVicar
T. R.
&
McCabe
M. F.
2008
Estimating land surface evaporation: a review of methods using remotely sensed surface temperature data
.
Surveys in Geophysics
29
(
4–5
),
421
469
.
Li
X.
,
Liang
S.
,
Yuan
W.
,
Yu
G.
,
Cheng
X.
,
Chen
Y.
, Zhao, T., Feng, J., Ma, Z., Ma, M. & Liu, S.
2014
Estimation of evapotranspiration over the terrestrial ecosystems in China
.
Ecohydrology
7
(
1
),
139
149
.
Lu
J.
,
Sun
G.
,
McNulty
S. G.
&
Amatya
D. M.
2003
Modeling actual evapotranspiration from forested watersheds across the southeastern United States
.
JAWRA Journal of the American Water Resources Association
39
(
4
),
886
896
.
Ma
N.
,
Szilagyi
J.
&
Zhang
Y.
2021
Calibration-free complementary relationship estimates terrestrial evapotranspiration globally
.
Water Resources Research
57
(
9
),
e2021WR029691
.
Neto
A. K.
,
Ribeiro
R. B.
&
Pruski
F. F.
2022
Assessment Water Balance Through Different Sources of Precipitation and Actual Evapotranspiration
.
Newcomb
D. J.
,
Stanuikynas
T. J.
&
Van Abs
D. J.
2000
Setting of the Raritan River Basin: A Technical Report for the Raritan Basin Watershed Management Project
.
PRISM Climate Group, Oregon State University
.
Available from: http://prism.oregonstate.edu (data created 4 Feb 2014, accessed 15 Nov 2020)
.
New Jersey Department of Environmental Protection
2017
New Jersey Water Supply Plan 2017–2022, 484 p. Available from
: http://www.nj.gov/dep/watersupply/wsp.html
(accessed 7 May 2020)
.
Paulson
R. W.
,
Chase
E. B.
,
Roberts
R. S.
&
Moody
D. W.
1991
National Water Summary, 1988-89 – Hydrologic Events and Floods and Droughts, U.S Geological Survey Paper 2375
.
United States Government Printing Office
,
Reston, VA
. https://pubs.er.usgs.gov/publication/wsp2375.
R Core Team
2020
R: A Language and Environment for Statistical Computing
.
R Foundation for Statistical Computing
,
Vienna
,
Austria
.
Reitz
M.
,
Sanford
W. E.
,
Senay
G. B.
&
Cazenas
J.
2017
Annual estimates of recharge, quick-flow runoff, and evapotranspiration for the contiguous US using empirical regression equations
.
JAWRA Journal of the American Water Resources Association
53
(
4
),
961
983
.
Ruhoff
A. L.
,
Paz
A. R.
,
Aragao
L. E. O. C.
,
Mu
Q.
,
Malhi
Y.
,
Collischonn
W.
, Rocha, H. R. D. & Running, S. W.
2013
Assessment of the MODIS global evapotranspiration algorithm using eddy covariance measurements and hydrological modelling in the Rio Grande basin
.
Hydrological Sciences Journal
58
(
8
),
1658
1676
.
Salhi
A.
,
Martin-Vide
J.
,
Benhamrouche
A.
,
Benabdelouahab
S.
,
Himi
M.
,
Benabdelouahab
T.
&
Ponsati
A. C.
2019
Rainfall distribution and trends of the daily precipitation concentration index in northern Morocco: a need for an adaptive environmental policy
.
SN Applied Sciences
1
(
3
),
277
.
Sanford
W. E.
&
Selnick
D. L.
2013
Estimation of evapotranspiration across the conterminous United States using a regression with climate and land-cover data
.
JAWRA Journal of the American Water Resources Association
49
(
1
),
217
230
.
Senay
G. B.
,
Leake
S.
,
Nagler
P. L.
,
Artan
G.
,
Dickinson
J.
,
Cordova
J. T.
&
Glenn
E. P.
2011
Estimating basin scale evapotranspiration (ET) by water balance and remote sensing methods
.
Hydrological Processes
25
(
26
),
4037
-
4049
.
Senay
G. B.
,
Friedrichs
M.
,
Singh
R. K.
&
Velpuri
N. M.
2016
Evaluating Landsat 8 evapotranspiration for water use mapping in the Colorado River Basin
.
Remote Sensing of Environment
185
,
171
185
.
Soleimani-Motlagh
M.
,
Soleimani-Sardo
M.
&
Mossivand
A. M.
2022
The efficiency of the Standardized Evapotranspiration Deficit Index (SEDI) in assessing the impact of drought on vegetation cover
.
Environmental Monitoring and Assessment
194
(
4
),
1
12
.
Sumner
D. M.
,
Nicholson
R. S.
&
Clark
K. L.
2012
Measurement and Simulation of Evapotranspiration at a Wetland Site in the New Jersey Pinelands
.
US Geological Survey Scientific Investigations Report 2012-5118, 30 p
.
Sun
G.
,
Alstad
K.
,
Chen
J.
,
Chen
S.
,
Ford
C. R.
,
Lin
G.
,
Liu
C.
,
Lu
N.
,
McNulty
S. G.
,
Miao
H.
&
Noormets
A.
2011a
A general predictive model for estimating monthly ecosystem evapotranspiration
.
Ecohydrology
4
(
2
),
245
255
.
Sun
G.
,
Caldwell
P.
,
Noormets
A.
,
McNulty
S. G.
,
Cohen
E.
,
Moore Myers
J.
,
Domec
J.C.
,
Treasure
E.
,
Mu
Q.
,
Xiao
J.
&
John
R.
2011b
Upscaling key ecosystem functions across the conterminous United States by a water-centric ecosystem model
.
Journal of Geophysical Research: Biogeosciences
116
(
G3
), 1–16.
Sun
G.
,
McNulty
S. G.
,
Amatya
D. M.
,
Skaggs
R. W.
,
Swift
L. W.
Jr.
,
Shepard
J. P.
&
Riekerk
H.
2002
A comparison of the watershed hydrology of coastal forested wetlands and the mountainous uplands in the Southern US
.
Journal of Hydrology
263
(
1–4
),
92
104
.
The Earth, Parts A/B/C 30(1–3), 69–79
.
Tobler
W. R.
1970
A computer movie simulating urban growth in the Detroit region
.
Economic Geography
46
(
supp1
),
234
240
.
US Army Corps of Engineers, New York District
1987
Flood Protection Feasibility Report for the Main Stem of the Passaic River, Volume I, Main Report and Environmental Impact Statement; Volume IV, Appendix C, Plan Formulation. NY: NY
.
USDA Natural Resources Conservation Service (NRCS)
2011
NEW YORK Rapid Watershed Assessment Profile, Hackensack-Passaic Watershed
.
United States Department of Agriculture: Syracuse State Office
.
Valipour
M.
2015
Temperature analysis of reference evapotranspiration models
.
Meteorological Applications
22
(
3
),
385
394
.
Watkins
M. M.
,
Wiese
D. N.
,
Yuan
D. N.
,
Boening
C.
&
Landerer
F. W.
2015
Improved methods for observing Earth's time variable mass distribution with GRACE using spherical cap mascons
.
Journal of Geophysical Research: Solid Earth
120
(
4
),
2648
2671
.
Wiese
D. N.
,
Yuan
D. N.
,
Boening
C.
,
Landerer
F. W.
&
Watkins
M. M.
2019
JPL GRACE Mascon Ocean, Ice, and Hydrology Equivalent Water Height RL06 CRI Filtered Version 02, version 2. PO. DAAC
.
Williams
D. G.
,
Cable
W.
,
Hultine
K.
,
Hoedjes
J. C. B.
,
Yepez
E. A.
,
Simonneaux
V.
,
Er-Raki
S.
,
Boulet
G.
,
De Bruin
H. A. R.
,
Chehbouni
A.
&
Hartogensis
O. K.
2004
Evapotranspiration components determined by stable isotope, sap flow and eddy covariance techniques
.
Agricultural and forest meteorology
125
(
3-4
),
241
258
.
Zhang
K.
,
Kimball
J. S.
,
Nemani
R. R.
&
Running
S. W.
2010
A continuous satellite-derived global record of land surface evapotranspiration from 1983 to 2006
.
Water Resources Research
46
(
9
),
1
21
.
Zhang
K.
,
Kimball
J. S.
&
Running
S. W.
2016
A review of remote sensing based actual evapotranspiration estimation
.
Wiley Interdisciplinary Reviews: Water
3
(
6
),
834
853
.
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY 4.0), which permits copying, adaptation and redistribution, provided the original work is properly cited (http://creativecommons.org/licenses/by/4.0/).