Recent Updates of the Calibration-Free Evapotranspiration Mapping (CREMAP) Method

1. Introduction

This study is the revision of an earlier evapotranspiration (ET) estimation technique [1], called Calibration-Free Evapotranspiration Mapping (CREMAP), applied for Nebraska (Figure 1). Major differences between the current and previous versions are

1. the ET rates of the winter months (December, January, February) are included in the current annual maps;

2. the method became largely automated and thus its application made simpler.

In principle, the recent modeling period could have been extended to include 2010 and 2011, however, the multi-institutional research project [2] that provided the monthly incoming global radiation values was terminated in 2010, thus no radiation data are available after June 2010 from that source. Rather than looking for alternative data sources, the original 10-year long modeling period, i.e., 2000-2009, was kept, thus ensuring that the same data types were employed throughout the study.

In the earlier version of the ET maps it was assumed that the ET rates in the winter time are negligible when one is concerned with the mean annual value. In the light of the present version of ET mapping, this assumption was found true only partially: there are regions within Nebraska for which winter ET indeed seems to be negligible (mostly the north-central and north-western parts) in comparison with water-balance data, while in other regions it is not so. These latter regions include parts of Nebraska (mostly the south and south-west portion) with the highest winter-time daily maximum temperatures and/or with most abundant precipitation (eastern, south-eastern portion of the state). As a result the precipitation (P) recycling ratio (i.e., ET / P) rose from a previously estimated mean annual value of 93% to 95%, leaving an estimated 5% of the precipitation to emerge as runoff (Ro) in the streams.

Naturally, as any estimation method, the current approach is not perfect. In the Pine Ridge Escarpment and in the Niobrara River Breaks regions (Figure 1) the ET estimates (similarly to the previous version of the ET map) had to be corrected via comparison with precipitation data because otherwise they would overestimate ET rates by about 10-20%. The reason is in the gross violation of the underlying hypotheses of the current ET estimation method in areas of very rough terrain. After the corrections in these distinct geomorphic regions, it is believed that the resulting ET rates are quite reasonable across the whole state. Overall, the method yields a state-representative ET rate value of 549 mm/yr, within 2% of the simplified water balance (P – Ro) derived rate of 538 mm/yr, employing the USGS [3] values of computed runoff for catchments with level-8 hydrologic unit codes (HUC), and explains 87% of the observed spatial variance of the water balance ET values among the HUC-8 catchments (there are 70 such watersheds within the state) for the same period.

2. Description of the current ET estimation method

An ET estimation method had been proposed by Bouchet [4], employing the complementary relationship (CR) of evaporation which was subsequently formulated for practical regional-scale ET applications by Brutsaert & Stricker [5] and Morton et al. [6]. In this study the WREVAP program of [6] was applied for the estimation of the regional-scale ET rates at monthly periods. Disaggregation of the regional ET value in space is based on the Moderate Resolution Imaging Spectroradiometer (MODIS) data [7] that have a nominal spatial resolution of about 1 km. The disaggregation is achieved by a linear transformation of the 8-day composited MODIS daytime surface temperature (T_s) values into actual ET rates on a monthly basis [1, 8] by first aggregating the composited T_s data into monthly mean values. Compositing is used for eliminating cloud effects in the resulting composite data by removing suspicious, low pixel-values in the averaging over each eight-day period. See [7] for more detail of data collection and characteristics.

Figure 1.

Stream network and selected geomorphic regions of Nebraska. MC: McConaughy Reservoir; LC: Lewis and Clark Reservoir; L: Lincoln; O: Omaha.

The transformation requires the specification of two anchor points in the T_s versus ET plane (Figure 2). The first anchor point is defined by the spatially averaged daytime surface temperature, <T_s>, and the corresponding regionally representative ET rate, E, from WREVAP. (The original FORTRAN source code can be downloaded from the personal website of the author: snr.unl.edu/ szilagyi/szilagyi.htm). The second anchor point comes from the surface temperature, T_sw, of a completely wet cell and the corresponding wet-environment evaporation, E_w, (defined by the Priestley-Taylor [9] equation with a coefficient value of 1.2). The two points identify the linear transformations of the T_s pixel values into ET rates for each month. The resulting line is extended to the right, since in about half the number of the pixels ET is less than the regional mean, E. A monthly time-step is ideal because most of the watershed- or large-scale hydrologic models work at this time-resolution, plus a monthly averaging further reduces any lingering cloud effect in the 8-day composited T_s values. Wet cells within Nebraska were identified over Lake McConaughy and the Lewis and Clark Lake on the Missouri River (Figure 1). An inverse-distance weighting method was subsequently used to calculate the T_sw value to be assigned to a given MODIS cell for the linear transformation.

