Agricultural Water Management 210 (2018) 340–353

Contents lists available at ScienceDirect

Agricultural Water Management

journal homepage: www.elsevier.com/locate/agwat

Accuracy of daily estimation of grass reference evapotranspiration using T

ERA-Interim reanalysis products with assessment of alternative bias
correction schemes
⁎
Paula Paredesa,b, , Diogo S. Martinsa,c, Luis Santos Pereiraa, Jorge Cadimad, Carlos Piresc
a
LEAF Linking Landscape, Environment, Agriculture and Food, Instituto Superior de Agronomia, Universidade de Lisboa, Portugal
b
Instituto Politécnico de Portalegre, Escola Superior Agrária de Elvas, Elvas, Portugal
c
Instituto Dom Luiz, Faculdade de Ciências, DEGGE, Universidade de Lisboa, Lisboa, Portugal
d
Centro de Estatística e Aplicações da Universidade de Lisboa (CEAUL) and Instituto Superior de Agronomia, Universidade de Lisboa, Portugal

A R T I C LE I N FO A B S T R A C T

Keywords: This study aims at assessing the accuracy of estimating daily grass reference evapotranspiration (PM-ETo)
Meteorological variables computed with ERA-Interim reanalysis products, as well as to assess the quality of reanalysis products as pre-
FAO Penman-Monteith ETo dictors of daily maximum and minimum temperature, net radiation, dew point temperature and wind speed,
Accuracy indicators which are used to compute PM-ETo. With this propose, ETo computed from local observations of weather
Raw reanalysis data
variables in 24 weather stations distributed across Continental Portugal were compared with reanalysis-based
Cross-validation
Bias correction methods
values of ETo (ETo REAN). Three different versions of these reanalysis-based ETo were computed: (i) an (un-
corrected) ETo based on the individual weather variables for the nearest grid point to the weather station; (ii) the
previously calculated ETo corrected for bias with a simple bias-correction rule based only on the nearest grid
point; and (iii) the ETo corrected for bias with a more complex rule involving all grid points in a 100 km radius of
the weather station. Both bias correction approaches were tested aggregating data on a monthly, quarterly and a
single overall basis. Cross-validation was used to allow evaluating the uncertainties that are modelled in-
dependently of any forcing; with this purpose, data sets were divided into two groups. Results show that ETo REAN
without bias correction is strongly correlated with PM-ETo (R2 > 0.80) but tends to over-estimate PM-ETo, with
the slope of the regression forced to the origin b0 ≥ 1.05, a mean RMSE of 0.79 mm day−1, and with EF gen-
erally above 0.70. Cross-validation results showed that using both bias correction methods improved the ac-
curacy of estimations, in particular when a monthly aggregation was used. In addition, results showed that using
the multiple regression correction method outperforms the additive bias correction leading to lower RMSE, with
mean RMSE of 0.57 and 0.64 mm day−1 respectively. The selection of the bias correction approach to be adopted
should balance the ease of use, the quality of results and the ability to capture the intra-annual seasonality of
ETo. Thus, for irrigation scheduling operational purposes, we propose the use of the additive bias correction with
a quarterly aggregation.

1. Introduction competition for water increases and water resources are gradually de-
pleted, the need to cope with water scarcity makes knowledge of in-
Evapotranspiration (ET) is a key variable in the hydrological cycle creasingly accurate values of ET more relevant (Allen et al., 1998;
and to quantify the water balance from the field to the hydrological Pereira et al., 2009; Pereira, 2017).
basin scales. Knowledge of ET is essential in water resource planning Agricultural water management practices require that crop water
and management, namely for multipurpose projects in irrigation, hy- and irrigation requirements be accurately estimated, which in turn
dropower, water transportation, flood control, and municipal and in- demands accurate knowledge of crop evapotranspiration (ETc). There
dustrial water uses, and in understanding the hydrological behavior of are various approaches for measuring and estimating ET at the various
natural and man-made landscapes (Jensen and Allen, 2016). As scales, depending upon the available data and the goals of each study

⁎
Corresponding author at: LEAF Linking Landscape, Environment, Agriculture and Food, Instituto Superior de Agronomia, Universidade de Lisboa, Portugal.
E-mail addresses: pparedes@isa.ulisboa.pt (P. Paredes), dmdmartins@fc.ul.pt (D.S. Martins), lspereira@isa.ulisboa.pt (L.S. Pereira),
jcadima@isa.ulisboa.pt (J. Cadima), clpires@fc.ul.pt (C. Pires).

https://doi.org/10.1016/j.agwat.2018.08.003
Received 1 May 2018; Received in revised form 1 August 2018; Accepted 3 August 2018
Available online 30 August 2018
0378-3774/ © 2018 Elsevier B.V. All rights reserved.
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

