Spatial and seasonal variations in evapotranspiration over Canada ’ s landmass

A 30 yr (1979–2008) dataset of actual evapotranspiration (ET) at 1 km resolution was generated over Canada’s landmass by integrating remote sensing land surface data and gridded climate data using the EALCO model run at a 30 min time step. This long-term high-resolution dataset was used to characterize the spatiotemporal variations in ET across Canada. The results show that annual ET varied from 600 mm yr −1 over several regions in the south to less than 100 mm yr −1 in the northern Arctic. Nationally, ET in summer (i.e., June to August) comprised 65 % of the annual total amount. ET in the cold season remained mostly below 10 mm month −1 over the country. Negative monthly ET was obtained over the Arctic region in winter, indicating EALCO simulated a larger amount of condensation than ET. Overall, the mean ET over the entire Canadian landmass for the 30 yr was 239 mm yr −1, or 44 % of its corresponding precipitation. Comparisons of available ET studies in Canada revealed large uncertainties in ET estimates associated with using different approaches. The scarcity of ET measurements for the diverse ecosystems in Canada remains a significant challenge for reducing the uncertainties; this gap needs to be addressed in future studies to improve capabilities in climate/weather modeling and water resource management.


Introduction
Evapotranspiration (ET) is the water lost from the land surfaces to the atmosphere.Annually and worldwide, ET returns about two-thirds of land-based precipitation (P ) to the atmosphere (Trenberth et al., 2007).Subtracting ET from P results in the available water to a hydrologic system for recharge to subsurface aquifers and to streamflows (Q).Determined by the climatic as well as land surface conditions, ET is a strong indicator of water availability in response to changes in climate and land use.ET is an important component in both the land surface energy and water budgets.The amount and timing of ET can strongly affect the atmosphere and surface/subsurface processes such as cloud development (Molders and Raabe, 1996), surface temperature (Taha, 1997;Wang, 2008), streamflow (Koster and Milly, 1997), groundwater recharge (Renger et al., 2007;Wang et al., 2011;Githui et al., 2012) and the ecosystem carbon cycle (Wang et al., 2002a).Thus, the characterization of the spatiotemporal variations of ET is critical to sustainable management of water resources, climate/weather modeling, and environmental applications (Shukla and Mintz, 1982;Cox et al., 1999;Wang et al., 2009b;Gerosa, 2011, Irmak, 2012).
Canada covers a land area of about 9.1 million km 2 , stretching from 42 • N in the south to 83 • N in the north, and from 52.5 • W in the east to 141 • W in the west.Quantifying ET for this large area has been challenged by data limitations and lack of understanding of ET processes.Traditionally, the surface water budget approach that calculates ET = P − Q is often used (e.g., Hare, 1980;Donohue et al., 2007).This approach is mainly used for long-term (multi-year or longer) ET estimates when total water storage change ( TWS, including surface water bodies, snow, glaciers, plant and soil water, and groundwater) involves negligible quantities compared to the accumulated P , ET and Q.It assumes no cross-boundary sub-surface water flow and treats one drainage basin as a single unit and recognizes no Published by Copernicus Publications on behalf of the European Geosciences Union.
ET spatial distribution within the basin.These limit its usefulness in understanding the detailed spatial and temporal dynamics of ET.Moreover, this approach is greatly challenged by the uncertainties in the basin-level P and Q and the longterm trends in TWS found over some Canadian basins.
Atmospheric models simulate ET consistently over broad spatial scales at high temporal resolutions.With advancement in data assimilation techniques, significant improvements in modeling accuracy have been achieved in the past decades.In Canada, ET from several major atmosphere models and reanalyses has been analyzed in the Water and Energy Budget Studies (WEBS, Szeto et al., 2008;Szeto, 2007).It should be noted that ET from atmosphere models and reanalyses is a purely forecasted variable that is more poorly constrained by observations than other variables such as air temperature (T a ), humidity and wind.Large-scale ET can also be estimated by an atmospheric moisture budget approach.For example, Serreze et al. (2003) estimated ET using the difference of P and the moisture flux convergence calculated from National Centers for Environmental Prediction (NCEP) reanalysis data over the terrestrial Arctic drainage system.Strong (1997) estimated ET over the Canadian Prairies based on sequential radiosonde soundings and moisture advection data from the CMC Data Assimilation System.This approach has the convenience of disregarding TWS, but closing the atmospheric moisture budget is quite difficult.The magnitudes of the residuals from balancing the budgets are often comparable to the budget terms themselves, and the differences of budget estimates from using different datasets are also large, suggesting that substantial improvements to models and observations are needed before the uncertainty in ET estimates can be significantly reduced (Szeto et al., 2008).Moreover, this approach can hardly characterize ET in a manner suitable for water resources management.
Remote sensing methods for ET estimation have been developed with varying degrees of success (Courault et al., 2005;Kalma et al., 2008;Li et al., 2009;Kustas and Anderson, 2009;Van Niel et al., 2011;Anderson et al., 2012).These methods differ greatly in the degree of physical realism, but they all inevitably face two major limitations.First, remote sensing itself can only provide partial information to constrain the ET physical process.Second, for high-latitude regions such as northern Canada, they have to rely on polar-orbiting satellite passive optical sensors (i.e., reflective and thermal), which, most of the time, can only image a given location once in several days, or even longer due to cloud and other contaminants and the polar night season.ET is a highly dynamic process that may have dramatic sub-daily changes with weather conditions.In order to obtain ET continuously over a period of time, these methods have to make some major assumptions and simplifications and rely heavily on data from additional sources as well as empirical relationships that require local tuning in applications.In fact, even instantaneous estimates of ET at satellite overpass time, remote sensing methods may still have greater variability and uncertainty when compared with physically based models (Kite and Droogers, 2000).However, in contrast to the water budget approaches discussed above, remote sensing methods are unique in providing high spatial resolution estimates of ET.
Physically based land surface models can give reliable ET estimates when they are fairly constrained (e.g., Amthor et al., 2001;Wang et al., 2002a;Hanson et al., 2004;Grant et al., 2005Grant et al., , 2006;;Wang, 2008), but this is often difficult for large-scale applications particularly in the high-latitude regions where observations are sparse.Thus, the integration of physically based models with remotely sensed data has become a promising approach to estimate large-scale ET (Mu et al., 2011;Sheffield et al., 2010;Cleugh et al., 2007).The first process-based model estimate of ET for Canada (Liu et al., 2003) used the Boreal Ecosystem Productivity Simulator (BEPS) model driven by Advanced Very High Resolution Radiometer (AVHRR)-derived vegetation data and one year of climate model forecast at NCEP.BEPS was a simplified bucket model with a daily time step and did not include the full water and energy closures of the land surface.The meteorological data used in their study were found to have large biases relative to observations.For example, a comparison at 96 climate stations across Canada showed that its daily total radiation was 20-40 % higher than the station measurements, and the correlation coefficient for daily P between the two datasets was only 0.62 (Liu et al., 2003).
In a later study, Fernandes et al. (2007) simulated Canadawide ET using the Ecological Assimilation of Land and Climate Observations (EALCO) model, land surface parameters derived from remote sensing, and observed long-term (30 yr) meteorological data at 101 stations where hourly surface climate observations are available.EALCO is a comprehensive land surface model that was developed with a focus on cold region processes and calibrated and validated using in situ measurements (e.g., tower fluxes) in Canada.EALCO simulates ET by solving the coupled radiation, energy, and water balance equations at a 30 min time step.Plant carbon and nitrogen processes are dynamically coupled with energy and water processes to simulate the plant physiological control on ET.The use of observed meteorological data in this study, as opposed to climate model outputs, also reduced the uncertainties in the ET estimates.However, this study only covered limited regions where climate observations are available.
In this study, the 30 yr  ET for the entire Canadian landmass was simulated at 1×1 km 2 spatial resolution by the EALCO model driven by remote sensing land surface data and gridded meteorological forcing products.The objectives are to characterize the spatial and seasonal variations of ET at Canada national scale and to provide water managers with long-term high-resolution ET datasets for understanding Canada's water resources and for policy development.