Figure 2.

Schematics of the linear transformation of the MODIS daytime surface temperature values into ET rates (after [1]), applied in CREMAP.

The core assumption of CREMAP is that the surface temperature of any MODIS cell is predominantly defined by the rate of evapotranspiration due to the large value of the latent heat of vaporization for water and that the energy (Q_n) available at the surface for sensible (i.e., heat convection) and latent heat (i.e., evapotranspiration) transfers are roughly even among the cells of a flat-to-rolling terrain. Heat conduction into the soil is typically negligible over a 24-hour period, and here considered negligible over the daytime hours as well. This last assumption is most likely true for fully vegetated surfaces where soil heat conduction is small throughout the day, and is less valid for bare soil and open water surfaces. While a spatially constant Q_n term at first seems to be an overly stringent requirement in practical applications due to spatial changes in vegetation cover as well as slope and aspect of the land surface, Q_n will change only negligibly in space provided the surface albedo (i.e., the ratio of in- and outgoing short-wave radiation) value also changes negligibly among the pixels over a flat or rolling terrain. For the study region, the MODIS pixel size of about 1 km may indeed ensure a largely constant Q_n value among the pixels since the observed standard deviation in the mean monthly (warm season) surface albedo value of 17% is only 1.2% among the MODIS cells.

A further assumption of the method is that the vertical gradient of the air temperature near the surface is linearly related to the surface temperature [10, 11], thus sensible heat (H) transfer across the land-atmosphere interface, provided changes in the aerodynamic resistance (r_a) among the MODIS pixels are moderate, can also be taken a linear function of T_s. This can be so because under neutral atmospheric conditions (attained for time-steps a day or longer) r_a depends linearly on the logarithm of the momentum roughness length, z_0m [11], thus any change in z_0m between pixels becomes significantly dampened in the r_a value due to the logarithm. Consequently, the latent heat (LE) transfer itself becomes a linear function of T_s under a spatially constant net energy (Q_n) term required by the CR, therefore Q_n = H + LE, from which LE = mT_s + c follows, m and c being constants for the computational time step, i.e., a month here, within a region.

8-day composited MODIS daytime surface temperature data were collected over the 2000 – 2009 period. The 8-day composited pixel values were averaged for each month to obtain one surface temperature per pixel per month, except for December, January, and February. The winter months were left out of the linear transformations because then the ground may have patchy snow cover which violates the constant Q_n assumption since the albedo of snow is markedly different from that of the land. Therefore in the wintertime the WREVAP-derived regional ET rates were employed without any further disaggregation by surface temperatures but, rather, with a subsequent correction, discussed later.

Mean annual precipitation, mean monthly maximum/minimum air and dew-point temperature values came from the PRISM database [12] at 2.5-min spatial resolution. Mean monthly incident global radiation data at half-degree resolution were downloaded from the GCIP/SRB site [2]. While previously [1] the regions were defined by subdividing the state into eight distinct areas (a largely arbitrary process) for the calculation of the regionally representative values of the mean monthly air temperature, humidity and radiation data, required by WREVAP, now such a subdivision is not necessary. Instead, a “radius of influence” is defined over which the regional values are calculated separately for each designated MODIS cell, very similar to a temporal moving-average process, but now in two dimensions of space. In principle, such a spatial averaging could be performed for each MODIS cell, in practice however, it becomes computationally overwhelming on the PC, and it is also unnecessary, since the spatial averages form a 2-D signal of small gradient, making possible that “sampling” (i.e., the actual calculation of the spatial mean values including the WREVAP-calculated ET rate) is performed only in a selected set of points, which was chosen as each tenth MODIS cell in space (both row-, and column-wise). The remaining cells were then filled up with spatial mean values, linearly interpolated first by row among the selected MODIS-cell values, and then by column, involving the already interpolated values in the rows as well. Near the eastern and southern boundaries of the state any necessary spatial extrapolation was done by the gradient method (i.e., keeping the first two terms of the Taylor-expansion). This “sampling” sped up calculations by at least two orders of magnitude.