(Farahani et al., 2007; Allen et al., 2011; Jensen and Allen, 2016). A Sentelhas et al., 2010) and is gradually being improved (Todorovic
commonly used ETc estimation method at the field scale consists of et al., 2013; Raziei and Pereira, 2013a; Ren et al., 2016a; Almorox
using the Kc-ETo approach combining the grass reference evapo- et al., 2018; Paredes et al., 2018b).
transpiration (ETo) with a crop coefficient (Kc) as proposed in FAO56 Recent approaches to estimating ETo combine the use of remotely
(Allen et al., 1998). ETo represents the evaporative demand of the at- sensed radiation with reanalysis data. De Bruin et al. (2010,2016),
mosphere, and is thus driven by the climate. The coefficient Kc is the proposed a methodology to estimate ETo from daily values of down-
ratio between crop ET and ETo (Kc = ETc/ETo), thus represents an in- ward solar radiation obtained with radiometers onboard the geosta-
tegration of the effects of three primary characteristics that distinguish tionary satellite Meteosat Second Generation (MSG) and air tempera-
the crop from the grass reference: crop height, the crop–soil surface ture at 2 m provided by operational weather forecasts from the
resistance and the albedo of the crop–soil surface. The theoretical bases European Centre for Medium-range Weather Forecasts (ECMWF). A
of the Kc-ETo approach are well founded despite the required empiri- similar approach was used by Cammalleri and Ciraolo (2013) but using
cism in obtaining Kc (Pereira et al., 1999), and applications worldwide radiation from MSG and observed air temperature. The accuracy of
are successful as recently overviewed (Pereira et al., 2015a,b; Jensen MSG estimates of solar radiation has been evaluated (Bojanowski et al.,
and Allen, 2016). 2014).
ETo is defined as the rate of evapotranspiration from a hypothetical Interpolation of observed weather data from nearby stations was
reference crop with an assumed crop height h = 0.12 m, a fixed daily assessed by Tomas-Burguera et al. (2017) for the Iberian Peninsula.
canopy resistance rs = 70 s m−1, and an albedo of 0.23, closely re- They concluded that the best solution was to first interpolate weather
sembling the evapotranspiration from an extensive surface of green data and then compute monthly ETo. This conclusion was also reported
grass of uniform height, actively growing, completely shading the by Raziei and Pereira (2013a,b) when comparing monthly PM-ETo
ground and not short of water (Allen et al., 1998). This definition computed from observations and gridded data for Iran, and by McVicar
corresponds to the application of the Penman-Monteith equation et al. (2007) in an application to the Loess Plateau of China. In addition,
(Monteith, 1965) to a grass crop cultivated in standard, optimal con- Raziei and Pereira (2013a,b) assessed the grid network by comparing
ditions, with ETo (mm day-1) described by the daily PM-ETo Eq. (Allen the gridded weather variables with those observed by the network of
et al., 1998): weather stations. The evaluation of gridded GDAS ETo against observed
CIMIS ETo for California was reported by Senay et al. (2008) and an
900
0.408Δ (Rn-G) + γ T u2 (es−ea ) assessment of the gridded NLDAS, comparing with ground-based
mean + 273
ETo = weather data measurements, is reported by Lewis et al. (2014). Good
Δ + γ(1 + 0.34 u2 ) (1)
examples of interpolation and mapping of air moisture variables using
where Rn is net radiation at the crop surface (MJ m−2 day-1), G is soil the gridded database PRISM are provided by Daly et al. (2015) and,
heat flux density (MJ m−2 day−1), Tmean is mean daily air temperature concerning ETo, by Strong et al. (2017).
at 2 m height (°C) computed from maximum and minimum air tem- The main advantages of the reanalysis products are their spatial and
perature (Tmax and Tmin, °C), u2 is wind speed at 2 m height (m s−1), es- temporal consistency over three or more decades, and free access
ea is vapour pressure deficit (kPa) computed from the saturation vapour (Sheffield et al., 2006). In contrast, limitations of reanalysis products
pressure (es, kPa), and the actual vapour pressure (ea, kPa), Δ is the include their coarse spatial resolution (typically 80–250 km) and in-
slope of vapour pressure curve (kPa °C−1), and γ is the psychrometric sufficient reliability which seems to vary with the climatic variable
constant (kPa °C−1). analysed, location, period of study and season (Boulard et al., 2016).
The PM-ETo Eq. (1) is not an empirical equation but a mechanistic Despite research showing the improvement of the new generation of
combination equation resulting from the parameterization of the reanalysis products there is the need to evaluate these products directly
Penman-Monteith equation to a grass crop as per its definition (Allen using near surface observations over different locations and different
et al., 1994; Pereira et al., 1999). ETo is often called potential ET (PET) climates.
but this concept is different as discussed, among others, by Pereira et al. Global atmospheric reanalysis datasets provide the weather vari-
(1999) and Jensen and Allen (2016). Moreover, PET is computed by a ables required for ETo estimation with high temporal resolution. There
variety of equations and when the PM-ETo equation is used the concepts are several available sources of global reanalysis data available for PM-
above referred are rarely assumed. For this reason, there are numerous ETo estimation, namely, the European Centre for Medium-range
empirical studies comparing dozens of equations or using heuristic Weather Forecasts (ECMWF), providing for the ERA-Interim reanalysis
mathematical approaches to compute ET. Results are very often com- products (Dee et al., 2011), the Japanese Meteorological Agency (JRA)
pared with ETo computed with full data sets, or are compared with (Onogi et al., 2007); the NASA-Modern Era Retrospective-analysis for
lysimeter data, which is an approach in contrast with the mechanistic Research and Applications (MERRA) (Rienecker et al., 2011); the Na-
approach used to derive the PM-ETo equation (Allen et al., 1994). tional Center for Environmental Prediction–National Center for Atmo-
The computation of the PM-ETo Eq. (1) requires data on Tmax and spheric Research (NCEP/NCAR) (Kanamitsu et al., 2002) and the Na-
Tmin, solar radiation, humidity and u2 (Allen et al., 1998). However, in tional Oceanic and Atmospheric Administration (NOAA) (Smith et al.,
many regions and locations these variables are not observed, are not 2008). As analyzed by Sheffield et al. (2006), the reanalysis products
freely available from the relevant meteorological services, or are of are constructed from numerical weather prediction and data assimila-
poor quality due to insufficient quality control. The unavailability of tion systems using a variety of atmospheric and sea surface observa-
such data is often overcome by using alternative ET equations or tions to provide for long-term, continuous fields in time and space of
heuristic approaches that use reduced weather data sets. However, atmospheric and land surface variables. Despite observed biases in
these solutions, mainly the heuristic ones, do not use the basic physics some reanalysis products (Berg et al., 2003; Sheffield et al., 2004), they
underlying the PM-ETo equation (Pereira et al., 2015a,b). Also, have already been successfully used by many authors in representing
adopting temperature-based methods may lead to biased results due to the spatio-temporal variability of surface climate variables over the
the effects of global warming as reported by Ren et al. (2016b). To globe (Sheffield et al., 2004; Mo et al., 2011; Hwang et al., 2014).
overcome the lack of complete weather data, Allen et al. (1998) pro- Few studies focused on ETo computation using reanalysis data,
posed estimators of the parameters of the PM-ETo equation using Tmax particularly when a daily time step is considered. When a fine space and
and Tmin data only, as recently revised by Paredes et al. (2018a). This time scale is needed such as in the study of a small catchment, then
PM-ETo temperature approach (PMT) has proved accurate for both appropriate downscale models are required (Ishak et al., 2010; Jaksa
monthly and daily ETo estimation for diverse climates and regions et al., 2013; Srivastava et al., 2013a). Otherwise gridded reanalysis
(Pereira et al., 2003; Popova et al., 2006; Jabloun and Sahli, 2008; weather data may be used directly to compute ETo, even though

341
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

RHmax RHmin
requiring the interpolation of weather variables to the target locations e o (Tmin ) + e o (Tmax )
100 100
before performing the computation. Martins et al. (2017) reported on ea =
2 (5)
the good accuracy of monthly ETo estimation using a blended reanalysis
product consisting of the combination of the NCEP/NCAR reanalysis When only the mean daily RH (RHmean, %) was available ea was
with observation-based datasets (Sheffield et al., 2006). Martins et al. computed as
also conducted a brief comparison of the performance of ETo compu- RHmean ⎡ e o (Tmax ) + e o (Tmin ) ⎤
tation with NCEP/NCAR and ERA-Interim reanalysis data, which con- ea =
100 ⎣ 2 ⎦ (6)
firmed the superiority of ERA-Interim as previously reported by
Srivastava et al. (2013). Boulard et al. (2016) used ERA-Interim re- To allow assessing the quality of Tdew estimates by ERA-Interim,
analysis data for a regional analysis of ETo at different times scales from Tdew for all weather stations data was computed from ea, i.e., solving
daily to interannual showing good accuracy of variables namely for ETo Eq. (4) in order to Tdew:
computation and for water balance computations. 116.91 + 237.3 ln (ea )
The ECMWF ERA-Interim reanalysis products (Dee et al., 2011) Tdew =
16.78− ln (ea ) (7)
were selected in the present study because of its good previous per-
−2 −1
formance (Srivastava et al., 2013; Boulard et al., 2016; Martins et al., The daily soil heat flux (G, MJ m day ) beneath a dense grass
2017). However, due to bias reported for ERA-Interim products, the use canopy is very small and it is assumed G = 0 (Allen et al., 1998).
of these data requires assessing the impact of bias correction (Berg Therefore, the energy available at the grass canopy reduces to the net
et al., 2003; Zhao et al., 2008; Hwang et al., 2014; Srivastava et al., radiation Rn, which represents the net balance of the short- and long-
2015) on the performance of estimations of reanalysis-based ETo when wave solar radiation:
compared with ETo computed with observation data. An additional Rn = Rns− Rnl (8)
benefit from using those reanalysis products is that the analysis uses
−2 −1
ECMWF short and medium term weather forecasts for supporting irri- where Rns is the net shortwave radiation (MJ m day ) and Rnl is the
gation management, with short term referring to real-time irrigation net balance of the outgoing and incoming longwave radiation (MJ m−2
scheduling and medium term to planning irrigation early in the season. day−1). Rns is computed as the difference between the incoming and the
In fact, the main goal of this research is to enable the application of ETo reflected shortwave solar radiation, thus
computed with reanalysis data as input to a soil water balance and ir-
Rns = Rs (1−α) (9)
rigation scheduling and management model for use in locations where
weather data are not available (Paredes et al., 2014, 2017; Pereira where Rs is the incoming solar radiation (MJ m−2 day−1) and α is the
et al., 2015a,b) as well as its use with ECMWF short and medium term albedo of the surface, with α = 0.23 for the grass surface as assumed in
weather forecasts (Paredes et al., 2015). the definition of ETo (Allen et al., 1998). When not observed, Rs is
Using daily meteorological observations from 24 weather stations in computed from the sunshine duration (n, h) using the Angström ra-
Continental Portugal and ERA-Interim reanalysis products, the objec- diation equation:
tives of the current study consist of: 1) assessing the performance of the
n
ERA-Interim weather variables when compared to local observations; 2) Rs = ⎛as + bs ⎞ Ra
⎝ N⎠ (10)
evaluating the performance of computing daily ETo with reanalysis
weather data relative to the PM-ETo computed with full observed data where n is the actual sunshine duration (h), N is maximum possible
sets; 3) assessing the quality of PM-ETo computed with raw reanalysis sunshine duration (h), n/N is the relative sunshine duration (di-
weather data compared with bias corrected ETo; and 4) assessing the mensionless), Ra is extraterrestrial radiation (MJ m−2 day−1), as is the
efficiency of alternative bias correction procedures using cross-valida- fraction of extraterrestrial radiation reaching the earth on overcast
tion. days, and as + bs is the fraction of extraterrestrial radiation reaching
the earth on clear-sky days. The adopted default values as = 0.25 and
2. Materials and methods bs = 0.50 are those recommended by Allen et al. (1998) and, for the
Iberian Peninsula, by Azorin-Molina et al. (2015).
2.1. Reference evapotranspiration The net long wave radiation (Rnl, MJ m−2 day−1) results from the
balance between the outgoing and incoming longwave radiation:
The computation of the reference evapotranspiration PM-ETo Eq. 4
⎡ Tmax , K + Tmin
4 ⎤ Rs
⎥ (0.34−0.14 ea ) ⎛1.35 −0.35⎞
,K
(1) was performed using the procedures proposed by Allen et al. (1998). Rnl = σ ⎢ ⎜ ⎟

The daily mean saturation vapour pressure (es, kPa) was computed from ⎢ 2 ⎥ ⎝ R so ⎠
⎣ ⎦ (11)
the mean of the saturation vapour pressures at the daily maximum and −9
minimum air temperatures e o (Tmax ) and e o (Tmin ) as where σ is the Stefan-Boltzmann constant = 4.903 10 MJ K m−2 −4
−1
day , Tmax,K and Tmin,K are respectively the maximum and minimum
e o (Tmax ) + e o (Tmin ) absolute temperature during the 24-hour period (K = °C + 273.16), ea
es =
2 (2) is the actual vapour pressure (kPa), Rs/Rso is the relative shortwave
where e o (Tmax ) and e o (Tmin ) were computed as radiation, Rs is the measured or computed solar radiation (MJ m−2
day−1) and Rso is the computed clear-sky radiation =Ra (as + bs ) (MJ
17.27 T ⎞ m−2 day−1). The parameters used in Eq. 11 are those proposed by Allen
e o (T ) = 0.6108 exp ⎛
⎝ T + 237.3 ⎠ (3) et al. (1998).
o All wind speed data was adjusted to the standard height of 2 m using
When the dewpoint temperature (Tdew, C) was available, as was the
the logarithmic wind speed profile proposed in FAO56 (Allen et al.,
case when using reanalysis products, the actual vapour pressure (ea,
1998), thus:
kPa) was computed as:
4.87
17.27Tdew ⎞ u2 = uz
ea = e o (Tdew ) = 0.6108 exp ⎛ ⎜ ⎟
ln (67.8 z−5.42) (12)
⎝ Tdew + 237.3 ⎠ (4) −1
where u2 is wind speed at 2 m height (m s ), uz is the measured wind
Differently, when observations of maximum and minimum relative speed at z m height (m s−1), and z is height of measurement above
humidity, RHmax and RHmin (%), were available, ea was computed as ground surface (m).