The EALCO model
EALCO is a physically based numerical model developed to simulate the ecological and hydrological processes of terrestrial ecosystems using earth observations.It runs at a 30 min time step and includes five major modules (Fig. 1) that are dynamically coupled with each other to simulate the land surface: (i) radiation transfer (Wang, 2005;Wang et al., 2007), (ii) energy balance (Wang et al., 2002a;Zhang et al., 2008), (iii) water dynamics including surface watergroundwater interactions (Wang, 2008(Wang, , 2012;;Wang et al., 2009a, b), and (iv) carbon and (v) nitrogen biogeochemical cycles (Wang et al., 2001(Wang et al., , 2002b)).EALCO has been calibrated and validated throughout using tower-flux measurements and watershed-level water budget approach, which were reported in the above literature as well as in Amthor et al. (2001), Potter et al. (2001), Hanson et al. (2004), Grant et al. (2005Grant et al. ( , 2006)), Fernandes et al. (2007), Mi et al. (2007Mi et al. ( , 2009)), Li et al. (2011), Widlowski et al. (2011), and Yang et al. (2012).In this study, we report the large-scale characterization of ET by EALCO across Canada's landmass and validation by comparing with other existing studies.Below is a brief overview of EALCO focusing on ET calculations.
The EALCO water module (Fig. 1, top panel) simulates the major water transfer processes in a terrestrial ecosystem, which include ET, canopy interception of rain or snow, dew and frost formations, infiltration at soil and snow surfaces, water transfer in the soil and snow profiles and the phase changes (freeze and thaw), surface runoff, snow cover dynamics, water transfer in the soil-root-stem-leafatmosphere continuum, and the diffuse exchange of soil water and groundwater.ET is obtained by simulating its components, including canopy transpiration through leaf stomata, canopy evaporation/sublimation of intercepted rain/snow or dew/frost, soil evaporation and snow sublimation at the ground surface, and evaporation of water from temporal puddles after rain events.Canopy interception of rain or snow is determined by the amount of P , canopy gap probability and canopy loading capacity.Dew and frost formations on the canopy or ground surface are simulated when the canopy or ground surface temperature drops below the dew point temperature of the air.Evaporation or sublimation of canopy intercepted water or snow occurs before transpiration and is simulated by solving the canopy energy balance equation.The transpiration algorithm is based on solving the governing equation system that represents the coupled canopy energywater-CO 2 transfer dynamics using the nested convergence approach.This algorithm dynamically integrates the roles of atmospheric demand, plant and soil hydraulic constraints, and the leaf stomata-related physiological activities in the ET process.When rain reaches the ground surface, infiltration is first calculated and if it is smaller than the P rate, the extra amount of water is temporally accumulated on the soil surface as a puddle until it exceeds the maximum water retention depth as surface runoff.Soil evaporation or snow sublimation on the ground surface is simulated by solving the soil/snow surface energy balance equations.The soil water movement in the soil profile is based on implicit solutions of the coupled soil heat and water transfer equations using a finite difference method for the multi-layered soils (7 layers were used in this study).The soil water is linked with groundwater (water table) through diffusive groundwater recharge and discharge.Water transfer in the soil-rootstem-leaf-atmosphere continuum is based on the hydraulic gradients and resistance in the system, which include the calculations of sunlit vs. shaded canopy stomatal resistances, multi-layer hydraulic resistances for soil and roots (both radial and axial), and plant water capacitance.
The other four modules in EALCO (Fig. 1, bottom panel) are dynamically coupled with the water module to simulate the land surface radiation, energy, C and N processes.Outputs from these modules directly or indirectly affect ET.For example, the radiation module not only provides driving energy for ET, but it also provides inputs for calculating other processes such as soil/snow thermal states and plant physiological activities, which indirectly affect ET.EALCO's radiation transfer scheme uses the gap probability-based successive orders of scattering approach to simulate the land surface albedo and radiation absorption by canopy, soil and snow at a user-specified number of wavelengths or wavebands (two broad wavebands of visible and near infrared were used in this study).Recent model validation studies have demonstrated the robustness of this algorithm (e.g., Widlowski et al., 2011).The energy module is also crucial to ET determination as the energy-related processes such as soil freeze and thaw, snow accumulation and melt, and plant phenological cycles play key roles in controlling the water processes such as evaporation and transpiration, surface infiltration and runoff, water movement in the soil profile, and groundwater recharge over high-latitude regions.EALCO simulates the energy states, fluxes, and storage changes of canopy, soil and snow through numerically solving the energy balance and transfer equations for each of them.The C and N modules simulate the plant physiological processes and their physiological control on ET.For example, the canopy stomatal conductance in EALCO is driven by leaf photosynthesis obtained in the C and N modules.The plant root hydraulic resistance and capacitance are determined by total root length and its distribution in the soil profile, which is also obtained in the C and N modules.
EALCO recognizes up to five surface types at sub-pixel level (Fig. 2): (i) vegetation with soil, (ii) vegetation with snow-covered ground, (iii) bare soil, (iv) snow-covered ground, and (v) surface water body.The snow cover fraction in a pixel is a dynamic variable simulated at each time step.The water body fraction in a pixel is prescribed and remains unchanged during the model run.Within a time step, EALCO runs separately for each of the sub-pixel surface types, after which the model outputs (e.g., ET) are scaled up to the pixel level using their areal fractions.

Land cover (LC)
LC maps are used to assign vegetation types in EALCO.The Canada national scale LC (Fig. 3) used in this study was extracted from a land cover map of North and Central America for the year 2000 (Latifovic et al., 2004).This product is at 1 km resolution and based on data acquired by the VEGETATION instrument on board the SPOT 4 satellite.It provides 28 classes based on modified Federal Geographic Data Committee/Vegetation Classification Standard classification system (Table 1).Vegetation in EALCO is categorized into six broad functional types: conifers, deciduous broadleaf, mixed conifers and deciduous, evergreen broadleaf, grassland, and crops.The vegetation type of mixed conifers and deciduous is specifically included in EALCO due to its wide occurrence in Canada.In this study, data from the LC map were regrouped into the broad vegetation types used in EALCO as shown in Table 1.

Leaf area index (LAI)
LAI maps are used together with LC for canopy parameterization in EALCO.The Canada national-scale LAI map (Fig. 4) used in this study was also retrieved from data acquired by the VEGETATION instrument (Fernandes et al., 2003).It represents the maximum LAI for each 1 × 1 km 2 pixel obtained from a series of 20 images covering the period of April-October in 2000.EALCO was developed to use remote sensing-derived LAI maps in alternate ways.When the time series of LAI maps (e.g., every 10 days) are available, EALCO can directly assimilate them to represent the seasonal variations of LAI.EALCO can also be configured to use the maximum LAI only, in which case the plant phenology and seasonal LAI changes will be simulated internally by  EALCO.In this study, the latter approach was adopted since LAI time series maps were not available for all the years in this study.

Soils
The soil texture dataset was obtained from the Soil Landscapes of Canada (SLC) database of Agriculture and Agri-Food Canada (1996).This dataset is based on existing soil survey maps that have been recompiled at a 1 : 1 million scale.Each area on the map is described by a standard set of attributes for a distinct type of soil.This SLC version 2.2 dataset is the only dataset available for the full coverage of Canada's landmass.It has undergone extensive manual review by regional pedologists and a formal set of validations.The soil texture dataset was used in EALCO for soil parameterizations such as porosity and saturated hydraulic conductivity based on literature such as Saxton and Rawls (2006).This soil-texture-based general estimate of parameter values is deemed sufficient for large-scale studies.

Digital elevation
The digital elevation was generated by combining three datasets: the Canadian Digital Elevation Data (CDED) from Natural Resources Canada, the Shuttle Radar Topography Mission (SRTM) data and the Global 30 Arc Second Elevation Data (GTOPO3) from USGS.CDED has resolution varying from 0.75 arc seconds to 3 arc seconds, depending on the geographic location.This dataset was extracted from the hypsographic and hydrographic elements of the National Topographic Data Base (NTDB) or various scaled positional data acquired from the provinces and territories, or from remotely sensed imagery.SRTM was flown aboard the space shuttle Endeavour to acquire radar data that were used to create detailed topographic maps.SRTM data cover 80 % of the earth's land surface between 60 • N and 56 • S latitude with data points posted every 3 arc seconds outside the US.Strictly, SRTM data are not a traditional digital elevation model (DEM), but rather a digital surface model (DSM).Given the scale of this study, the possible difference between SRTM and DEM is ignored.The GTOPO30 is a global dataset with a resolution of 30 arc seconds.These three datasets were combined to produce the full coverage DEM data for Canada's entire landmass based on optimum data resolution and quality (Fig. 5).The terrain slopes were then calculated from this DEM using the Geomatica software.The final product was up-scaled to 1 km resolution by averaging all the sub-pixels in the original data and used to parameterize the maximum water retention depth in EALCO as stated in Sect.2.1.

Ecozones
The Canadian ecozones (Fig. 6) were used to aggregate the results for analysis.There are 15 ecozones recognized over Canada's landmass.These ecozones represent the highestlevel ecological units recognized on a sub-continental scale.They are characterized by major abiotic and biotic factors including climate, vegetation, soil, and terrain, which are also the most relevant factors in determining ET.More detailed classifications of ecological regions in Canada are also available, but this paper will focus on the ecozone level.

Meteorological forcing data
EALCO is driven by the meteorological variables of downward shortwave and longwave radiations, T a , P , humidity, wind speed, atmospheric pressure, and CO 2 concentration.
In this study, the first seven variables were obtained from Sheffield et al. (2006).This dataset is constructed by combining a suite of global observation-based datasets with atmospheric model reanalysis and available globally at 1.0 degree and 3 h resolutions for 1948-2008.Biases in the reanalysis P and near-surface meteorology are corrected using observation-based datasets of P , T a and radiation.Corrections are also made to the rain day statistics of the reanalysis P , which have been found to exhibit a spurious wavelike pattern in high-latitude wintertime.Downward shortwave and longwave radiation have undergone trend correction and probability-weighted scaling for biases.This dataset provides a long-term, spatially consistent near-surface meteorology and has been noted to be most suitable for studying land surface hydrology in the long term and broad scales.In this study, these data were disaggregated to 10×10 km 2 grids using bilinear interpolation and to 30 min time step using the method developed in Global Soil Wetness Project (http: //www.iges.org/gswp2/util/drv_finterp.f90).The time series of CO 2 data is obtained from the Global View CO 2 (2010) data product.Spatial variations of atmosphere CO 2 concentration were not considered in the model run.

Results
The main controlling factors of ET -climate, vegetation, soil, and terrain topography -all vary substantially over this vast landmass.The overall annual mean T a is mostly under 10 • C (Fig. 7).Areas with T a > 5 • C are small and mainly distributed over the west and east coastal regions, southern Mixed Wood Plains close to the Great Lakes, and the southern edge of the Prairies.While areas with T a > 0 • C extend significantly northward, the central strip of the landmass, together with the high-elevation regions in the west Montane Cordillera and the southeast corner of the Arctic, has T a between 0 and −10 • C. The three Arctic ecozones have very low T a , which reaches below −20 • C in the northern Arctic.
The annual P varies dramatically in space and time across the landmass.It generally decreases along with latitude, ranging from over 1000 mm yr −1 in the south to less than 100 mm yr −1 in the high Arctic (Fig. 8).The Pacific and Atlantic coasts and some windward mountain slops in the cordillera receive the greatest amounts of P , which could exceed 2000 mm yr −1 .In sharp contrast, the Interior Plains region, corresponding to the Prairies and Mackenzie Lowlands, is the "rain shadow" of the cordillera, rarely receiving P over 500 mm yr −1 and frequently much less than 400 mm yr −1 .For most of the country, the maximum amounts of P occurred in summer.On the coastal and east of the Great Lakes regions, a relatively large amount of P falls in winter.April was a transitional month across the basins in southern Canada and rainfall began to replace snow as dominant form.October marked the transition from mainly rain to snowfall across northern Canada.
ET obtained for the landmass varies from over 600 mm yr −1 in the south to under 100 mm yr −1 in the north (Fig. 9).The regions with ET > 500 mm yr −1 are mainly distributed in the western Pacific Maritime, southern Montane Cordillera, southern Mixed Wood Plains and Boreal Shield in southern Ontario, and southeastern Atlantic Maritime.Some regions in the south were modeled with a depressed ET of < 300 mm yr −1 , which include the mountainous regions in Pacific Maritime and Montane Cordillera and the central part of the Prairies.ET in the central-northern part was mostly < 300 mm yr −1 , and the Arctic Archipelago region has ET < 100 mm yr −1 .Some regions in the Taiga Shield also show noticeably lower ET than the surrounding regions.
After the ET was aggregated over the 15 ecozones (Fig. 10), the Mixed Wood Plain and the Atlantic Maritime showed the highest ET at around 450 mm yr −1 , followed by the Prairie, Pacific Maritime, Boreal Plain, Boreal Shield, and Montane Cordillera, which all have ET over 300 mm yr −1 .The Prairies is the driest region in Canada, but its ET was relatively high among the 15 ecozones.This is mainly due to the high ET values in some areas such as the Lake Manitoba Plain and the Aspen Parkland.These areas have relatively high P and atmosphere evaporation demand.The Montane Cordillera also has relatively high ET due to the pronounced values in its southern part.The ecozones with ET < 200 mm yr −1 include the three Arctic ecozones, Taiga Shield, Taiga Cordillera, and Boreal Cordillera.The Arctic Cordillera has ET under 100 mm yr −1 .A large number of pixels in the northern Arctic Cordillera were not included in the model run due to the gaps in input data, so the actual value of ET for this ecozone may be even lower.
The error bars in Fig. 10 indicate the degree of spatial variation (one standard deviation, σ ) or heterogeneity of ET within an ecozone.The Pacific Maritime shows the largest value of σ = 159 mm among the 15 ecozones, with coefficient of variation (CV = σ /mean) of 45 %, largely due to the high ET values in the west coastal regions and low ET values in the east mountainous region.The Boreal Shield, which stretches across a vast geographical region with very different climates, and Montane Cordillera, which has high ET values in its south and low ET values in its north, also show high spatial variations.The Atlantic Maritime and Mixed Wood Plain also have relatively high σ , but their CVs are comparatively low as they have high mean ET values.In contrast, the Taiga Shield has a high CV but a moderate σ of ET.
Climatologically, the P -ET determines the renewable water availability over a region.The most noticeable pattern of the ET / P ratio (Fig. 11) is the high values over the rain shadow regions of the cordillera.These regions have relatively low P and high T a .In the central-south Prairies, where the driest region in Canada is located, this ratio reaches close to 1.0 and highly correlates with P , indicating water constraints on ET.On the other hand, the ET / P ratios are as low as 0.2 or less in some regions such as the west coast due to the very high P amount.Overall, the vast Taiga Shield, Boreal Shield, Hudson Plains, Atlantic Maritime, Pacific Maritime, and a large part of Montane Cordillera and Boreal Cordillera  have ET / P ratios below 0.5, indicating energy control on ET.These regions contain most of the Canadian lakes and wetlands.
The seasonal variations of ET over the landmass are dramatic, typical for high-latitude regions where climate and vegetation both have pronounced seasonal changes.ET in the cold season (October-April) is very low over the major part of the landmass (Fig. 12).Negative ET was obtained in the Arctic regions in some winter months (mostly October-December), indicating that EALCO simulated condensation > evapotranspiration.Areas with relatively high ET during the cold season are mainly found in the west and east maritime, the southern part of Montane Cordillera, and the vast boreal forest regions.It is worth noting that ET over the agriculture-and grassland-dominated regions in winter, such as in the Prairies and the Mixed Wood Plains, is lower than the forest regions to their north, although their climate conditions are more favorable for ET than those over the northern forest regions.This is mainly due to the fact that snow-covered bare ground over the agricultural/grassland regions has high surface albedo, which reduces energy absorption from solar radiation for ET and increases the chance for condensation, while the forested region has low albedo resulting in the higher radiation absorption.The snow interception by the forest canopy and its sublimation is another important process in contributing to a higher ET than the agricultural/grassland regions.ET in April and May generally increases all over the landmass.The differences of ET between the above-mentioned agricultural/grassland regions and their northern forest regions diminish, mainly due to the wet soil after snowmelt over the agricultural/grassland re-  the central-south Prairies have much lower ET than the boreal forest region to its north due to the soil water limitation.
The long-term average monthly ET for the 15 ecozones (Fig. 13) all show the maximum values in July but with large difference in the seasonal variation ranges.The Mixed Wood Plains, which has the highest annual total ET (Fig. 10), also shows the highest peak season ET in July.The Prairies have the second highest peak season ET, followed by the Atlantic Maritime (in opposite order to the annual total ET of these two ecozones).The Boreal Plain and Boreal Shield also have a very pronounced peak season ET.In contrast, the Pacific Maritime, which has relatively high annual total ET, has relatively low peak season ET in summer but the highest ET in winter months.The winter month ET for the Atlantic Maritime, Boreal Shield, and Montane Cordillera is also relatively high.The three Arctic ecozones, together with the Taiga shield, have the lowest peak season ET in summer and the smallest seasonal variations.
Differences in skewness of the monthly ET distributions can be observed among the ecozones (Fig. 13).For example, the Prairie shows a significant positive skewness, which is the result of better soil water conditions in May-June than that in August-September.In contrast, the ecozones of Atlantic Maritime, Pacific Maritime, and Hudson Plain show significant negative skewness, the result of less favorable atmospheric conditions for ET in May-June than after July.
Overall, the ET over the entire landmass is very low for about half the year due to the long cold season with dormant vegetation (Fig. 14).ET in summer has a dramatic increase and peaks in July due to the more favorable atmospheric conditions and vegetation growth.The sum of ET in June, July, and August takes about 65 % of the annual total ET at the national level.In contrast, the seasonal distribution of P is less pronounced than ET and about 60 % of the annual total amount of P occurring in the second half of the year (July through December).The ET amount is lower than the corresponding quantity of P in all of the 12 months.The minimum P -ET occurs in July, and it is about 6.5 mm month −1 , or 10 % of the P of 67 mm month −1 .The maximum P -ET occurs in September through November, which is above 30 mm month −1 .Overall, the modeled average annual ET for the entire landmass during 1979-2008 was 239 mm yr −1 .It accounted for 44 % of the corresponding P of 545 mm yr −1 .

Comparison with other studies
The ET obtained in this study is compared with existing large-scale studies in Canada.The Hydrological Atlas of Canada (1978) published the baseline ET map of Canada where ET was calculated using water surface evaporation (E) and empirical coefficients representing the ratio of actual to potential ET for various land surface types.Hare (1980) delineated the pattern of ET over Canada south of 60 • N using a surface water budget approach based on the observed fields of P and Q.Although these studies reported at very coarse scales and were based on data from very different time periods, it is found that results from this study are fairly consistent with these two products.For example, both the quantities and spatial patterns of the highest/lowest ET values from these products agree well across the landmass.The northward-protruding patterns of high ET in the rain shadow of the cordillera are also evident in these products.The ET reported in Liu et al. (2003) for Canada's landmass was 228 mm yr −1 (Fig. 15).It is close to the estimate of 239 mm yr −1 from this study despite the substantial differences in both models and meteorological input data.However, some notable differences in the spatial distributions of ET can be found.For example, Liu et al. (2003) did not show the areas with pronounced ET on Vancouver Island, the coastal region of BC, and the southern cordillera as found in this study as well as in Hare (1980) and the Hydrological Atlas of Canada (1978).On the other hand, they reported higher ET values over the Boreal Plains and Taiga Shield than the three other studies.
The Mackenzie River basin (MRB) is one of the most studied large basins in Canada.It is located in northwest Canada and is one of the major river systems of the world (Fig. 5).It covers about 1.8 million km 2 and crosses six ecozones including Taiga Plains, Taiga Cordillera, Boreal Cordillera, Taiga Shield, Boreal Plains and Boreal Shield.Szeto et al. (2008) calculated its ET from five different atmosphere models and reanalysis datasets of NCEP, ERA-40, NARR, CMC and CRCM.In general, the differences in ET resulting from the above datasets were substantial (Fig. 15).The reported average annual ET during 1997-2002 varied from 270 mm yr −1 with the CRCM to as high as 639 mm yr −1 with the NCEP.The ET suggested by NCEP was probably overestimated as its value was even higher than the estimate for its water surface evaporation (Wang and Yang, 2011).
ET for the MRB was also estimated using surface water budget approach (i.e., ET = P − Q) in three of the studies shown in Fig. 15 (i.e., Serreze et al., 2003;Szeto et al., 2008 Vinukollu et al., 2011).Despite the difference in their study periods, the discrepancy of ET among them is likely due to the difference in the P datasets.Szeto et al. (2008), who reported an ET of 292 mm yr −1 , used a gridded monthly P dataset known as CANGRID, which was based on climate station data.The P dataset used in Serreze et al. (2003), who reported an ET of 241 mm yr −1 , was also based on climate station data but processed and interpolated differently.The P dataset used in Vinukollu et al. (2011), who reported an ET of 233 mm yr −1 , was based on the Global Precipitation Climatology Center (GPCC) data.Figure 15 also shows two studies that obtained ET for the MRB using the atmospheric moisture budget approach (i.e., Serreze et al., 2003;Su et al., 2006).In this approach, ET is calculated as the difference of P and atmospheric moisture flux convergence.Again, despite the difference in the study periods, the discrepancy of ET (269 mm yr −1 in Serreze et al. (2003) and 327 mm yr −1 in Su et al., 2006) between these two studies could be due to the difference in the datasets they used.Serreze et al. (2003)  the PT-Fi model obtained similar values of ET (231 mm yr −1 and 248 mm yr −1 , respectively), and they were significantly higher than that estimated by the PM-Mu model (155 mm yr −1 ).The hydrological model Variable Infiltration Capacity (VIC) was used in two independent studies in Fig. 15, and the estimated ET was 285 mm yr −1 in Su et al. (2006) and 198 mm yr −1 in Sheffield and Wood (2007).In another study, Louie et al. (2002) used the model of Morton (1983) and reported an average basin-level ET of 277 mm yr −1 .
ET was also estimated for the Saskatchewan River basin (SRB) using atmosphere models and reanalysis datasets by Szeto (2007).The SRB is located to the south of MRB and is bounded by the Rocky Mountains to the west and the boreal forest to the north and east (Fig. 5).It covers over 400 000 km 2 , which drains most of the Prairies and a portion of Boreal Plains in the northeast.The SRB is dominated by cropland and grassland, and it is one of the most important regions for grain crop productions in Canada.The reported ET for SRB during 1997-2002 varied from 383 mm yr −1 in CRCM to as high as 763 mm yr −1 in NCEP (Fig. 15).The observed ET, averaged from three flux tower sites located in the southern boreal forest within the SRB, was 339 mm yr −1 .Again, the value suggested by NCEP was likely overestimated as its value was even higher than the estimate for its water surface evaporation (Wang and Yang, 2011).
The average ET obtained in this study was 226 mm yr −1 for MRB and 352 mm yr −1 for SRB during the 30 yr of 1979-2008.The water surface (rivers and lakes) within the basins covers a total area of about 11 % in the MRB and 6 % in the SRB.With the adjustment of water surface evaporation of 454 mm yr −1 for the MRB and 681 mm yr −1 for the SRB (Wang and Yang, 2011), the basin-level composite ET was 251 mm yr −1 for the MRB and 372 mm yr −1 for the SRB.They are very comparable to the majority of estimates from other studies (Fig. 15).Note that the above water surface fractions do not include the seasonal water ponds formed by the spring snowmelt, which is a significant phenomenon particularly in SRB.In addition, agricultural irrigation is widely S. Wang et al.: Evapotranspiration over Canada's landmass used in SRB, which would significantly increase ET but is not included in this study.The four atmosphere models and reanalysis of ERA-40, NARR, CMC and NCEP are positively biased significantly over EALCO for both MRB and SRB.As indicated earlier, ET from these datasets is a purely forecasted variable that was poorly constrained by observations.Nevertheless, given that regional ET measurements are not available for rigorous validation, the spread of the ET magnitudes among the different studies gives a measure of the degree of uncertainty in ET estimates.
Negative ET (condensation > evapotranspiration) was obtained in EALCO over the Arctic region during the cold season.The Arctic region in winter was dominated by loss of radiative energy, which results in lower land surface temperatures than the lower atmosphere.The simulated condensation by EALCO is not unreasonable; in fact, ice crystal and hoarfrost deposits on the cold snow surface are frequently reported in winter (e.g., Burns, 1974).Interestingly, negative ET was also reported in the Canadian Arctic region during the winter months by Serreze et al. (2003), Su et al. (2006), and Dai and Trenberth (2002) using the atmospheric moisture budget approach, despite the use of different datasets for calculating the atmospheric moisture convergence and regional P .In an independent study, Louie et al. (2002) also reported negative ET in winter.However, Serreze et al. (2003) and Su et al. (2006) suggested, without providing adequate evidence, that the negative ET they obtained could be due to the estimation error in either atmospheric moisture convergence or P .None of the studies included rigorous estimate of the blowing snow sublimation, which could be a significant process contributing to ET (Pomeroy and Essery, 1999).Nevertheless, the ET during the cold season is poorly validated and needs to be addressed in further studies.

Discussion and conclusion
This study characterized the spatiotemporal variations in ET over Canada's landmass for 1979Canada's landmass for -2008.The results show that ET varied remarkably both in time and space.Annually, ET exceeds 500 mm yr −1 over several regions in the south, and is below 100 mm yr −1 in the northern Arctic.For the major part of the country, ET is mainly controlled by energy and it is significantly lower than P , except some central-west regions such as the Prairies where water limitation plays a significant role in the ET process.The ET in July is the highest for all of the 15 ecozones.Nationally, ET in June, July, and August comprises about 65 % of the annual total.In the cold season, which lasts about half a year, ET remains below 10 mm month −1 for most of the ecozones.Monthly ET was found lower than the corresponding P yearround for all ecozones.Overall, the modeled average ET over the entire Canadian landmass for the 30 yr of 1979-2008 is 239 mm yr −1 , or 44 % of the corresponding P .Comparisons of regional ET with other existing studies (Fig. 15) revealed Morton's model ( 5) EALCO ( 6) CRCM ( 7) ERA-40 ( 7) CMC ( 7) NCEP ( 7) EALCO ( 6) BEPS ( 8) EALCO ( 6) Evapotranspiration (mm year  Su et al. (2006Su et al. ( ), 1979Su et al. ( -1999;; (5) Louie et al. (2002Louie et al. ( ), 1961Louie et al. ( -1990;;(6) this study, 1979-2008;(7) Szeto (2007(7) Szeto ( ), 1997(7) Szeto ( -2002;;and (8) Liu et al. (2003), 1996.For consistency, the EALCO values include both land surface evapotranspiration (ET) and water surface evaporation (E) for MRB and SRB, but land surface ET only for Canada.
large uncertainties in ET estimates associated with using different approaches.In general, our ET was found close to the results obtained using surface water budget approach, which is deemed to involve less uncertainty in ET estimates than other approaches.
Modeling ET on a large scale is often challenged by a variety of factors including the availability and quality of input data, and the robustness and validity of modeling algorithms and model parameters.Information on vegetation and soils over high-latitude regions in northern Canada is still very limited.Remote sensing has largely been used to fill the data gaps.This study used the best available national-scale data for vegetation retrieved from satellite observations.One limitation of the land cover product is that plant species detail is not present.Different plant species in the same plant functional type may have significantly different hydraulic characteristics such as leaf stomatal conductance.This modeling study is unable to account for these effects on ET.
The soil texture data have much coarser spatial resolution than other land surface datasets.Validations of this dataset have been more focused on the agricultural regions, and the uncertainties may be high in non-agricultural areas.Nevertheless, this is the only Canada-wide soil dataset available.Soil texture affects ET in the model mainly through its water-holding capacity and hydraulic conductivity, which control the soil water dynamics and root water uptake.For regions where P is significantly higher than ET, the impact of soil texture on ET estimates is small.Since most of nonagricultural areas in Canada are relatively wet, the overall impact of uncertainties in the soil texture data on the nationalscale ET estimates is probably small.Meteorological forcing fields are critical in large-scale ET modeling but difficult to obtain.Both the quantity and quality of meteorological observations deteriorate in the northern cold regions, which are a major constraint in the production of high-quality atmospheric forcing fields.The meteorological data used in this study were constructed by combining a number of global observation-based datasets with the NCEP/NCAR reanalysis.This dataset provides a valuable data source for large-scale long-term hydrological modeling studies.One weakness of this data product is the low spatial resolution of 1.0 degree.The bilinear spatial interpolation used in this study further smoothed the original data representing the grid averages and may lead to more frequent P events with lower intensity.Its direct impact on ET modeling would include the overestimation of canopy interception, thereafter canopy evaporation, which, at the same time, depresses transpiration.The indirect impact would include overestimation of surface water infiltration when heavy rain storms occur.The fact that ET is mostly energy constrained over the landmass reduces the above impacts.It is important to note that the ET from this study is based on gridded forcing fields and does not reflect micro-scale variations.The site-level meteorological conditions could have substantial variations within a grid of the forcing fields.As such, the ET measured at a specific site might not be directly comparable with our mapped ET.As an example, the three Fluxnet Canada Research Network (FCRN) tower flux sites of mature boreal aspen, black spruce, and jack pine forests in the Saskatchewan province are all located within 100 km of each other and are exposed to the same mesoscale weather systems, but considerable difference in summer rainfall can happen over the three sites.Specifically, during the severe drought in 2001-2003 that occurred in mid-west Canada, the black spruce and jack pine sites received significantly more P than the aspen site.This resulted in the large difference in the impact of drought on their ecosystem fluxes (Kljun et al., 2006).
Large gaps in data and knowledge on ET over Canada's landmass exist, particularly for the north.The ET process over the Canadian landmass is highly complicated, involving important changes over time and space such as snow accumulation and melt, soil freeze and thaw, distinct phenological cycle of vegetation, and highly diversified land cover types.All these ensure the importance of robust modeling algorithms and accurate model calibration.The limitations of using simple models for ET estimation as compared to process-based models such as EALCO were analyzed and discussed in Fernandes et al. (2007).EALCO is a comprehensive model that includes all the major physical and physiological processes related to ET. Special effort has been made to develop rigorous cold region processes.However, process-based models, like any other models, may suffer from inadequate model calibrations and validations.Eddy correlation tower flux measurements have been a major data source for model calibration and validation in recent studies (e.g., Liu et al., 2003;Cleugh et al., 2007;Zhang et al., 2010;Mu et al., 2011;Vinukollu et al., 2011).Similarly, EALCO has been calibrated and validated throughout using flux measurements particularly data from Fluxnet Canada Research Network (FCRN) (see references cited in Sect.2.1).However, it is worth noting that the FCRN sites (measurements have been ceased) only represent a tiny fraction of Canada's landmass (one over a million), and these sites were mostly located in the south.The observational data and knowledge gaps on ET need to be filled to reduce further the uncertainty in large-scale model applications as in this study.

Fig. 2 .
Fig. 2. The EALCO model scheme for a heterogeneous land surface pixel.

Fig. 3 .
Fig. 3. Land cover of Canada's landmass retrieved from the SPOT 4 VEGETATION data for year 2000.

Fig. 4 .
Fig. 4. The maximum leaf area index (LAI, m 2 m −2 ) of Canada's landmass retrieved from the SPOT 4 VEGETATION data for year 2000.

Fig. 5 .
Fig. 5. Terrain elevation (meters above sea level) of Canada's landmass generated by combining three datasets of CDED, SRTM and GTOPO3 (see Sect. 2.2.4 for details).The two polygons represent the Mackenzie River basin (MRB, delineated by the black line) and the Saskatchewan River basin (SRB, delineated by the white line).

Fig. 10 .
Fig. 10.Mean annual evapotranspiration (1979-2008) modeled for the 15 ecozones in Canada.The error bars represent one standard deviation of ET among the pixels in each ecozone.

Fig. 14 .
Fig. 14.Seasonal variations of the average evapotranspiration and precipitation over the entire Canadian landmass for 1979-2008.

Table 1 .
Legends for the land cover classes in Fig.1and broad vegetation types in the EALCO model.
gions during this time period.ET in the summer months of June through August increases dramatically particularly in the southern part of the landmass.The hot spots with high monthly ET have similar spatial patterns to those of the annual ET discussed earlier.The agricultural regions in the southern Mixed Wood Plains reached higher ET than the boreal forest region to its north due to the high atmospheric evaporative demand.In contrast, the agricultural regions in

Apr May Jun Jul Aug Sep Oct Nov Dec Evapotranspiration (mm)
Fig. 13.Seasonal variations of the average evapotranspiration over each of the 15 ecozones in Canada.