Potential evapotranspiration method influence on climate change impacts on river flow: a mid-latitude case study

Projected changes in 21st century climate are likely to impact global water resources substantially, in part due to an anticipated intensification of the hydrological cycle (for example, following the Clausius-Clapeyron relationship). A general pattern of wet areas becoming wetter and dry areas drier is emerging as a robust thermodynamic response to global warming (e.g. Held & Soden 2006). However, the extent to which this relationship holds at the local scale has been questioned (Roderick et al. 2014) and is further complicated by simultaneous changes in plant stomatal conductance under higher atmospheric carbon dioxide concentrations (e.g. Kay et al. 2013). Given these uncertainties, representation of evapotranspiration (ET) from the land surface is a critical aspect of hydrological modelling under climate change scenarios.

In hydrological models, actual ET (AET) is rarely measured explicitly, with ET typically quantified using a conceptual variable, potential ET (PET). PET can be derived using a variety of meteorological variables and empirical relationships (e.g. Lu et al. 2005). One of the most complete methods is the Penman (1948) equation (and subsequent versions, e.g. Allen et al. 1998), which directly incorporates all of the controlling meteorological variables: net radiation, temperature, saturation vapour pressure deficit (VPD) and wind speed. In many instances, however, there are insufficient data for all four variables, leading to the use of a number of alternative equations that are based on fewer meteorological variables and which make greater use of empirical relationships (e.g. Priestley-Taylor, Hargreaves-Samani).

Whilst a number of empirical relationships have been found to represent PET satisfactorily in different locations (e.g. Vörösmarty et al. 1998; Lu et al. 2005), it has been demonstrated that net radiation, temperature, VPD and wind speed all need to be considered to enable a full understanding of historical trends and variations in PET (e.g. Roderick et al. 2007; Donohoe et al. 2010). Consideration of all four climatological variables is expected to be equally important under climate change, given the multitude of possible combinations of changes in these variables. Indeed, this is evidenced by the work of a number of past studies, at the station scale (e.g. McKenney & Rosenberg 1993), national scale (e.g. Prudhomme & Williamson 2013), and global scale (e.g. Kingston et al. 2009), which showed that different PET methods produce markedly different estimates of the PET climate change signal (even to the extent of different directions of change).

Previous studies have shown that differences in the PET climate change signal originating from choice of PET method can go on to influence the climate change signal for catchment-scale water resources. For example, Bae et al. (2011) found variation between seven PET methods could cause up to 14% variation in annual runoff in the Chungju basin in Korea for the 2071–2100 period, but over 100% variation in the monthly climate change signal. Similarly, Kay & Davies (2008) found that for three catchments across Britain, the annual and seasonal uncertainty in median, high and low flow climate change signal could vary by up to 20% between just two PET methods. For the Mekong river in southeast Asia, a comparison of six PET methods resulted in a maximum variation in runoff climate change signal of about 10% (Thompson et al. 2014). Notwithstanding these important contributions, uncertainty associated with PET method remains a relatively under-researched subject – it is difficult to draw parallels between studies based on tropical continental-scale rivers and fully distributed hydrological models (e.g. Thompson et al. 2014), and those based on relatively simple models in either highly seasonal (Bae et al. 2011) or temperate catchments (Kay & Davies 2008). The possible magnitude of PET method-related uncertainty is further complicated by consideration of the work of Oudin et al. (2005), which indicated that choice of PET method actually had very little impact on hydrological model performance across 308 catchments in Australia, France and the USA. Although not focussed on climate change impacts on runoff, these results suggest that the sensitivity to PET method indicated for the Mekong and catchments in the UK and Korea may not be present across all physiographic settings or for all model representations of catchment hydrology.

The current study adds to this important area of research through a case study of a small catchment (the Waikaia) in the South Island of New Zealand that generally has energy-limited evaporation. A conceptual hydrological model of the Waikaia is developed using HBV-Light, and is forced with a series of climate change scenarios based on five general circulation models (GCMs) and six different methods of calculating PET. As such, an estimate is provided of the possible magnitude of PET-related uncertainty in simulated runoff, within the context of the more commonly examined GCM-related uncertainty in climate change impacts.