342
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Table 1
Weather stations coordinates, anemometers height, recorded relative humidity (RH), recorded solar radiation (Rs, MJ m−2 day−1) or sunshine duration (n, h),
periods of observation and number of daily data records.
Weather stations Latitude* Longitude** Elevation (m) Anemometer height (m) Recorded RH data Recorded radiation Data records

Period Number

Alvalade do Sado 37° 55' 44” 08° 20' 45’’ 78.5 2 RHmax, RHmin Rs 2001–2013 3419
Aveiro 40° 38’ 7.44’’ 08° 39’ 34.6’’ 5 10 RHmean Rs 2003–2013 3392
Beja 38° 1’ 32.62’’ 07° 52’ 2.35’’ 246 10 RHmean n 1979–2013 11915
Braga 41° 34’ 32.13’’ 08° 27’ 3.99’’ 65 10 RHmean Rs 2003–2013 2339
Bragança 41° 48’ 14.86’’ 06° 44’ 34.69’’ 690 10 RHmean Rs 2003–2013 3637
Castelo Branco 39° 50’ 20.60’’ 07° 28’ 43.35’’ 449 10 RHmean Rs 2003–2013 3662
Coimbra 40° 12’ 48.49’’ 08° 27’ 18.55’’ 35 10 RHmean Rs 1998–2013 5557
Elvas 38° 53’ 26.88’’ 07° 08’ 23.14’’ 208 10 RHmean n 1979–2000 7529
Estremoz 38° 52' 20” 07° 35' 49” 404 2 RHmax, RHmin Rs 2006–2013 2511
Évora 38° 32’ 9.78’’ 07° 53’ 15.09’’ 246 10 RHmean n 1979–2013 12082
Faro 37° 0’ 59.68’’ 07° 58’ 19.03’’ 8 10 RHmean Rs 2003–2013 3752
Guarda 40° 31’ 42.81’’ 07° 16’ 43.23’’ 1020 10 RHmean Rs 2003–2013 3455
Leiria 39° 46’ 49.99’’ 08° 49’ 15.48’’ 45 10 RHmean Rs 2003–2013 1665
Miranda do Douro 41° 29’ 55.76’’ 06° 16’ 17.49’’ 693 10 RHmean Rs 2009–2013 1572
Montalegre 41° 49' 12” 07° 46' 48” 1005 10 RHmean Rs 2003–2013 3107
Odemira 37° 30' 06” 08° 45' 12” 91.5 2 RHmax, RHmin Rs 2002–2013 3862
Portalegre 39° 17’ 39.06’’ 07° 25’ 16.74’’ 597 10 RHmean Rs 2003–2013 3847
Portimão 37° 08’ 50.93’’ 08° 34’ 59.82’’ 38 10 RHmean Rs 2003–2013 2572
Porto 41° 14’ 0.61’’ 08° 40’ 52.80’’ 63 10 RHmean Rs 2003–2013 3627
Santarém 39° 14’ 10.54’’ 08° 41’ 23.07’’ 99 10 RHmean Rs 2003–2013 3562
Setúbal 38° 31’ 25.86’’ 08° 52’ 51.75’’ 32 10 RHmean Rs 2003–2013 3134
Sines 37° 57’ 14.50’’ 08° 50’18.61’’ 103 10 RHmean Rs 2003–2013 3866
Vila Real 41° 16’ 26.71’’ 07° 43’ 2.61’’ 561 10 RHmean Rs 1998–2013 3869
Viseu 40° 42’ 53.24’’ 07° 53’ 46.30’’ 636 4 RHmean Rs 1998–2013 5423

*Latitude - North; ** Longitude – West.

2.2. Data km−1 was adopted (Berg et al., 2003; Zhao et al., 2008; Grouillet et al.,
2016), which is an intermediate value between the dry and wet adia-
Daily weather data sets were collected across Continental Portugal batic lapse rate. This simplified approach has highly improved the re-
in 24 weather stations, 21 from the Instituto Português do Mar e da lationships between reanalysis and local Tmax, Tmin and Tdew.
Atmosfera (IPMA) and 3 from the Centro Operativo e de Tecnologia de The use of the ERA-Interim reanalysis products may be carried out,
Regadio (COTR). They are listed in Table 1 and their locations are with or without bias correction of the calculated ETo REAN (Fig. 2). The
presented in Fig. 1. The weather stations were selected so as to ensure simplest bias correction approach consists of using the data from the
regular distribution throughout the country with each station having a nearest grid point of the target location, while a possibly more precise
record of at least 4 years for comparing with reanalysis data. Daily approach considers multiple nearby grid points. Usable grid points are
weather data refers to maximum and minimum air temperatures (Tmax selected: (a) at a distance < 100 km of the target location, that concerns
and Tmin, oC) measured at 2 m height, relative humidity (RH, %), wind limited changes of the weather variables; and (b) when the forced to the
speed (u2, m s−1), and solar radiation (Rs, MJ m−2 day−1) or sunshine origin (FTO) regression coefficient b0 and the determination coefficient
duration (n, h). The quality of observed data was previously assessed R2 of the ordinary regression between Tmax REAN and Tmax OBS, as well as
using the techniques described in FAO56 (Allen et al., 1998). between Tmin REAN and Tmin OBS simultaneously satisfy respectively the
The ERA-Interim reanalysis products selected for the current study conditions 0.70 ≤ b0 ≤ 1.30 and R2 ≥ 0.70. For all locations the re-
cover the period from 1979 to present on a regular grid with a spatial analysis nearest grid point overcome these conditions.
resolution of 0.75° × 0.75° latitude-longitude, corresponding to an ap-
proximately uniform spacing of 79 km (Dee et al., 2011). Data referred 2.3. Reanalysis ETo bias correction
to eight 3-h forecasts for every day. All grid points located in Con-
tinental Portugal (Fig. 1) were considered for retrieving daily weather ETo may be computed with the referred reanalysis data (ETo REAN)
products on maximum and minimum air temperature (Tmax, Tmin, K), with or without bias correction of the computed ETo REAN. Bias cor-
dew point temperature (Tdew, K), wind speed at 10 m height (u10, rection seeks to reduce the differences between ETo computed with
m s−1), and solar radiation (Rs, W m−2). reanalysis data and with observation data, since reanalysis products are
To obtain the daily reanalysis variables for Eq. 1 (identified with the often biased due to errors in the host weather forecast models as dis-
subscript REAN) the following procedures were used: (a) the daily Tmax cussed by Berg et al. (2003) for a hydrologic application of ECMWF
REAN and Tmin REAN were selected, respectively, as the maximum and reanalysis products. Hwang et al. (2014) reported on bias correction of
minimum values among the eight daily available 3-h values of Tmax and reanalysis precipitation and temperature data used in hydrologic si-
Tmin series; (b) the daily Tdew REAN were obtained by computing the 24- mulations. Various bias correction methods have been used and as-
h average of the eight 3-h values of Tdew; (c) the daily Rs REAN values sessed in several studies (Maraun, 2013; Fang et al., 2015).
were obtained by computing the 24-h cumulative value of the 12-h Rs With a view to finding an accurate but simple bias correction pro-
values; and (d) the wind speed u2 REAN were computed first as the 24-h cedure, it was applied to the computed ETo REAN and not to the in-
average of the eight 3-h values of u10, and then these values were dividual reanalysis variables used to compute it. Bias correction applied
converted to 2 m height using Eq. 12. In what follows, all units of to reanalysis computed ETo is yet rarely reported in literature
variables were converted into the units used in Eq. 1, that is, K were (Srivastava et al., 2015) since very few applications are reported; al-
converted into °C and W m−2 into MJ m−2 day-1. In addition, all tem- though, bias correction of individual weather variables used in model
perature data were corrected for the elevation of the considered target computations has been reported (Baigorria et al., 2007; Maurer et al.,
location. Following Soares et al. (2012), a fixed lapse rate of 6.5 °C 2013).

343
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