Care had to be exercised with the choice of the radius of influence. Rather than applying a constant radius, a spatially changing one was required because near the boundary of the state the “window” becomes asymmetrical around the MODIS cells, therefore the radius changed linearly with distance to these boundaries from a starting value of 25 cells up to a maximum of 125 cells (at a rate of 4/5 cell by each line or column) in the central portion of the domain. It was simpler to define a rectangular region around each designated MODIS cell, rather than a circular one, therefore the radius of influence is half the side-length of the resulting square.

Once the spatial mean values were available for each MODIS cell, the actual linear transformation of the T_s to ET values was performed for each month (except the winter months). The linear transformation of the T_s values into ET rates assumes a negligible change in the r_a value among the cells. As was mentioned above, r_a is directly proportional –up to a constant and with a negative slope— to the logarithm of z_0m under neutral stability conditions of the atmosphere, provided the wind speed at the blending height (about 200 m above the ground) is near constant within the region [11]. The momentum roughness height, z_0m, of each MODIS cell has been estimated over the state (Figure 3) with the help of a 1-km digital elevation model, as the natural logarithm of the standard deviation in the elevation values among the 25 neighboring cells surrounding a given cell, including the chosen cell itself. The minimum value of z_0m has been set to 0.4 m, so when the estimate became smaller than this lower limit (possible for flat regions), the value was replaced by 0.4 m. Note that the z_0m = 0.4 m value is the upper interval value for a “prairie or short crops with scattered bushes and tree clumps” in Table 2.6 of [13]. The rugged hills regions of Nebraska (e.g., the Pine Ridge and the Pine Bluffs, just to name a few) are characterized (Figure 3) by the largest z_0m values (larger than 3 m), covering the 3-4 m range for “Fore-Alpine terrain (200-300 m) with scattered tree stands” of [13]. Since the r_a estimates are proportional (up to a constant) to the logarithm of the z_0m values, their change among the MODIS cells is much subdued: about 67% of the time they are within 5% of their spatial mean and more than 94% of the time they remain within 15% (Figure 4), supporting the original assumption of the present ET mapping method.

Cells that had r_a values smaller than 95% of the mean r_a value (involved about 20% of all cells) were identified, and the corresponding ET values corrected by the relative change in r_a, considering that the sum of the latent (LE) and sensible heat (H) values are assumed to be constant among the cells (equaling Q_n) and that H is proportional to dT_z / r_a [11], where dT_z is the vertical gradient of the air temperature above the surface, itself taken proportional to T_s. The reason that only the “overestimates” of ET are corrected is that the linear transformation of the T_s values into ET rates seems to be more sensitive to more rugged-than-average terrain than to smoother one. That is why the most rugged part of Nebraska, i.e., the Pine Ridge, required an additional (to the above) 10% ET adjustment if no Ponderosa Pine was detected in the 3x3-cell region of the land use-land cover map around a given cell, and a 20% cut if it was. The assumption is that in this extremely rugged region cells with other than Ponderosa pine designation, may still contain scattered trees, if in the vicinity there are pine-forested areas plus the air turbulence, enhanced by the rugged terrain, may have a wake with a characteristic length of about a km. Within the Niobrara River Breaks region (less rugged than the Pine Ridge) only a 10% additional ET adjustment was applied without regard if the cell-surroundings are pine-covered or not. The underlying reasons of these deviations may be (after accepting that the PRISM precipitation values are correct) the way z_0m is estimated, perhaps a DEM with a finer resolution would yield better results.

Or maybe the type of vegetation, even at a 1-km resolution, has relevance (similar to plot-scale applications), in addition to the primary elevation variance. Or even, due to the enlarged surface area of the rugged terrain, the global radiation value should be reduced, which would lower ET. This topic certainly requires further research.

Figure 3.

Estimated values of the a) momentum roughness height (z_0m), and; b) relative change in the aerodynamic resistance (r_a) around the state-wide mean. The numbers along the left and bottom edge of the panels are the MODIS cell numbers.