River discharge data for model calibration were obtained for the Piano Flat gauging station in the upper Waikaia (operated by Environment Southland), giving a modelled catchment area of 489 km². The Piano Flat river flow record extends from 1979 to the present day. Meteorological data (precipitation, temperature, shortwave radiation, relative humidity and wind speed) were sourced from a combination of observed records at Piano Flat (recorded by Environment Southland) and the Virtual Climate Station network (VCSN) (Tait et al. 2006). Daily observations of precipitation at Piano Flat began in 1977, but the temperature record only goes back to 1998, which would leave a relatively short period for hydrological model calibration and validation. Furthermore, relative humidity data only exist from 2005, with no radiation or wind data recorded within the catchment. To overcome these data limitations, VCSN data corresponding to the location of Piano Flat were used for temperature, humidity, wind and radiation. The VCSN is a gridded product at approximately 5 km² resolution, based on a thin plate spline interpolation of station data. Past studies have shown that this data set is generally reliable at elevations below 500 m (Tait et al. 2012).

Climate change scenario data for precipitation, temperature, humidity and radiation were obtained from the QUEST-GSI scenario set (Todd et al. 2011; http://www.met.reading.ac.uk/research/quest-gsi/). Scenario wind data were not available, so baseline values were used throughout. The QUEST-GSI data set was developed using the ClimGen pattern scaling technique, which is based on the assumption that the spatial pattern of change in a given climate variable is constant for a given GCM, and that climate variables have a linear response to changing global mean annual temperature (Arnell & Osborn 2006; Todd et al. 2011). In this study, climate scenarios were used corresponding to global mean temperature 2 °C above the 1961–1990 baseline, the presumed threshold for dangerous climate change. This approach has been used successfully in a number of previous studies focussed on uncertainty in climate change impacts (e.g. Gosling et al. 2011; Todd et al. 2011; Thompson et al. 2013).

For the purposes of this study five GCMs were selected from the QUEST-GSI subset: CCCMA-CGCM31, IPSL CM4, MPI-ECHAM5, NCAR-CCSM30, and UKMO-HadCM3. These GCMs have been previously identified both as providing a reasonable spread of possible future climates (Todd et al. 2011), and that satisfactorily simulate the major present day features of climate in the New Zealand region (MfE 2008).

The QUEST-GSI scenario data are at a spatial resolution of 0.5° latitude/longitude, and are calculated as a change in climate from the CRU TS 3.0 climate baseline (1961–1990; Todd et al. 2011). Given that the Waikaia catchment is substantially smaller than this spatial resolution, and is situated within a larger area of land surface complexity (i.e. it is situated on a relatively small island – the South Island – which contains marked topographic variation), the change in climate from CRU TS 3.0 baseline to scenario was applied to the observed baseline data using a delta factor approach to generate scenario data for input to HBV-Light.

This investigation focuses on the three common PET method types (i.e. physically, radiation- and temperature-based). The analysis allows for the comparison of the uncertainty between the individual PET methods and different types of PET groups. Six PET methods were selected as a representative sample of the many different methods in existence (reported at over 50: Lu et al. 2005). The specific methods were selected because of their past use, particularly in other climate change impact assessment studies. The two physical methods, from Penman (1948) and Granger (1991), incorporate the four key variables associated with PET, i.e. net radiation, humidity, windspeed and temperature. The two radiative methods (Jensen & Haise 1963; Priestley & Taylor 1972) omit direct humidity and windspeed data from their calculations, and the final two methods (Hamon 1963; Hargreaves & Samani 1985) are based on empirical relationships that enable temperature alone to model the evaporative power of the atmosphere.

Following the approach of Kay & Davies (2008) and Bae et al. (2011), one PET method was selected for the baseline simulation of Waikaia river flow (Priestley-Taylor, the method that resulted in the best Nash-Sutcliffe value during calibration). The baseline-to-scenario percentage change in each of the PET methods was calculated, and then applied to the baseline Priestley-Taylor PET values to provide the scenario PET input data for HBV-Light.

The HBV-Light model is a widely used conceptual hydrologic model (Seibert & Vis 2012). The model is operated here on a daily timestep, and is driven by daily temperature, and precipitation. Furthermore, we use the mode of HBV-Light involving direct input and use of daily PET data. Temperature and precipitation data are lapsed with altitude (part of the HBV-Light calibration process); there is no corresponding function for PET so a single value has to be applied across the catchment.

A series of routines, each representing a key hydrologic function, determine the eventual runoff within HBV-Light. A brief summary follows here, for more details see Seibert & Vis (2012). Snow accumulation and storage is modelled using a freezing threshold and degree-day method. Actual evaporation and release of water to runoff are calculated according to the current level of water storage. Following rainfall or snowmelt, water is either stored in the ‘soil box’ and eventually evaporated, or enters the recharge routine. Evaporation from the soil box occurs at the potential rate until the water in this box lowers to a (user-defined) threshold, after which a linear reduction is applied. Water entering the recharge routine is released for runoff either relatively quickly, or ‘percolates’ into a slower responding water store. Finally, a user-defined smoothing function is applied to output from the quick and slow-responding water stores before the final model output is produced.

Both manual and automatic calibration of the HBV-Light model parameters were undertaken. Manual calibration was performed to identify a physically reasonable range within which parameter values were likely to vary, before Monte-Carlo based automatic calibration (part of the HBV-Light software) was used to sample across this range. The simulation began with 50 randomly generated parameter sets that are specified to certain ranges. The sets were then evaluated by running the model, and the goodness of fit calculated for each parameter set. Parameter sets with high values were given higher probability to generate new sets than those sets that gave poorer results (Seibert 2000).

The spilt sample calibration/validation procedure was adopted, following Klemes (1986). Accordingly, model parameter values were varied to achieve the best possible fit between modelled and observed discharge for the period 1980–1998. In the second half of the analysis period (1999–2013), the parameter values used in the first period were run unaltered, and model performance compared to the first period.

Model performance is at a similar level to other New Zealand-based studies. For example, the TopNet model performance for the Lindis and Matukituki catchments in the South Island resulted in monthly Nash-Sutcliffe values of 0.69 and 0.68, respectively (Gawith et al. 2012). Soil and Water Assessment Tool (SWAT) model performance within the Motueka catchment at the top of the South Island ranged from 0.31 to 0.67 between subcatchments, with a value of 0.78 for the catchment as a whole during calibration (Cao et al. 2006). Adopting the performance thresholds of Henriksen et al. (2003) for the monthly resolution (following Thompson et al. 2013), performance for the Waikaia can be classified as ‘very good’.

The annual baseline-to-scenario change in the key meteorological variables used in the various PET methods is shown in Table 1. NCAR precipitation increases the least (1.4%), while IPSL has the largest increase in precipitation of 12.3% (Table 1). CCCMA had the largest increase in mean temperature, with a value of 1.7 °C, and MPI having the smallest temperature increase (1.2 °C). Vapour pressure increases by between 0.8% (MPI) and 1.3 hPa (HadCM3, CCCMA), whilst net radiation decreases for all GCMs, by 0.3% (HadCM3) to 3.5% (CCCMA).

The annual PET percentage change for the PET methods under each GCM is shown in Table 2. With the exception of the HadCM3 scenario, the two physical PET methods (Penman, Granger) have the smallest increase or a decrease in annual PET. Granger decreases for three of the five GCM scenarios. The radiation and temperature-based methods have increases in annual PET for all GCMs. The Jensen-Haise and the Priestley-Taylor methods give relatively similar results for each GCM, whereas the Hamon and Hargreaves are frequently quite different. This is largely due to the variation in the Hargreaves method between GCMs, which ranges from a relatively small increase (MPI, 4.1% or 12.1 mm from the 294 mm baseline) to the largest increase across all method and scenarios (NCAR, 17.1% or 50.3 mm).

Most PET method-GCM combinations result in increasing annual runoff, with the notable exception of the NCAR GCM, for which all PET methods lead to a decrease in runoff (Table 3). The 0.9–6.1% decreases for NCAR represent an absolute decrease of 7.5–50.4 mm from the 832.8 baseline mean annual total. IPSL results in the largest increases in annual runoff (9.4–12.6%, or 78.2–105.2 mm). The range in annual runoff change between PET methods under each GCM is between 1.7% (HadCM3) and 3.6% (CCCMA) (Table 3). The range between GCMs for each PET method is far larger: 13.9% (Hamon) to 15.5% (Hargreaves and Penman). The two physical PET methods are generally associated with relatively high increases in annual runoff for each GCM (or relatively low decreases in the case of NCAR), but there are no obvious similarities between the two radiation- or temperature-based methods. The lowest increases in runoff generally result from the Hamon method.

In contrast to the differences in magnitude of change between GCM scenarios, the uncertainty associated with PET method is relatively small (Figure 5). Differences in climate change runoff signal between PET methods are largest in summer in both percentage and absolute terms, at up to 8% (NCAR), equating to a range of PET uncertainty in summer runoff of 6 mm for this GCM scenario. Variation in climate change signal due to PET method is generally below 5% (3 mm) for other GCMs and other times of year for NCAR. During summer and autumn, when the climate change signal in runoff is close to zero, uncertainty associated with PET method can influence the direction of change in runoff. This situation occurs in at least one month of the year for all GCMs except for HadCM3.

Differences in annual PET climate change signal of up to 16.6 percentage points are observed between six PET methods under a climate scenario of a 2°C increase in global mean temperature (Table 2). Although the scenarios are based on a relatively simplistic delta-change scenarios, we consider these results robust in terms of differences at the 30-year mean resolution. This data set has also been previously widely used and tested (e.g. Todd et al. (2011)). These differences are approximately similar to the latitudinally averaged results of Kingston et al. (2009) for the same GCM scenarios. The greater uncertainty in PET climate change signal between methods (in comparison to between GCMs) also mirrors the findings of Kingston et al. (2009).

The tendency for the two physically based methods to lead to lower increases (or decreases) in annual PET (Table 2) appears linked to changes in net radiation (Table 1), particularly for Granger. The largest decrease for this method corresponds to the GCM with the largest decrease in net radiation (CCCMA). Similarly, the largest increase in Granger, and second largest for Penman, occurs with the GCM with the smallest decrease in net radiation (HadCM3). The finding that physically based methods provide a relatively moderate response to climate change again follows previous work (e.g. Kingston et al. 2009; Shaw & Riha 2011; Prudhomme & Williamson 2013). However, it should also be noted that the absence of scenario wind data may have led to a somewhat muted response of the two physical methods to climate change. Wind speed data can be rather important for PET variability (McVicar et al. 2012), and projections indicate increases in the strength of the prevailing westerly circulation in the New Zealand region under climate change (MfE 2008).

Interestingly, the two radiation-based methods are not as closely linked to changes in radiation as the physical methods are, suggesting that changes in temperature also strongly influence Priestley-Taylor and Jensen-Haise. For example, although the HadCM3 GCM leads to the greatest increases in both methods, CCCMA does not lead to the smallest increase. Furthermore, PET for both of these methods increases for all GCMs, despite the decrease in net radiation.

The tendency for Hamon to result in the largest increase in PET (and lowest increase in runoff) is consistent with previous findings (e.g. Kingston et al. 2009; Shaw & Riha 2011). The difference in change in Hamon PET between GCMs follows the difference in Tmean change between GCMs (i.e. the greatest increase in both occurs for CCCMA, and smallest for MPI; Tables 1 and 2). In contrast, the other temperature based method (Hargreaves) does not follow GCM changes in T_mean to the same extent, indicating that the empirical relationships contained in this method do not result in a simple function of temperature.

Similar to Kay & Davies (2008), this study found that for the Waikaia catchment there is only a consistent positive percentage change in monthly PET from the temperature and radiation-based methods (Figure 4). In the winter months under the two physical methods there is mostly a decrease in PET. As with the changes at the annual resolution, this decrease appears linked to the decrease in net radiation (Table 1). Although the magnitude of the percent changes in PET appear relatively large during winter (Figure 4), this is a function of the low magnitude of baseline PET during winter.

The impact of choice of PET method on scenario runoff is relatively small, at under 5% on the monthly scale and at most 5.2% at the annual scale (Figure 5, Table 3). Similar runoff ranges were found by Thompson et al. (2014), at 5.6% across eight gauging stations within the Mekong, and the seven GCMs. Bae et al. (2011) and Kay & Davies (2008) both reported greater variation in runoff response between different PET methods (>20% for some GCM-PET method combinations), but under a scenario of greater overall climate change (SRES A2, 2071–2100) than the 2°C scenario used in the present study.

Given the different emission scenarios and GCMs used across past studies, a more informative comparison (than comparing PET percent changes) is the range in PET climate change signal between methods versus the range in runoff climate change signal between PET methods. All three of Thompson et al. (2014), Bae et al. (2011) and Kay & Davies (2008) record approximately equal ranges of climate change signal from PET method to PET influence on runoff. In contrast, for the Waikaia the mean PET related range in runoff climate change signal (3.2%, Table 3) is low compared to the mean range in signal for just PET (9.9%, Table 2). Even for the GCM with the greatest range between PET methods (NCAR), the range in annual runoff change is just 5.2%, which compares with a 16.6% variation in scenario PET for this GCM.

Aside from the general decrease in variation when moving from PET to AET, soil moisture and modelled runoff (Sperna-Weiland et al. 2012), the relatively low sensitivity of runoff climate change signal to PET method could indicate that the hydrological model as used here is somewhat insensitive to PET. Subsequent analysis of water in the HBV-Light soil box is consistent with this insensitivity, as mean monthly soil box moisture varies by 4% at most in all months for all GCMs except NCAR (where variations of up to 6% are seen in the summer months). It is possible that this apparent insensitivity could be a function of model calibration, although the model fit to observed river flow during the calibration and validation periods suggests this is not the case (Figure 2). In support of a more general insensitivity of HBV-Light to PET, this model was one of those used by Oudin et al. (2005), a study that suggested that the PET method was relatively unimportant for hydrological model performance. Bae et al. (2011) also found that the contribution of the PET method to uncertainty in scenario runoff varied between three different hydrological models.

This study examined the impact of a 2°C climate change scenario on river runoff for the upper Waikaia catchment. The primary changes comprise an increase in winter runoff, followed by a decrease in spring which is likely due to reduced proportion of precipitation falling as snow during winter, and resultant reduced snowmelt influence in spring simulated runoff. The choice of the PET method was found to lead to a relatively large difference in scenario PET. In many instances, there was a greater difference in scenario PET between methods for one GCM than between GCMs for one PET method. Such results run contrary to the usual finding that GCMs are the major source of uncertainty in climate change scenarios, and may have important implications for agricultural use of water within the catchment, particularly given recent expansion of irrigation for dairy farming. However, the changes in PET produced by the different PET methods appear to have only a minor impact on the eventual catchment runoff. The reasons for this disparity in sensitivity to PET method are not fully understood, although it seems likely that a lack of available surface water is not the cause: at the monthly resolution, ET is not water limited during the summer half-year. Limited sensitivity of the HBV-Light hydrological model to PET may be a cause of these results, but this possibility requires further research to confirm.

This work was partly supported by a University of Otago Masters postgraduate publishing bursary. Tracy Connolly (Department of Geography, University of Otago) produced Figure 1. The comments of two anonymous reviewers led to a number of improvements in this manuscript.