i=n
ETo OBS (t ) = β0, t + ∑ βi, t ETo REAN , i (t ) + εt
i=1 (14)
where the multiple regression coefficients β0, β1, β2,… βn are estimated
by least squares, from the ETo OBS and the ETo REAN,i computed for each
of n reanalysis grid points i. Using the estimated multiple regression
coefficients, β̂j , it is possible to obtain new values for ETo REAN:
i=n
ˆ (t ) = βˆ +
ETo REAN ∑ βˆi, t ETo REAN , i (t )
0, t
i=1 (15)
As with the bias correction scheme a), Eqs. (14) and (15) were fitted
for monthly, quarterly and overall aggregated data.
A major uncertainty of bias correction refers to how well it performs
for conditions different from those used at calibration. Thus, a cross-
validation procedure was applied (e.g. Arlot and Celisse, 2010; Hofer
et al., 2012). The cross-validation procedure consisted in providing a
validation of model fit with a set of data that is independent of the
model fitting set. In the present study the models were fitted in-
dividually for each weather station location and validated on in-
dependent data sets for the same location. The cross-validation proce-
dure consisted in dividing each data set into two subsets of the same
size, separated in time (the more recent and the more distant data). For
each cross-validation iteration, one group was used for training and the
other was reserved for validation (Fig. 2). In a first step the correction
parameters were obtained from the first set of data (training) and the
same correction parameters were used with the second set of data
(validation). Subsequently, the second set of data was used for deriving
the correction parameters (training) and the first one for verification.
The assessment of the performance of each bias correction procedure
was performed on the validation/verification sets. The performance
results of the cross-validation were averaged over the two validation
sets.

2.4. Statistical accuracy indicators

Fig. 1. Spatial distribution of the ERA-Interim reanalysis grid points (at each
The accuracy assessment focused on the pairwise comparison be-
0.75°) in Continental Portugal and location of the 24 weather stations used in
tween observed values (Oi) and the corresponding reanalysis values (Pi)
the current analysis.
for each variable used in ETo computations and for the ETo values,
whose means are respectively O and P . Several statistical indicators
Two bias corrections schemes were applied (Fig. 2): were used to assess the performance of the reanalysis datasets in re-
presenting the amplitude and time variation of the diverse variables. To
a) a simple correction which adds a constant c(t) to the uncorrected verify the similarity between Oi and Pi the following set of indicators
ETo REAN unc computed for the nearest grid point, where t denotes a were used:
time period and i) The regression coefficient (b0) of a linear regression forced to the
origin (FTO) computed as
c (t ) = ETo REAN unc (t )− ETo OBS (t ) (13)
n
∑i = 1 (Oi Pi )
b0 = n
Thus, c(t) is the difference between the mean daily values of the un- ∑i = 1 Oi 2 (16)
corrected ETo REAN unc and of the ETo OBS computed with the observed
variables at the target location. These mean values were computed The FTO slope b0, being an overall constant of proportionality between
grouping the available daily data (for all the years used to estimate the the reanalysis and observed values, is often interpreted as a measure of
bias) for different periods of time (t): (1) twelve monthly averages bias, with b0 < 1 suggesting underestimation and b0 > 1 over-
(January to December); (2) four quarterly averages (JFM, AMJ, JAS and estimation.
OND); and (3) a single overall average. A similar approach was used in ii) The coefficient of determination (R2) of the of the ordinary least
various studies (Terink et al., 2010; Hofer et al., 2012; Hempel et al., squares regression (OLS) used to assess the dispersion of pairs of Oi and
2013). Pi values along the regression line, computed as
n
∑i = 1 (Oˆi −O)2
• a combined approach, also for the same time periods t, which uses R2 = n
∑i = 1 (Oi−O)2 (17)
the multiple regression relating ETo OBS computed with observed
data at the target location with the nearby ETo REAN unc computed at where Ôi is the fitted regression line (Oˆi = a + bPi ). Large R indicate 2

n nearby grid points that satisfy the conditions on b0 and R2 that that a large fraction of the variance of observations is explained by the
were mentioned above. model;
iii) The percent bias (PBIAS, %) measures the average tendency of
The multiple linear regression equation may be expressed as the reanalysis data to be larger or smaller than the corresponding ob-
served data with positive values indicating an over-estimation bias, and
negative values indicate an under-estimation bias.

344
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Fig. 2. Flow chart presenting the procedures to estimate daily PM-ETo from reanalysis data and comparing with observations using two bias correction methods.

345
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Fig. 3. Frequency (%) distributions of the performance indicators comparing observations and reanalysis predicted maximum and minimum temperature corrected
for elevation.

n n
∑i (Pi−Oi ) P −O ∑i = 1 (Oi−Pi )2
PBIAS = 100 n = 100 EF = 1− n
∑i = 1 Oi O (18) ∑i = 1 (Oi−O )2 (21)

iv) The root mean square error (RMSE), whose units are those of the EF values close to 1.0 indicate excellent accuracy of reanalysis esti-
considered variable, is an overall measurement of the differences be- mates. EF ≤ 0 indicate very poor estimates with MSE > σ.
tween observed and the corresponding reanalysis values. It was com- Further descriptions of these indicators and their use are reported in
puted as previous applications (Pereira et al., 2015a,b; Martins et al., 2017;
n
Paredes et al., 2018a,b).
∑i = 1 (Oi−Pi )2
RMSE =
n (19)
3. Results and discussion
v) The normalized RMSE (NRMSE, %), defined as the ratio between
RMSE and the mean of observations O, is an adimensional measure of 3.1. Assessing the accuracy of reanalysis weather variables
the relative error of estimate and was computed as follows:
The accuracy of reanalysis estimated Tmax REAN and Tmin REAN when
RMSE
NRMSE = 100 compared to observed Tmax OBS and Tmin OBS is summarized in Fig. 3
O (20)
where the frequencies of the computed indicators b0, R2, PBIAS, NRMSE
vi) In addition, the Nash and Sutcliffe (1970) model efficiency (EF) and EF are presented. The coefficients of determination (R2) are high,
was used. EF measures the relative magnitude of the mean square error mainly for Tmax, indicating that a large fraction of the variance of ob-
(MSE = RMSE2) relative to the observed data variance (σ) (Legates and servations is explained by the reanalysis estimated variables. Co-
McCabe Jr., 1999) and is computed as herently, EF are generally high for Tmax but low EF were obtained for

346
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Fig. 4. Frequency (%) distributions of the performance indicators comparing observations and reanalysis predicted shortwave incoming radiation and dewpoint
temperature.

Tmin. Consequently, errors are relatively small for Tmax, with more
frequent NRMSE ranging from 5 to 15%, while errors for Tmin are
larger, mainly in the interval 20 to 50%.
Results for b0 and PBIAS reveal under- and over-estimation ten-
dencies respectively for Tmax and Tmin, however with larger PBIAS in
case of Tmin. These results agree with those reported by Soares et al.
(2012) when estimating daily Tmax and Tmin using ERA-Interim re-
analysis for Portugal and by Martins et al. (2017) using monthly
blended NCEP/NCAR reanalysis for the Iberian Peninsula. The same
tendencies were reported for ERA-Interim estimates of Tmax and Tmin
across Australia (Fu et al., 2016) and Canada prairies (Betts and
Beljaars, 2017). Other studies also report a tendency for over-estima-
tion of Tmin by ERA-Interim reanalysis (Simmons et al., 2010; Mooney
et al., 2011) and by NCEP/NCAR products (Srivastava et al., 2015). A
study applied to the globe presented a generally good adherence of
ERA-40 reanalysis Tmax to observations but a tendency for over-esti-
mation of Tmin (Sillmann et al., 2014). Less good accuracy observed for
Tmin may be due to the influence of global warming in the models used
Fig. 5. Comparing wind speed from reanalysis and from observations through to produce the reanalysis temperature estimates as analyzed by
the linear regression between u2 REAN and u2 OBS.
Simmons et al. (2010) for ERA-40 and ERA-Interim products. These