A final correction was applied for cells of “extreme” elevation. Namely, when the elevation of a cell differed from the regional mean value by more than 100 m, its surface temperature was corrected by 0.01 Kelvin per meter, reflecting the dry-adiabatic cooling rate of the air.

Figure 4.

Relative histograms of the momentum roughness height (z_0m) and the relative change in the aerodynamic resistance (r_a) around its spatial mean value [m(r_a)] across Nebraska, estimated from a 1-km resolution digital elevation model.

The WREVAP model is based on the complementary relationship [4] which performs the worst in the cold winter months [8, 14, 15], thus the resulting WREVAP-obtained winter ET rates become the most uncertain. A yet unpublished study by the present author, involving the Republican River basin, indicated that inclusion of the winter ET rates of WREVAP improved the mean annual ET estimates in comparison with water balance derived [3, 12] data. Other studies [1, 16] also indicated that WREVAP somewhat overestimates ET rates in the Nebraska Sand Hills region even without inclusion of the winter months. Finally, a water balance based [3, 12] verification of the current ET estimates indicated that the WREVAP-provided winter ET rates are necessary in the most humid eastern, south-eastern part of the state. Based on these comparisons, the WREVAP winter months were fully included in the mean annual ET rates [besides the wettest part of the state, defined by (P_sm – ET_WREVAP) > 50 mm] for the Republican River basin, and for areas where the mean monthly daytime maximum temperature values exceeded 5 ºC. The latter area almost fully covers the Republican River basin, plus the south and south-western part of the panhandle region. P_sm designates the spatially smoothed precipitation values of PRISM, applying a 30-by-30-cell window, to filter out the unrealistic grainy structure of the PRISM precipitation field (Figure 5) due probably to its spatial interpolation method. No winter ET rates were included in the mean annual ET values wherever (P_sm – ET_WREVAP) < 10 mm; and a 50% reduction of the WREVAP winter ET rates were used for areas where [10 mm < (P_sm – ET_WREVAP) < 50 mm] held true.

Figure 5.

Distribution of the mean annual precipitation (mm) values across Nebraska from the PRISM data, 2000-2009. The state-wide mean annual precipitation rate is 577 mm.

3. Results and conclusion

The mean annual ET rates across Nebraska are displayed in Figure 6. By and large, ET follows the distribution of precipitation, as expected. Most of the ET values are between 250-500 mm in the panhandle, around

Figure 6.

Estimated mean annual ET rates (mm) in Nebraska (2000-2009). The state-wide mean ET value is 549 mm/yr.

500-650 mm in the middle of the state and near 650 mm in the eastern portion of it. Locally, however, there are large differences due to land use and land cover variance. The sizeable reservoirs (McConaughy, Lewis and Clarke, Harlan County, Swanson, Calamus, etc.) large enough to fully accommodate a MODIS cell, display the largest ET rates, around 1000 mm annually. The reservoirs/lakes are followed by the wider rivers (i.e., Platte, Missouri, Loups, Elkhorn, and partly the Republican) and their valleys in ET rates. The reason, beside the presence of the open water surface, is in the relatively small distance to the groundwater table in these river valleys, enabling the root system of the vegetation to tap into it, plus in the accompanying large-scale irrigation within the valleys. The river valleys on the ET scale are followed by areas of intensive irrigation, reaching 750 mm a year. The driest regions, with the smallest rate of ET in eastern Nebraska are the urban areas of Omaha and Lincoln (Figure 1), where the built in surfaces enhance surface runoff. The eastern outline of the Sand Hills is clearly visible, as well as the sandy areas (the elongated green-colored features) between the Loup and the Platte Rivers. The sandy soil, due to its large porosity favors deep percolation of the water often out of reach of the vegetation.

Figure 7 depicts the monthly ET rates from January through December. In July and August the irrigated plots in south-central Nebraska can evaporate as intensively as the open water surfaces. Even in September, when most of the produce has been harvested, the soil through its enhanced moisture due to summer irrigation, evaporates more than the surrounding, non-irrigated land. In November the distribution of ET rates becomes zonal and follows the precipitation distribution.

Figure 7.

Estimated mean monthly ET rates (mm) in Nebraska., 2000-2009.

While in absolute numbers the south-central portion of the state produces the highest ET rates, the picture changes significantly, when one looks at the ET to precipitation (so-called precipitation recycling) ratios of Figure 8.

Figure 8.

Estimated mean annual ET to precipitation ratios (%), 2000-2009. The state-wide mean ratio is 95%.

Figure 9.

Distribution of areas with the largest observed groundwater decline (at least 3 m, to up to 8 m) over the 2000-2009 period, overlain the ET / P map. For the correct color to ratio correspondence, please, use the colors in Figure 8 instead of the current ones.

Lake McConaughy is the clear winner (followed by smaller lakes in the vicinity), evaporating about twice as much as it receives from precipitation. It does not mean, of course, that the other small lakes in the area would not evaporate as much as Lake McConaughy per unit area, they probably evaporate even more (the smaller the lake the larger typically its evaporation rate, provided other environmental factors are equal), but their size inhibits MODIS to detect their surface temperatures without “contamination” from the surrounding land. Note again the eastern outline of the Sand Hills and the elongated sandy areas between the Loup and Platte Rivers as areas of relatively low ET rates. The two urban areas of Omaha and Lincoln are clearly visible again.

Two large irrigated areas stand out clearly as the most intensive water users (relative to precipitation), one in the Republican River basin and the other in the North-Platte River valley of the panhandle. In these areas ET rates significantly exceed (up to 50%) precipitation rates. Another significant irrigated area in Box Butte County (at the western edge of the Sand Hills) plus the one in the Republican River basin coincide largely with regions of extensive groundwater depletions, displayed in Figure 9.

The lines in Figure 9 designate areas (after [17]) where groundwater decline was at least 3, 5, 8 m over the 2000-2009 period. Naturally, in heavily irrigated areas close to major streams (e.g., North- and South Platte, Platte River), such groundwater depletions are absent (but not around Lake McConaughy, where reservoir water levels have been below normal most of 2000-2009) since the chief source of the irrigation water is the stream itself. Figure 10 displays the distribution of irrigated land, overlain the ET / P values.

Figure 10.

Irrigated land (marked in black) distribution [after 18, 19] in Nebraska, 2005.

As seen in Figure 10, not all land areas with larger than unity ratios are connected to irrigation, good examples are the Sand Hills wetlands. Similarly, not all areas that come up with values larger than 100% do actually evaporate more than they receive from precipitation. Such an area is the table-land just south of the western edge of Lake McConaughy, between the North- and South-Platte Rivers (please, refer to Figure 8 for corresponding precipitation recycling ratios, the colors of Figure 10 are slightly off because it was produced by a different software that enabled marking the irrigated areas on top of the ratios). In this area irrigation is largely absent (or at least it was in 2005, the date of the irrigation data), yet the ratio is between 100-120%. The error may be caused by several factors, namely a) the well-known, often 10% underestimation of precipitation; b) inaccuracy of the ET estimates, and; c) problems with the spatial interpolation of the measured precipitation values. The latter is well demonstrated in Figure 5, which shows that in the southern panhandle region there can be a difference of 125 mm (about 25-30% of the annual value) in the precipitation values within a distance of 30 km or less. Added to this uncertainty is the wide-spread underestimation of precipitation, especially in windy areas where a measurable portion of the raindrops (and especially snowflakes) is swept away from the rain gage. Finally, the present ET estimates have an error term (discussed further later) of about plus/minus (±) 5-10%. Employing a ±5% error in the latter, another -5% underestimation in the measured precipitation values together with yet another independent ±5% error in the interpolated values, the resulting ET / P ratio may contain an error of -5% to 16%, coinciding well with the error extent found in the table-land area. Therefore, the ratios in Figure 8 must be treated with this uncertainty in mind.

A comparison with the previous version [1] of the mean annual ET map (Figure 11) reveals that the largest differences are found in the Republican River basin and the southern panhandle region, where now the full values of the WREVAP-estimated winter-time ET rates were added to the warm-season values (March-November). Note that the procedures used for preparing the two maps are different (application of a radius of influence around each MODIS cell versus distinct geographic regions) as was explained above. The perceptible diagonal and level straight lines suggest some problems with the interpolation method employed in the previous ET map.

Figure 11.

Differences in the present and previous [1] mean annual ET maps (mm). Mean is 18.5 mm.