347
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Table 2
Comparing raw (non bias-corrected) reanalysis-based ETo with ETo computed with observed weather data: frequency (%) distribution of the statistical accuracy
indicators relative to all locations.
Intervals (%) Intervals (%) Intervals (%) Intervals (%) Intervals (%)
of b0 of R2 of PBIAS (%) of RMSE (mm day−1) of EF

]1.30, [ 0 ]0.95, 1.00] 0 ] ,−20.0] 0 [0.00, 0.20] 0 ]0.95, 1.00] 0

]1.15, 1.30] 21 ]0.90, 0.95] 29 ]−20.0, −10.0] 0 ]0.20, 0.35] 0 ]0.90, 0.95] 8
]1.05, 1.15] 17 ]0.80, 0.90] 71 ]−10.0, −2.5] 8 ]0.35, 0.50] 0 ]0.80, 0.90] 63
]0.95, 1.05] 46 ]0.70, 0.80] 0 ]−2.5, 0.0] 4 ]0.50, 0.60] 17 ]0.70, 0.80] 13
]0.85, 0.95] 17 ]0.60, 0.70] 0 ]0.0, 2.5] 4 ]0.60, 0.75] 25 ]0.50, 0.70] 8
]0.70 ,0.85] 0 ]0.50, 0.60] 0 ]2.5, 10.0] 33 ]0.75, 1.00] 50 ]0.00, 0.50] 8
], 0.70] 0 [0.00, 0.50] 0 ]10.0, 20.0] 29 ]1.00, [ 8 ] , 0.00] 0
]20.0, [ 21

b0 is the regression coefficient of the FTO; R2 is the coefficient of determination of the OLS regression; PBIAS is the Percent Bias; RMSE is the Root Mean Square Error;
EF is the model efficiency.

Fig. 6. Performance of daily ETo REAN computed with bias correction, comparing RMSE and absolute PBIAS for all 24 weather stations when performing cross-
validation for: (a) an additive bias correction and (b) multiple regression comparing results to different aggregation periods: monthly (Δ), quarterly (•), and the entire
period (*).

Table 3
Performance of computing daily ETo REAN with bias correction: mean performance indicators and respective range of values for all 24 weather stations resulting from
the cross validation of using additive bias correction and multiple regression approaches computed with different aggregation periods.
Bias correction approach Aggregating period (t) Mean accuracy indicators (and range)

R2 Absolute PBIAS (%) RMSE EF

(mm day−1)

Additive Monthly 0.89 [0.81 to 0.94] 3.5 [0.1 to 12.1] 0.64 [0.48 to 0.88] 0.88 [0.78 to 0.94]
Quarterly 0.89 [0.80 to 0.94] 3.5 [0.2 to 12.3] 0.65 [0.48 to 0.89] 0.87 [0.74 to 0.94]
Entire period 0.89 [0.83 to 0.94] 3.5 [0.1 to 12.4] 0.67 [0.48 to 0.93] 0.86 [0.64 to 0.93]
Multiple regression Monthly 0.91 [0.87 to 0.95] 3.3 [0.3 to 13.3] 0.57 [0.44 to 0.89] 0.91 [0.84 to 0.95]
Quarterly 0.91 [0.87 to 0.96] 3.4 [0.1 to 16.1] 0.58 [0.44 to 0.95] 0.90 [0.83 to 0.96]
Entire period 0.90 [0.85 to 0.96] 3.5 [0.0 to 20.3] 0.61 [0.44 to 1.03] 0.89 [0.79 to 0.96]

R2 is the coefficient of determination of the OLS regression; PBIAS is the Percent Bias; RMSE is the Root Mean Square Error; EF is the model efficiency.

authors hypothesized that reanalysis is capturing warming over land in 96% of locations. These results suggest that Rs were well estimated
more than over sea, which affects night-time temperatures more, thus by Rs REAN.
making that Tmin REAN larger than Tmin OBS data. Comparable results for Rs REAN were reported in a study applied to
Results in Fig. 4 show that Rs REAN agrees well with Rs OBS, with b0 Europe using ERA-Interim reanalysis data, which outlined for
ranging 0.95 to 1.05 in 79% of cases, and most PBIAS between 2.5 and Portuguese stations a NRMSE ranging from 15% to 40% (Bojanowski
10%, thus with a slight tendency for under–estimation. R2 > 0.70 for et al., 2014). Urraca et al. (2017) reported a tendency for ERA-Interim
all cases and ranging from 0.80 to 0.90 in 79% of locations, thus in- reanalysis to slightly over-estimate daily Rs in Spain with an average
dicating that the variability of Rs REAN reflects the time variability of Rs NRMSE of 13.9%. These authors hypothesize that this over-estimation
OBS in the majority of locations. Despite the reported good agreement, may be due to lack of accuracy of ERA-Interim in estimating aerosols
the NRMSE ranged from 15 to 30% in most cases, but with EF > 0.80 loadings or the effects of water vapour when estimating the

348
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Fig. 7. Frequency (%) distribution of the statistical indicators measuring the performance of ETo estimation with reanalysis data without bias correction compared
with adopting an additive bias correction and a multiple regression correction with quarterly aggregation of data.

atmospheric transmissivity. Using blended NCEP/NCAR reanalysis in 46% of cases, however, EF < 0.5 were observed in 29% of the lo-
products, Martins et al. (2017) reported good results for monthly esti- cations. Srivastava et al. (2016) reported a general tendency for NCEP
mates of Rs for the Iberian Peninsula. Other studies relative to daily and ERA-Interim reanalysis products to under-estimate Tdew. Studies
solar radiation reported a general tendency for reanalysis products to comparing ea derived from reanalysis RH products also reported a
over-estimate observations (Ishak et al., 2010; Vidal et al., 2010; Decker tendency for under-estimation of ea Eq. (4) (Berg et al., 2003; Simmons
et al., 2010), including when adopting short time steps for the analysis et al., 2010). Differently, Martins et al. (2017) reported for Iberia no
(Srivastava et al., 2015). Dee et al. (2011) reported that the tendency tendency for over- or under-estimation of RH when using the blended
for reanalysis global solar radiation in ERA-Interim to overestimate NCEP/NCAR reanalysis monthly products.
observations is partly due to programming errors in the calculation of Daily u2 REAN have shown a quite low accuracy. Differently from the
incident solar radiation (2 W m−2) and to considering a constant value previous analyzed variables, u2 REAN values reveals a very high dis-
of 1370 W m−2 for the solar irradiance without taking into account the persion along the regression line (Fig. 5), with low R2 = 0.26. Results
solar cycle. also show a general tendency for over-estimation of u2 (b0 = 1.08),
The Tdew REAN results (Fig. 4) indicate a tendency for under-esti- with only 8% of locations with b0 in the range 0.95-1.05. Estimation
mation, with b0 ≤ 0.95 in 75% of locations. Most cases present high R2 errors were large, with NRMSE > 30% for all cases, resulting in
(0.70–0.90), with only 4% of cases having R2 < 0.60. PBIAS ranged MSE > σ for 75% of cases. The low accuracy of wind speed derived
from −10.0 to +10.0% in 25% of cases but in 20% of locations from reanalysis was reported in several other studies applied to Por-
PBIAS > 20%. Estimation errors were relatively high, with NRMSE tugal (Carvalho et al., 2014), and the Iberian Peninsula (Lorente-Plazas
> 30% for 46% of locations. Good EF results (EF > 0.70) were found et al., 2015; Martins et al., 2017).

349
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Table 4 NCEP/NCAR blended reanalysis for Portugal (Martins et al., 2017).

Mean and range values of the accuracy indicators relative to the validation ETo Moreover, other studies using ERA-Interim reanalysis data with small
REAN computation with and without bias correction aimed at operational use time steps also reported an overall over-estimation bias for ETo
with calculations performed with quarterly time aggregation. (Srivastava et al., 2016). Thus, although computing ETo without bias
Performance Uncorrected bias Additive bias Multiple correction leads to acceptable results, it is appropriate to assess the
indicators correction, regression accuracy of ETo REAN computation with bias correction, as analyzed in
quarterly correction, the next section.
aggregation quarterly
aggregation

b0 Mean 1.06 0.98 0.98 3.3. Statistical evaluation of the bias correction approaches
Range 0.88 to 1.31 0.89 to 1.05 0.85 to 1.05
R2 Mean 0.89 0.89 0.91 The mean performance indicators (RMSE and absolute values of
Range 0.83 to 0.95 0.89 to 0.95 0.87 to 0.95
PBIAS (%) Mean 12.4 1.4 −0.6
PBIAS) relative to the cross-validation applied to the bias correction of
Range −5.4 to 37.2 −12.0 to 9.6 −9.5 to 14.7 the daily computation of ETo using reanalysis weather data of the
RMSE (mm Mean 0.78 0.63 0.57 nearest grid point are presented in Fig. 6a. The additive bias correction
day−1) were computed using the three aggregation periods t previously re-
Range 0.53 to 1.18 0.48 to 0.87 0.44 to 0.89
ferred (Fig. 2): monthly, quarterly, and the entire period. As expected,
EF Mean 0.80 0.88 0.90
Range 0.33 to 0.93 0.75 to 0.95 0.85 to 0.95 results in Fig. 6a show that computing the bias correction monthly
outperforms the other considered periods, and that quarterly compu-
b0 is the regression coefficient of the FTO; R2 is the coefficient of determination tations lead to better results than using the entire period. It is likely that
of the OLS regression; PBIAS is the Percent Bias; RMSE is the Root Mean Square RMSE and PBIAS result smaller when the additive correction is com-
Error; EF is the efficiency of modelling. puted for a shorter period that better considers the seasonality of cli-
mate. Therefore, both period lengths and seasonality combine. Similar
3.2. Assessing the accuracy of ETo REAN estimates using reanalysis data cross-validation results (Fig. 6b) were obtained when the bias correc-
without bias correction tion was performed using the multiple regression relating various grid
points with the target location (Fig. 2), i.e., computations considering
PM-ETo (Eq. 1) was computed with reanalysis data relative to the the monthly aggregation performed better than quarterly and the entire
nearest grid point (Fig. 2). Though most variables have shown various period. These results confirm the assumption that a bias correction
levels of accuracy, as described in the previous sections, the perfor- produce better effects when the seasonality of climate is considered.
mance of ETo REAN estimates when compared with ETo OBS is generally Comparing Figs. 6a and b show that using the multiple regression to
good as shown in Table 2 where the frequencies of various statistical perform bias correction leads to reduced RMSE and PBIAS relative to
accuracy indicators are presented. A strong correlation between ETo OBS using an additive bias correction. Those results are explained by the fact
and ETo REAN was observed for all locations (R2 > 0.80). Overall, these that the former correction refers to various nearby locations while the
results show that ETo REAN explain well the variability of ETo OBS. additive correction bases upon the nearest point only. However,
However, a tendency for over-estimating ETo was observed, with adopting the multiple regression is much more demanding in terms of
b0 > 1.15 in 21% of locations and PBIAS > 10% in 50% of cases. This computation as it also happens when using computations relative to the
tendency for over-estimation of ETo is likely due to the overestimation monthly aggregation. Requirements for accuracy have to be carefully
of Rs and wind speed for various stations (Fig. 4 and 5) as analyzed in considered by the user when selecting the most appropriate bias cor-
the previous section. It results an acceptable but high RMSE, larger than rection approach.
0.75 mm day−1 in 58% of locations and an overall average RMSE of The accuracy of computing daily ETo REAN using both bias correction
0.79 mm day−1. EF are generally high, with EF > 0.70 in 84% of lo- approaches is summarized in Table 3, where mean performance in-
cations, thus indicating that the mean square errors were generally dicators and the respective range of values are presented. Accordingly,
much smaller than the ETo OBS variance. (see also Fig. 6), it can be concluded that a bias correction through the
Results for daily uncorrected ETo REAN are more accurate than those multiple regression approach outperforms the use of an additive bias
previously obtained also using ERA-Interim without bias correction correction computed from the nearest grid point. However, the accu-
(Paredes et al., 2017b), which is likely due to the altitude correction racy results are not very different: all R2 averages are close to 0.90, all
used in this study. Results in Table 2 also agree with those of a previous absolute PBIAS average about 3.5%, all RMSE averages are around
study assessing the accuracy of monthly ETo REAN computed with 0.60 mm day−1, and all EF averages are close to 0.90. The ranges of

Fig. 8. Spatial distribution of the root mean squared error (RMSE) measuring the performance of ETo REAN estimation with uncorrected bias, additive bias correction
and a multiple regression correction with quarterly aggregation of data.

350
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

these indicators are shorter and closer to the optimal targets in case of approach (PMT) and the Hargreaves-Samani equation (HS). This is the
the multiple regression, when compared to the additive correction. case of a study by Paredes and Rodrigues (2010) that assessed the ac-
Also, the indicators computed when data are aggregated to the entire curacy of PMT for 8 locations in Portugal, as well as various studies
period are less good than those using quarterly or monthly aggregations focusing on the use of the HS equation in several Mediterranean loca-
(Table 3). tions (e.g., Gavilán et al., 2006; Jabloun and Sahli, 2008; Berti et al.,
Despite the improvement of results when aggregating data at the 2014). The accuracy of the current approaches revealed particularly
monthly time scale, the accuracy indicators are not appreciably better better when referring to arid, semi-arid and sub-humid climates.
than those obtained with a quarterly aggregation. Thus, the quarterly
aggregation seems to adequately capture the intra-annual seasonality of 4. Conclusions
climate and of ETo and consequently may be used for bias correction in
applications to Portugal. The accurate estimation of daily ETo is essential but data sets re-
Moreover, bias correction computed from a multiple regression as- ferring to all weather variables required to compute ETo are often not
sociating various nearby grid points may provide better results. available. Then, evaluating the potential of using daily ERA-Interim
However, due to its simple usage and since related accuracy indicators reanalysis data to compute the PM-ETo in alternative to the use of re-
compare well with those due to the multiple regression calculation, the duced data sets is needed. This evaluation refers to assess the ability of
use of the additive bias correction computed with the quarterly ag- reanalysis products to reproduce the temporal variability of the ob-
gregation may also be appropriate. A more focused comparison is served climatic variables used to compute ETo, as well as the accuracy
therefore performed in the next section. To do so, considering that cross of daily ETo REAN in reproducing the PM-ETo computed with full data
validation produces different bias corrections according to the data set sets across continental Portugal. In addition, two approaches for bias
used for training, that analysis was performed computing the bias correcting ETo REAN were assessed.
correction for the most recent half set of data because it may provide for Results show that Rs, the main driving variable of ETo, can be ac-
climate variability closer to present climate and related climate curately estimated from reanalysis with only a slight trend for over-
warming (Cardoso et al., 2018). estimation of observations. High correlations between reanalysis and
observed temperature variables - Tmax, Tmin and Tdew – were obtained
3.4. Operational selection of the bias correction methodology after their correction for altitude. However, a tendency for under-esti-
mation of both Tmax and Tdew was observed, together with a tendency
The selection of the bias correction approach and the computation for over-estimation of Tmin. Nevertheless, errors were not high.
of the related bias correction parameters was based upon the accuracy Differently, the reanalysis wind speed estimation has shown to be less
indicators computed for the validation sets when using the most recent accurate as it is often reported in literature relative to this variable.
base-period sets for calibration/training. Results comparing the accu- Despite the limitations in representing those weather variables the
racy indicators relative to computing ETo REAN without bias correction resulting errors for computing daily PM-ETo are quite small when bias
and with both the additive correction and multiple regression calcu- corrections were adopted, with RMSE not exceeding 0.87 mm day−1
lated using the quarterly aggregation are presented in Fig. 7. Results when using an additive bias correction computed from the nearest grid
show that both ETo-bias correction approaches highly improved the ETo point and data aggregated quarterly to consider seasonality effects.
REAN accuracies. Notably, after bias correction, most values of b0 are RMSE does not exceed 0.89 mm day−1 when bias correction adopted a
concentrated in the range 0.95–1.05, thus eliminating the over-esti- multiple regression with nearby grid points and data were also ag-
mation tendency previously observed. All other indicators also im- gregated quarterly. When data are aggregated monthly errors are re-
proved but less evidently. Errors for both bias corrected computations duced but computation requirements highly increase. However, cross-
are relatively small with RMSE < 0.75 mm day−1 in 88% and 92% of validation results did not show noticeable accuracy differences between
locations when using, respectively, the additive bias correction and the using quarterly or monthly aggregation periods.
multiple regression. Selection between approaches was based upon easiness of use, ac-
The average values and range of the accuracy indicators are pre- curacy of ETo REAN estimates and ability of capturing the intra-annual
sented in Table 4, which clearly shows the advantages of the bias cor- seasonality of ETo. Thus, aiming at local irrigation scheduling opera-
rection in reducing the over-estimation bias and reducing RMSE and the tional purposes, adopting an additive bias correction based upon re-
respective range. As previously referred, results of both bias correction analysis data relative to the nearest grid point and using a quarterly
approaches do not differ much. Therefore, because the additive bias aggregation may be an adequate approach.
correction is more straightforward in its application to a single location,
it may be selected for operational irrigation scheduling purposes. Acknowledgements
The spatial variability of the RMSE over Portugal (Fig. 8) show that
higher values of RMSE of the ETo REAN uncorrected for bias are mostly The support of the Fundação para a Ciência e a Tecnologia,
concentrated in the coastal areas. These results agree with those re- Portugal, through the Post-Doc grant SFRH/BPD/102478/2014 to the
ported by Martins et al. (2017) when using the NCEP/NCAR blended first author, the PhD grant SFRH/BD/92880/2013 to the second author,
reanalysis data. As already pointed out, both bias correction methods and the grants attributed to the research units LEAF (UID/AGR/04129/
allowed to improve RMSE results for most locations particularly in the 2013) and CEAUL (UID /MAT/00006/2013) are acknowledged.
northern part of the country and in the coastal areas (Fig. 8). Higher
estimation errors were found in interior locations influenced by aridity. References
However, RMSE differences between bias correction methods were re-
vealed small. Allen, R.G., Smith, M., Perrier, A., Pereira, L.S., 1994. An update for the definition of
Both bias correction approaches improve daily ETo estimations, thus reference evapotranspiration. ICID Bull. 43 (2), 1–34.
Allen, R.G., Pereira, L.S., Raes, D., Smith, M., 1998. Crop Evapotranspiration. Guidelines
resulting a mean RMSE of 0.63 and 0.57 mm day−1 respectively for the for Computing Crop Water Requirements, FAO Irrigation and Drainage Paper 56.
additive bias correction and the multiple regression correction. These FAO, Rome, Italy, pp. 300.
results are generally better than those reported by Martins et al. (2017) Allen, R.G., Pereira, L.S., Howell, T.A., Jensen, M.E., 2011. Evapotranspiration informa-
tion reporting: I. Factors governing measurement accuracy. Agric. Water Manage. 98
for several Portuguese stations when assessing the accuracy of monthly (6), 899–920.
ETo estimations using NCEP/NCAR blended reanalysis data. The results Almorox, J., Senatore, A., Quej, V.H., Mendicino, G., 2018. Worldwide assessment of the
in the current study are better or comparable to those of several studies Penman–Monteith temperature approach for the estimation of monthly reference
evapotranspiration. Theor. Appl. Climatol. 131, 693–703.
using temperature only for computing ETo, the PM temperature

351
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

Arlot, S., Celisse, A., 2010. A survey of cross-validation procedures for model selection. Legates, D.R., McCabe Jr., G.J., 1999. Evaluating the use of “goodness-of-fit” measures in
Stat. Surv. 4, 40–79. hydrologic and hydroclimatic model validation. Water Resour. Res. 35, 233–241.
Azorin-Molina, C., Vicente-Serrano, S.M., Sanchez-Lorenzo, A., McVicar, T.R., Morán- Lewis, C.S., Geli, H.M.E., Neale, C.M.U., 2014. Comparison of the NLDAS weather forcing
Tejeda, E., Revuelto, J., El Kenawy, A., Martín-Hernández, N., Tomas-Burguera, M., model to agrometeorological measurements in the western United States. J. Hydrol.
2015. Atmospheric evaporative demand observations, estimates and driving factors 510, 385–392.
in Spain (1961–2011). J. Hydrol. 523, 262–277. Lorente-Plazas, R., Montávez, J.P., Jerez, S., Gómez-Navarro, J.J., Jiménez-Guerrero, P.,
Baigorria, G.A., Jones, J.W., Shin, D.-W., Mishra, A., O’Brien, J.J., 2007. Assessing un- Jiménez, P.A., 2015. A 49 year hindcast of surface winds over the Iberian Peninsula.
certainties in crop model simulations using daily bias-corrected Regional Circulation Int. J. Climatol. 35, 3007–3023.
Model outputs. Clim. Chang. Res. Lett. 34, 211–222. Maraun, D., 2013. Bias correction, quantile mapping, and downscaling: revisiting the
Berg, A.A., Famiglietti, J.S., Walker, J.P., Houser, P.R., 2003. Impact of bias correction to inflation issue. J. Climate 26, 2137–2143.
reanalysis products on simulations of North American soil moisture and hydrological Martins, D.S., Paredes, P., Raziei, T., Pires, C., Cadima, J., Pereira, L.S., 2017. Assessing
fluxes. J. Geophys. Res. 108 (D16), 4490. https://doi.org/10.1029/2002JD003334. reference evapotranspiration estimation from reanalysis weather products. An ap-
Berti, A., Tardivo, G., Chiaudani, A., Rech, F., Borin, M., 2014. Assessing reference eva- plication to the Iberian Peninsula. Int. J. Climatol. 37, 2378–2397.
potranspiration by the Hargreaves method in north-eastern Italy. Agric. Water Maurer, E.P., Das, T., Cayan, D.R., 2013. Errors in climate model daily precipitation and
Manage. 140, 20–25. temperature output: time invariance and implications for bias correction. Hydrol.
Betts, A.K., Beljaars, A.C.M., 2017. Analysis of near-surface biases in ERA-Interim over the Earth Syst. Sci. Discuss. 17, 2147–2159.
Canadian Prairies. J. Adv. Model. Earth Syst. 9, 2158–2173. McVicar, T.R., Van Niel, T.G., Li, L.-T., Hutchinson, M.F., Mu, X.-M., Liu, Z.-H., 2007.
Bojanowski, J.S., Vrieling, A., Skidmore, A.K., 2014. A comparison of data sources for Spatially distributing monthly reference evapotranspiration and pan evaporation
creating a long-term time series of daily gridded solar radiation for Europe. Sol. considering topographic influences. J. Hydrol. 338, 196–220.
Energy 99, 152–171. Mo, K.C., Long, L.N., Xia, Y., Yang, S.K., Schemm, J.E., Ek, M., 2011. Drought indices
Boulard, D., Castel, T., Camberlin, P., Sergent, A.-S., Bréda, N., Badeau, V., Rossi, A., Pohl, based on the climate forecast system reanalysis and ensemble NLDAS. J.
B., 2016. Capability of a regional climate model to simulate climate variables re- Hydrometeorol. 12, 181–205.
quested for water balance computation: a case study over northeastern France. Clim. Monteith, J.L., 1965. Evaporation and environment. The State and Movement of Water in
Dyn. 46, 2689–2716. Living Organisms, 19th Symp. Soc. Exp. Biol. 205–234.
Cammalleri, C., Ciraolo, G., 2013. A simple method to directly retrieve reference eva- Mooney, P.A., Mulligan, F.J., Fealy, R., 2011. Comparison of ERA-40, ERA-Interim and
potranspiration from geostationary satellite images. Int. J. Appl. Earth Obs. Geoinf. NCEP/NCAR reanalysis data with observed surface air temperatures over Ireland. Int.
21, 149–158. J. Climatol. 31, 545–557.
Cardoso, R.M., Soares, P.M.M., Lima, D.C.A., Miranda, P.M.A., 2018. Mean and extreme Nash, J.E., Sutcliffe, J.V., 1970. River flow forecasting through conceptual models: part 1.
temperatures in a warming climate: EURO CORDEX and WRF regional climate high- A discussion of principles. J. Hydrol. 10, 282–290.
resolution projections for Portugal. Clim. Dyn. https://doi.org/10.1007/s00382-018- Onogi, K., et al., 2007. The JRA-25 reanalysis. J. Meteorol. Soc. Jpn. 85, 369–432.
4124-4. Paredes, P., Rodrigues, G.C., 2010. Necessidades de água para a rega de milho em
Carvalho, D., Rocha, A., Gómez-Gesteira, M., Silva Santos, C., 2014. WRF wind simulation Portugal continental considerando condições de seca. In: Pereira, L.S., Mexia, J.T.,
and wind energy production estimates forced by different reanalyses: comparison Pires, C.A.L. (Eds.), Gestão do Risco em Secas. Métodos, tecnologias e desafios.
with observed data for Portugal. Appl. Energy 117, 116–126. Edições Colibri e CEER, Lisboa, pp. 301–319.
Daly, C., Smith, J.I., Olson, K.V., 2015. Mapping atmospheric moisture climatologies Paredes, P., Rodrigues, G.C., Alves, I., Pereira, L.S., 2014. Partitioning evapotranspiration,
across the conterminous United States. PLoS ONE 10 (10), e0141140. https://doi. yield prediction and economic returns of maize under various irrigation management
org/10.1371/journal.pone.0141140. strategies. Agric. Water Manage 135, 27–39.
De Bruin, H.A.R., Trigo, I.F., Jitan, M.A., Temesgen Enku, N., van der Tol, C., Gieske, Paredes, P., Ribeiro, A.F.S., Pires, C.A.L., Pereira, L.S., 2015. Usabilidade de previsões
A.S.M., 2010. Reference crop evapotranspiration derived from geo-stationary satellite sazonais de ensemble da temperatura e precipitação na modelação da rega e da
imagery: a case study for the Fogera flood plain, NW-Ethiopia and the Jordan Valley. produção. Aplicação a cevada e milho. In: Pires, C.A.L., Pereira, L.S. (Eds.),
Jordan. Hydrol. Earth Syst. Sci. 14, 2219–2228. Predictabilidade Sazonal de Secas - Avaliação ao nível regional e agrícola. ISA Press,
De Bruin, H.A.R., Trigo, I.F., Bosveld, F.C., Meirinka, J.F., 2016. Thermodynamically Lisboa, pp. 307–332.
based model for actual evapotranspiration of an extensive grass field close to FAO Paredes, P., Pereira, L.S., Rodrigues, G.C., Botelho, N., Torres, M.O., 2017a. Using the
reference, suitable for remote sensing application. J. Hydrometeorol. 17, 1373–1382. FAO dual crop coefficient approach to model water use and productivity of proces-
Decker, M., Brunke, M.A., Wang, Z., Sakaguchi, K., Zeng, X., Bosilovich, M.G., 2010. sing pea (Pisum sativum L.) as influenced by irrigation strategies. Agric. Water
Evaluation of the reanalysis products from GSFC, NCEP, and ECMWF using flux tower Manage. 189, 5–18.
observations. J. Climate 25, 1916–1944. Paredes, P., Martins, D.S., Cadima, J., Pires, C., Pereira, L.S., 2017b. Accuracy of daily
Dee, D.P., et al., 2011. The ERA-Interim reanalysis: configuration and performance of the PM-ETo estimations with ERA-Interim reanalysis products. Europ. Water 59, 239–246.
data assimilation system. Q. J. R. Meteorol. Soc. 137, 553–597. Paredes, P., Fontes, J.C., Azevedo, E.B., Pereira, L.S., 2018a. Daily reference crop eva-
Fang, G.H., Yang, J., Chen, Y.N., Zammit, C., 2015. Comparing bias correction methods in potranspiration with reduced data sets in the humid environments of Azores islands
downscaling meteorological variables for a hydrologic impact study in an arid area in using estimates of actual vapour pressure, solar radiation and wind speed. Theor.
China. Hydrol. Earth Syst. Sci. 19, 2547–2559. Appl. Climatol. https://doi.org/10.1007/s00704-017-2329-9.
Farahani, H.J., Howell, T.A., Shuttleworth, W.J., Bausch, W.C., 2007. Evapotranspiration: Paredes, P., Fontes, J.C., Azevedo, E.B., Pereira, L.S., 2018b. Daily reference crop eva-
progress in measurement and modeling in agriculture. Trans. ASABE 50 (5), potranspiration in the humid environments of Azores islands using reduced data sets.
1627–1638. Accuracy of FAO PM temperature and Hargreaves-Samani methods. Theor. Appl.
Fu, G., Charles, S.P., Timbal, B., Jovanovic, B., Ouyang, F., 2016. Comparison of Climatol. https://doi.org/10.1007/s00704-017-2295-2.
NCEP‐NCAR and ERA‐Interim over Australia. Int. J. Climatol. 36, 2345–2367. Pereira, L.S., 2017. Water, agriculture and food: challenges and issues. Water Resour.
Gavilán, P., Lorite, I.J., Tornero, S., Berengena, J., 2006. Regional calibration of Manage. 31, 2985–2999.
Hargreaves equation for estimating reference ET in a semiarid environment. Agric. Pereira, L.S., Perrier, A., Allen, R.G., Alves, I., 1999. Evapotranspiration: review of con-
Water Manage. 81, 257–281. cepts and future trends. J. Irrig. Drain. Eng. 125 (2), 45–51.
Grouillet, B., Ruelland, D., Ayar, P.V., Vrac, M., 2016. Sensitivity analysis of runoff Pereira, L.S., Cai, L.G., Hann, M.J., 2003. Farm water and soil management for improved
modeling to statistical downscaling models in the western Mediterranean. Hydrol. water use in the North China Plain. Irrig. Drain. 52, 299–317.
Earth Syst. Sci. 20, 1031–1047. Pereira, L.S., Cordery, I., Iacovides, I., 2009. Coping With Water Scarcity. Addressing the
Hempel, S., Frieler, K., Warszawski, L., Schewe, J., Piontek, F., 2013. A trend-preserving Challenges. Springer, Dordrecht, pp. 382.
bias correction – the ISI-MIP approach. Earth Syst. Dyn. 4, 219–236. Pereira, L.S., Allen, R.G., Smith, M., Raes, D., 2015a. Crop evapotranspiration estimation
Hofer, M., Marzeion, B., Mölg, T., 2012. Comparing the skill of different reanalyses and with FAO56: past and future. Agric. Water Manage. 147, 4–20.
their ensembles as predictors for daily air temperature on a glaciated mountain Pereira, L.S., Paredes, P., Rodrigues, G.C., Neves, M., 2015b. Modeling malt barley water
(Peru). Clim. Dyn. 39, 1969–1980. use and evapotranspiration partitioning in two contrasting rainfall years. Assessing
Hwang, S., Graham, W.D., Geurink, J.S., Adams, A., 2014. Hydrologic implications of AquaCrop and SIMDualKc models. Agric. Water Manage. 159, 239–254.
errors in bias-corrected regional reanalysis data for west central Florida. J. Hydrol. Popova, Z., Kercheva, M., Pereira, L.S., 2006. Validation of the FAO methodology for
510, 513–529. computing ETo with missing climatic data. Application to South Bulgaria. Irrig. Drain.
Ishak, A.M., Bray, M., Remesan, R., Han, D., 2010. Estimating reference evapo- 55, 201–215.
transpiration using numerical weather modelling. Hydrol. Process. 24, 3490–3509. Raziei, T., Pereira, L.S., 2013a. Spatial variability analysis of reference evapotranspiration
Jabloun, M., Sahli, A., 2008. Evaluation of FAO-56 methodology for estimating reference in Iran utilizing fine resolution gridded datasets. Agric. Water Manage. 126, 104–118.
evapotranspiration using limited climatic data: applications to Tunisia. Agric. Water Raziei, T., Pereira, L.S., 2013b. Estimation of ETo with Hargreaves-Samani and FAO-PM
Manage. 95, 707–715. temperature methods for a wide range of climates in Iran. Agric. Water Manage. 121,
Jaksa, W.T., Sridhar, V., Huntington, J.L., Khanal, M., 2013. Evaluation of the com- 1–18.
plementary relationship using Noah Land Surface Model and North American Ren, X., Qu, Z., Martins, D.S., Paredes, P., Pereira, L.S., 2016a. Daily reference evapo-
Regional Reanalysis (NARR) data to estimate evapotranspiration in semiarid eco- transpiration for hyper-arid to moist sub-humid climates in Inner Mongolia, China: I.
systems. J. Hydrometeorol 14, 345–359. Assessing temperature methods and spatial variability. Water Resour. Manage. 30,
Jensen, M.E., Allen, R.G. (Eds.), 2016. Evaporation, Evapotranspiration, and Irrigation 3769–3791.
Water Requirements. American Society of Civil Engineers, Reston, Virginia. https:// Ren, X., Martins, D.S., Qu, Z., Paredes, P., Pereira, L.S., 2016b. Daily reference evapo-
doi.org/10.1061/9780784414057.fm. transpiration for hyper-arid to moist sub-humid climates in Inner Mongolia, China: II.
Kanamitsu, M., Ebisuzaki, W., Woollen, J., Yang, S.H., Hnilo, J.J., Fiorino, M., Potter, Trends of ETo and weather variables and related spatial patterns. Water Resour.
G.L., 2002. NCEP– DOE AMIP-II reanalysis (R-2). Bull. Am. Meteor. Soc. 83, Manage. 30, 3793–3814.
1631–1643. Rienecker, M.M., et al., 2011. MERRA: NASA’s modern-era retrospective analysis for

352
P. Paredes et al. Agricultural Water Management 210 (2018) 340–353

research and applications. J. Climate 24, 3624–3648. Research and Forecasting model. Atmos. Sci. Lett. 28, 4419–4432.
Senay, G.B., Verdin, J.P., Lietzow, R., Melesse, A.M., 2008. Global daily reference eva- Srivastava, P.K., Islam, T., Gupta, M., Petropoulos, G., Dai, Q., 2015. WRF dynamical
potranspiration modeling and evaluation. J. Am. Water Resour. Assoc. 44, 969–979. downscaling and bias correction schemes for NCEP estimated hydro-meteorological
Sentelhas, P.C., Gillespie, T.J., Santos, E.A., 2010. Evaluation of FAO Penman–Monteith variables. Water Resour. Manage. 29, 2267–2284.
and alternative methods for estimating reference evapotranspiration with missing Srivastava, P.K., Han, D., Islam, T., Petropoulos, G.P., Gupta, M., Dai, Q., 2016. Seasonal
data in Southern Ontario, Canada. Agric. Water Manage. 97, 635–644. evaluation of evapotranspiration fluxes from MODIS satellite and mesoscale model
Sheffield, J., Ziegler, A.D., Wood, E.F., Chen, Y., 2004. Correction of the high-latitude rain downscaled global reanalysis datasets. Theor. Appl. Climatol. 124, 461–473.
day anomaly in the NCEP–NCAR reanalysis for land surface hydrological modeling. J. Strong, C., Khatri, K.B., Kochanski, A.K., Lewis, C.S., Allen, L.N., 2017. Reference eva-
Climate 17, 3814–3828. potranspiration from coarse-scale and dynamically downscaled data in complex
Sheffield, J., Goteti, G., Wood, E.F., 2006. Development of a 50-year high-resolution terrain: sensitivity to interpolation and resolution. J. Hydrol. 548, 406–418.
global dataset of meteorological forcings for land surface modeling. J. Climate 19, Terink, W., Hurkmans, R.T.W.L., Torfs, P.J.J.F., Uijlenhoet, R., 2010. Evaluation of a bias
3088–3111. correction method applied to downscaled precipitation and temperature reanalysis
Sillmann, J., Kharin, V.V., Zwiers, F.W., Zhang, X., Bronaugh, D., Donate, M.G., 2014. data for the Rhine basin. Hydrol. Earth Syst. Sci. 14, 687–703.
Evaluating model-simulated variability in temperature extremes using modified Todorovic, M., Karic, B., Pereira, L.S., 2013. Reference evapotranspiration estimate with
percentile indices. Int. J. Climatol. 34, 3304–3311. limited weather data across a range of Mediterranean climates. J. Hydrol. 481,
Simmons, A.J., Willett, K.M., Jones, P.D., Thorne, P.W., Dee, D.P., 2010. Low-frequency 166–176.
variations in surface atmospheric humidity, temperature, and precipitation: Tomas-Burguera, M., Vicente-Serrano, S.M., Grimalt, M., Beguería, S., 2017. Accuracy of
Inferences from reanalyses and monthly gridded observational data sets. J. Geophys. reference evapotranspiration (ETo) estimates under data scarcity scenarios in the
Res. 115, D01110. https://doi.org/10.1029/2009JD012442. Iberian Peninsula. Agric. Water Manage. 18, 103–116.
Smith, T.M., Reynolds, R.W., Peterson, T.C., Lawrimore, J., 2008. Improvements to Urraca, R., Martinez-de-Pison, E., Sanz-Garcia, A., Antonanzas, J., Antonanzas-Torres, F.,
NOAA’s historical merged land-ocean surface temperature analysis (1880–2006). J. 2017. Estimation methods for global solar radiation: case study evaluation of five
Climate 21, 2283–2296. different approaches in central Spain. Renew. Sust. Energ. Rev. 77, 1098–1113.
Soares, P.M.M., Cardoso, R.M., Miranda, P.M.A., de Medeiros, J., Belo-Pereira, M., Vidal, J.-P., Martin, E., Franchistéguy, L., Baillon, M., Soubeyroux, J.-M.A., 2010. 50-year
Espirito-Santo, F., 2012. WRF high resolution dynamical downscaling of ERA-Interim high-resolution atmospheric reanalysis over France with the Safran system. Int. J.
for Portugal. Clim. Dyn. 39, 2497–2522. Climatol. 30, 1627–1644.
Srivastava, P.K., Han, D., Rico-Ramirez, M.A., Islam, T., 2013. Comparative assessment of Zhao, T., Guo, W., Fu, C., 2008. Calibrating and evaluating reanalysis surface temperature
evapotranspiration derived from NCEP and ECMWF global datasets through Weather error by topographic correction. J. Climate 21, 1140–1446.

353