Verification of the estimated mean annual ET rates can be best performed on a watershed-by-watershed basis by subtracting the stream discharge values (expressed in mm) from the mean precipitation values of the catchments, assuming that groundwater level changes are negligible over the study period, i.e., 2000-2009. As seen above, the latter is not true in many regions within Nebraska, but a transformation of these groundwater-level changes into water depth values would require the state-wide distribution of the specific yield value (also called drainable porosity, defined by dividing the drained water volume value with that of the control volume, fully saturated with water at the start of drainage) of the water bearing aquifer, a hydro-geological parameter not available for the whole state. Figure 12 displays the water-balance derived ET rates employing the PRISM precipitation values [12] and USGS-derived watershed-representative runoff values [3] for the HUC-8 watersheds within Nebraska, while Figure 13 displays the spatial distribution of estimation error (predicted minus water-balance derived) of the CREMAP ET values among the same catchments. As seen, in the majority of the watersheds the estimation error is within 30 mm of the “observed” value. The largest overestimation takes place for watersheds within the Missouri River basin, within the Sand Hills, and within and just north of the Republican River basin. The latter area corresponds to one of the most severe groundwater depletion regions within the state, where the “missing” water certainly contributes to elevated ET rates, not detectable by the simplified water balance approach.

Figure 12.

Water balance derived (P – Ro) mean annual ET rates (mm) of the USGS HUC-8 watersheds in Nebraska, 2000-2009.

Note also that a systematic underestimation of the precipitation rates automatically leads to a virtual overestimation of the ET rates by the present method during verification. Another problem with the verification is that the watershed area employed for the transformation of the discharge values into water depth, may also be somewhat uncertain, since the groundwater catchment does not always overlap perfectly with the surface-water catchment delineated by the help of surface elevation values. Probably that is why the largest over- and under-estimation of ET is found within the Sand Hills, in catchments very close to each other. Also, the USGS watershed runoff values employ simplifications that may cause serious errors in the estimated watershed runoff rate, such as found for the Lower Republican basin in Kansas (not shown here), where USGS reports a mean runoff rate of 5 mm/yr for 2000-2009, while a Kansas Geological Survey study [20] found a mean annual runoff rate of 106 mm/yr, an almost twenty-fold difference.

Figure 13.

Distribution of the estimation error (mm) in the mean annual ET rates among the USGS HUC-8 watersheds within Nebraska, 2000-2009.

Figure 14 summarizes the ET verification results. It can be seen that in the vast majority of the USGS HUC-8 watersheds the estimated values are within 10% of the simplified water balance derived values. Five of the seven overestimates (above the upper intermittent line) of Figure 14, found between 400 and 500 mm, correspond to the large groundwater depletion area in and around the Republican River basin, displayed in Figure 9, so in those cases the CREMAP ET estimates may better represent reality than the simplified water balance derived values. The explained variance, R², has a value of 0.87, meaning that 87% of the spatial variance found in the HUC-8 water-balance derived ET values is explained by the CREMAP estimates. In summary, the annual and monthly ET maps are recommended for use in future regional-scale water-balance calculations with the resolution and accuracy of the estimates kept in mind. The maps are certainly not recommended for reading off individual cell values, because the exact cell coordinates maybe slightly off due to the geographically referenced data manipulations necessary to produce those maps. For example, the author found some problem with coordinate referencing when cells are extracted from a grid employing another grid with differing cell size. The maps are best suited for studies of spatial scale larger than one km.

Figure 14.

14. Regression plot of the water-balance derived and CREMAP-estimated mean annual ET rates (mm) among the USGS HUC-8 catchments. R² is the portion of the spatial variance of the water balance ET rates that is explained by the CREMAP estimates. The upper and lower envelope lines designate the P – Ro value plus or minus 10%.

Acknowledgments

This work has been supported by the Hungarian Scientific Research Fund (OTKA, #83376) and the Agricultural Research Division of the University of Nebraska. This work is connected to the scientific program of the "Development of quality-oriented and harmonized R+D+I strategy and functional model at BME" project. This project is supported by the New Szechenyi Plan (Project ID: TAMOP-4.2.1/B-09/1/KMR-2010-0002).

Disclaimer: The views, conclusions, and opinions expressed in this study are solely those of the writer and not the University of Nebraska, state of Nebraska, or any political subdivision thereof.