Gridded climate datasets are produced in many parts of the world by applying various interpolation methods to weather observations, to which are sometimes added secondary information (in addition to geographic location) such as topography and radar or atmospheric model outputs. For a region of interest, the choice of a dataset for a given study can be a significant challenge given the lack of information on the similarities and differences that exist between datasets, or about the benefits that one dataset may present relative to another. This study aims to provide information on the spatial and temporal differences between gridded precipitation datasets and their implication for hydrological modeling. Three gridded datasets for the province of Quebec are considered: the Natural Resources Canada (NRCan) dataset, the Canadian Precipitation Analysis (CaPA) dataset, and the dataset from the Ministère du Développement Durable, de l’Environnement et de la Lutte contre les Changements Climatiques du Québec (MDDELCC). Using statistical metrics and diagrams, these precipitation datasets are compared with each other. Hydrological responses of 181 Quebec watersheds with respect to each gridded precipitation dataset are also analyzed using the hydrological model HSAMI. The results indicate strong similarities in the southern parts and disparities in the central and northern parts of the province of Quebec. Analysis of hydrological simulations indicates that the CaPA dataset offers the best results, particularly for watersheds located in the central and northern parts of the province. MDDELCC shows the best performance in watersheds located on the south shore of the St. Lawrence River and comes out as the overall second-best option.
Gridded precipitation datasets are popular in hydrology and climatology because of the low spatial and temporal density of weather stations in many parts of the world. They are also convenient for watershed modeling because they provide continuous values for input data, uniformly distributed in space and time. Precipitation is one of the key inputs to hydrological models and, along with temperature, constitutes the bare minimum in terms of inputs to run hydrological models. The various hydrological processes are often simulated by such models using from subdaily to daily time steps. Therefore, the model needs precipitation data at the same frequency. The spatial interpolation of daily precipitation grids from a limited number of ground observations presents a major challenge, in part because the spatial density of ground stations might not be high enough to fully grasp the spatial structure of the real underlying precipitation field. For instance, in kriging, the spatial structure of data is represented by a variogram. The optimal equation for this semivariogram is dependent on regional weather and geography, but also on the time step of interest (e.g., Grimes and Pardo‐Igúzquiza 2010; New et al. 2001). As stressed by Grimes and Pardo‐Igúzquiza (2010), a daily variogram for precipitation accumulation is most influenced by local weather while a monthly variogram should also account for large-scale patterns of storms and more general climate patterns. It is well established that spatial and temporal variability in precipitation affects several aspects of hydrological modeling such as the watershed response, the timing of peak runoff, the estimation of model parameters, and the hydrological model outputs, as reviewed by Ly et al. (2013).
There are many methods to obtain gridded precipitation datasets from weather stations, ranging from simple techniques [e.g., arithmetic mean, Thiessen polygons (Thiessen 1911), and inverse distance weighting (Shepard 1968)] to more sophisticated techniques like kriging (Goovaerts 2000) and thin plate splines (Hutchinson and Bischof 1983). The last two interpolation techniques can incorporate additional information such as elevation (Goovaerts 2000; Hutchinson et al. 2009; Tapsoba et al. 2005) and radar (Haberlandt 2007). Choosing the most appropriate gridded precipitation dataset for a given hydroclimatic study in a given region may prove crucial and is limited by the available information regarding the similarities and distinctions between different types of gridded datasets and the advantages of one type of dataset relative to another. Therefore, comparative analyses of different precipitation gridded datasets for regions of interest are necessary to analyze their respective benefits.
Analyses of gridded precipitation datasets for regions such as the Athabasca watershed in Canada (Eum et al. 2014), the midwestern states in the United States (Ensor and Robeson 2008), the Sierra Nevada of California in the United States (Mizukami and Smith 2012), the continental United States (Essou et al. 2016), Australia (Tozer et al. 2012), and Brazil (Silva et al. 2007) can be found in the literature. Eum et al. (2014) conducted an intercomparison of three gridded datasets: the North American Regional Reanalysis (NARR), the Canadian Precipitation Analysis (CaPA), and the thin plate smoothing splines [named for Natural Resources Canada (NRCan) in this study]. They found systematic differences between the gridded datasets and significant differences in hydrologic model simulations made with input data from the three datasets. Differences between gridded datasets were also found across Australia by Tozer et al. (2012), as well as marked differences between the gauged precipitation and the precipitation grid cell within which the corresponding gauge sits. They also reported markedly different runoff responses associated with each gridded dataset. Ensor and Robeson (2008), in their study of gridded precipitation datasets for midwestern states in the United States, found a significant increase of the frequency of low precipitation and a great reduction of the frequency of heavy precipitation in gridded datasets compared to observations from weather stations. They attributed those differences to the interpolation processes. The same discrepancies are also described by Silva et al. (2007). Mizukami and Smith (2012) analyzed the temporal inconsistencies that occur in gridded precipitation time series for the Sierra Nevada, due to various factors such as gauge relocation and changes in data processing methods and their impact on hydrologic model output. They proposed an approach to uncover and correct the inconsistencies and found that river flow simulations based on the corrected input data were more consistent and precise than those based on the inconsistent gridded data. Essou et al. (2016) compared temperature and precipitation of three gridded datasets over the continental United States in order to evaluate their differences and their impact on lumped hydrological modeling. They found significant differences between the datasets, but their performances regarding the hydrological modeling were overall similar. A summary of the studies mentioned in this paragraph along with their main findings is given in Table 1. To the author’s knowledge, a specific intercomparison of CaPA with other popular interpolated precipitation grids with a focus on hydrological applications for the easternmost portion of Canada has not yet been reported in the literature. The province of Quebec is especially appropriate to carry such a study because of the importance of hydrological modeling for this part of the world. In fact, Quebec is the most important hydropower producer in Canada, which is among the top world hydropower producers.
This study presents an intercomparison of three different gridded precipitation datasets over the province of Quebec. The precipitation gridded datasets are compared with each other and the comparison is carried out for the entire study area. Then, to go one step further, the datasets are also used for the calibration of a lumped conceptual hydrological model on 181 watersheds with different physical and hydroclimatic characteristics. This allows for further comparisons between the gridded datasets in terms of their hydrologic response for a variety of watersheds. The three gridded precipitation datasets available for this research were produced by various interpolation methods in order to provide spatialized information in Canada, even in areas where there are very few weather stations. The first two datasets were interpolated from those weather stations. Details regarding the available stations, the interpolation methods, and the spatial coverage of the datasets are given in section 2b. The third gridded precipitation dataset combines information from weather stations with the outputs of an atmospheric model.
The paper is organized as follows. The study area, the gridded datasets, and the observed data are presented and described in section 2. In section 3, the various metrics and diagrams used to conduct the intercomparison are detailed, as is the hydrological model calibration technique. The results regarding the comparisons between the different gridded precipitation datasets are presented in section 4 along with the model simulations results. A discussion and concluding remarks are given in section 5.
2. Study area and gridded and observed datasets
a. Study area
The province of Quebec is located in the eastern part of Canada and covers a surface of over 1.6 million km2 (Fig. 1). The climate is characterized by significant temperature and precipitation variations from region to region, depending on latitude, topography, and maritime influence. The province of Quebec comprises 430 major watersheds, 100 of which drain an area larger than 4000 km2 (Environment Québec 2002). Of these 430 watersheds, 181 were selected for the study presented here. The selection was based on the availability of discharge data for the period of interest (2002–10). More details are given in the methodology section (section 3).
b. Gridded and observed data
Three different gridded datasets are used for this study, with information accumulated on a daily basis. Those datasets are produced by different institutions with different interpolation methods (see Table 2 for a summary of some of the datasets characteristics and Fig. 2 for the geographical location of stations used by each gridded dataset). The 2002–10 time period, which is common to all three datasets, is considered in this study.
Ministère du Développement Durable, de l’Environnement et de la Lutte contre les Changements Climatiques du Québec (MDDELCC) grids are based on an exact interpolation method, namely, ordinary kriging, while NRCan and CaPA are based on inexact interpolation methods, respectively, splines and optimum interpolation (OI; Table 2). Exact interpolation methods ensure that the predicted values at points for which measurements are available will be exactly equal to those values. This is clearly an advantage for datasets where the measurements are very reliable and associated to low uncertainty. On the contrary, with inexact interpolation methods, predicted values at points for which measurements are available will not necessarily be equal to measured value. For precipitation, especially when the measurement network is sparse, this characteristic of inexact interpolation methods is an advantage over exact interpolation methods. For example, precipitation measurements at ground stations do not necessarily correspond to the true precipitation at this location for a number of reasons. In particular, there is the well-known phenomenon of precipitation undercatch (e.g., Pollock 2012). In addition, some inexact interpolation methods allow for the fusion of two or more sources of information regarding the variable of interest, which is the case here for CaPA. In the case of gridded precipitation products, the location expressed by geographic coordinates represents the gridcell centers. The precipitation value is given for the gridcell area defined by the product resolution and not a single location, as it is the case for the ground station measurements.
1) Gridded dataset from Natural Resources Canada
NRCan used the Australian National University Splines package (ANUSPLIN) algorithm developed by Hutchinson (1995) to apply the thin plate smoothing splines interpolation method to Canadian weather observations. Initially, Hutchinson et al. (2009) interpolated daily minimum and maximum temperature and daily precipitation from 1961 to 2003 for all of Canada on a high-resolution grid (10 km). These daily data were updated by Hopkinson et al. (2011). They extended the data up to 2011 and reduced the residuals caused by the existence of different climatological days. This updated version took into account 7514 observation stations from Environment and Climate Change Canada (ECCC) for 62 years starting from 1950. This dataset is called NRCan throughout this paper (Table 2).
2) Gridded dataset from ECCC
To further reduce the uncertainty that arises from the low density of weather stations in the northern part of Canada, ECCC developed an innovative precipitation analysis (CaPA) that combines the data available from weather stations with additional (or secondary) information from short-range precipitation forecasts produced by the Global Environmental Multiscale (GEM) model. Since 2014, CaPA also assimilates radar observations and satellite data, but this is not the case for the 2002–10 period on which this study focuses.
More precisely, the CaPA algorithm is based on OI (Table 2), which is similar to residual kriging (e.g., Baillargeon 2005). OI involves a spatial covariance function, while kriging involves a variogram. Precipitation quantities can be influenced by explanatory variables such as altitude, latitude, or others. Both OI and residual kriging can account for additional information from such explanatory variables during the interpolation. OI produces an analyzed precipitation field, which is a weighted combination of background field and observation field. For CaPA under the 2002–10 period, the background field is the precipitation forecast from GEM. At ground station locations, those stations are given high credibility (a weight of one or close to one) while at locations far away from ground stations the background (GEM) is given more credibility (high weight). The values of the weights in OI are obtained by minimizing the variance of the analysis, under the assumption that background and observations are uncorrelated and that each is unbiased.
Similarly, in residual kriging, a regression is first performed between precipitation and an external variable. The external variable may well be the background field from a model. In any case, the external variable should have a higher spatial density than the ground stations network. At points where the true value of precipitation is known, the error (or residual) of the aforementioned regression can be computed and kriged to obtain a map of corrections. The precipitation at any point in space can then be obtained by adding the correction from this map to the result of the regression between precipitation and the external variable.
Brasnett (1999) used OI at the Canadian Meteorological Centre (CMC) to generate screen-level analyses of several meteorological variables. The short-range precipitation forecasts are used as background values. The method performs the analysis on innovations, that is, the difference between an observation and the corresponding background value. The innovations are then weighted according to their error statistics to provide the precipitation estimates (Fortin and Roy 2011; Mahfouf et al. 2007). It should be noted that the version of CaPA used for this study does not include radar and satellite data as only recent CaPA products include those data (from 2014 to present).
The observation database used by CaPA consists of 6-h surface precipitation accumulations available in real time at the Meteorological Service of Canada (MSC). These archives include North American data from four networks: the surface synoptic observations (SYNOP), the aviation routine weather report (METAR), the Cooperative Observer Program (COOP), and the Réseau Météorologique Coopératif du Québec (RMCQ). The SYNOP network contains about 300 manual stations and 750 automated weather stations maintained by ECCC and partner organizations. The METAR network contains more than 1300 stations generally located at airports. The COOP is a network of more than 11 000 volunteers in the United States who collect observations in U.S. territory. The RMCQ is a cooperative network of private companies and provincial agencies in Quebec comprising about 400 stations (Haghnegahdar et al. 2014).
Finally, note that the spatial resolution of CaPA changes with time. It follows the improvement of the numerical weather prediction model (GEM) on which it is based. On 18 May 2004, the horizontal resolution of the GEM regional deterministic model went from 0.22° (approximately 24 km) to 0.1375° (approximately 15 km). Further successive refinements of the horizontal resolution were also performed after 2010. The 2002–10 period therefore includes two different horizontal resolutions. For this study, CaPA grids were retrieved through an interactive data portal that allows the user to choose from different predefined horizontal resolutions (here 0.1°). Data retrieved in this fashion are not really interpolated in the strict sense but rather downscaled using the nearest neighbor method.
3) Gridded dataset from the government of Quebec
Gridded climate data are also produced by the MDDELCC. The grid covers the province of Quebec with a spatial resolution of 0.1° (approximately 11 km; Table 2). The data are interpolated from 329 weather stations of the Programme de Surveillance du Climat du Québec (PSC) and 41 ECCC stations located north of the 49th parallel to ensure spatial coverage of the northern part of the province of Quebec. The interpolation is performed by ordinary kriging using a daily spherical variogram. The maximal distance of separation between points for the construction of this variogram is limited to 200 km because it is inferior to the synoptic scale and provides the best results in terms of agreement with ground station measurements (Bergeron 2015). The data cover the period from 1961 to the present. This dataset is hereafter referred to as MDDELCC.
c. Observed streamflow data
a. Regridding for standardization of the three gridded datasets
Several evaluation metrics require a gridcell-by-gridcell comparison of datasets. Therefore, a regridding of NRCan to the same grid as MDDELCC was necessary to allow for comparisons. The regridding of CaPA was performed during data retrieval as mentioned in section 2b(2). The following methodology was used for regridding of NRCan: either one or two NRCan grid points (cell centers) fit into one CaPA/MDDELCC grid cell. When only one NRCan cell fit into the CaPA/MDDELCC cell area, the value of the NRCan cell was simply assigned to the standardized NRCan grid cell. When there were two NRCan grid points within a CaPA/MDDELCC gridcell area, the two values were averaged to form the standardized grid cell.
b. General intercomparison of gridded datasets
This intercomparison of the gridded precipitation datasets is conducted using daily accumulated precipitation. Special care was taken regarding the time of validity of each gridded dataset. The NRCan grids are based on daily station data from ECCC, which are valid either at 0600 UTC (automatic and synoptic stations) or at 1300 UTC (for manned climatic stations). All stations are pooled together before interpolation, and it is considered that the grid is valid at 1300 UTC. The 6-h CaPA product that was used in this study was cumulated on 24-h periods ending at 1200 UTC. Similarly, for the MDDELCC grids, hourly precipitation data are collected from ground stations and cumulated. The MDDELCC grids are valid at 1300 UTC (see Table 2).
The first comparison between grids involves the analysis of annual, seasonal [winter (DJF), spring (MAM), summer (JJA), and autumn (SON)], and monthly metrics computed with, or between, the different datasets. The results are either presented as maps for the entire provincial territory or as curves as a function of time. The following metrics are used:
Precipitation annual cycle: annual cycle of precipitation provides relevant information about the distribution of the precipitation throughout the year (variation, peaks, etc.). The mean total amounts of precipitation are first computed for each month of the year in each grid cell (12 values per grid cell) for the 2002–10 period. The monthly values are then spatially averaged over the entire province, for each dataset (one curve per dataset).
Relative bias: the gridcell-by-gridcell cumulative relative differences (RDs) are computed using the daily precipitation time series of two datasets. Seasonal time series are considered (one map per season for every combination of two datasets).
Ratio of variances: the gridcell-by-gridcell ratio of variances (RV) is computed between the daily precipitation time series of two datasets. Seasonal time series are considered (one map per season for every combination of two datasets).
c. Hydrological modeling: Input data and model calibration
When comparing different gridded precipitation datasets, it might be difficult to obtain data from ground stations that are completely independent from the gridded products. In fact, when designing a new gridded precipitation product, one typically wishes to use all the available information so very few ground stations, if any at all, are left out. To avoid this problem and allow for an independent evaluation of each gridded dataset for hydrological purposes, here watersheds are considered as large independent “gauges” capable of measuring a transformation of precipitation (e.g., Oudin et al. 2006). Consequently, a hydrological model can be fed by gridded precipitation dataset A or by gridded dataset B. Corresponding simulated streamflow series can each be compared to observed series. If series simulated from A match the observed streamflow series better than simulated series from B, then gridded dataset A is of superior quality compared to gridded dataset B. Of course, this framework for comparison rests on a series of assumptions: the hydrological model must be an adequate representation of the watershed, groundwater and surface water exchanges are almost time invariant or very well represented by the model, and all measurements (for instance streamflow) are free of error.
Here the lumped conceptual HSAMI hydrological model (Fortin 2000) was chosen to simulate daily discharges, using daily precipitation from each of the three gridded datasets in turn. HSAMI is a reservoir-based model that simulates the main hydrological processes: snowmelt, evapotranspiration, infiltration, and runoff, using simple water balance principles and empirical equations for snow processes and potential evapotranspiration. HSAMI is utilized to simulate natural inflows on 181 watersheds with surface areas ranging from 301 to 69 191 km2. As mentioned above, the streamflows are valid at 0000 LST and the precipitation from all gridded datasets are valid at 0700 or 0800 LST (1200 or 1300 UTC). To avoid the noncausal response in small watersheds where the streamflow would respond before the rain starts, the watersheds with an area less than 300 km2 were removed.
Since it is a lumped model, precipitation and temperature inputs must be spatially averaged over each watershed before being fed into the model. The Thiessen’s polygons method (Thiessen 1911) was used to calculate the average data (precipitation and temperature) on watersheds for each dataset. It should be recalled that the aim is to evaluate the similarities and discrepancies in the streamflow simulations from different precipitation inputs. The CaPA dataset only provides precipitation, while MDDELCC and NRCan provide both precipitation and temperature. The temperature from the NRCan datasets was used with all three precipitation datasets, as NRCan is the one most widely used in the literature (e.g., Chen et al. 2011; Shrestha et al. 2012; Troin et al. 2015).
HSAMI has been used in an operational framework for short-term flow forecasting by Hydro-Québec for over 30 years. It has also been used in climate change impact studies (Chen 2012; Minville et al. 2009; Poulin et al. 2011). HSAMI comprises a total of 23 free parameters that need to be calibrated. Model calibration is conducted automatically using the covariance matrix adaptation evolution strategy (CMA-ES) algorithm (Hansen and Ostermeier 1996, 2001). Arsenault et al. (2014) tested 10 stochastic optimization methods and found that the CMA-ES algorithm was an optimal choice for calibrating the HSAMI model. The Nash–Sutcliffe efficiency (NSE; Nash and Sutcliffe 1970) is used to evaluate the performance of the model simulations during calibration. First the even years are used for calibration while the odd years are used for the validation. The calibration is then performed using the odd years while validation is performed on the even years. For each precipitation dataset and for each watershed, five parameter sets are obtained using the even/odd approach and five parameter sets are obtained using the odd/even approach, for a total of 10 different parameter sets. The overall best parameter set according to NSE in the validation period is selected.
The following metrics are used to compare simulated and observed streamflows:
differences in NSE criterion between series of simulated streamflow based on different gridded precipitation datasets for the validation period;
maps showing performance in terms of NSE for watersheds and for each gridded precipitation dataset;
relative bias (RB), correlation coefficient (CC), and RV between simulated and observed streamflow time series; and
mean absolute error between simulated and observed streamflows below the 10th (MAE10) and above the 90th (MAE90) percentile of the observed flows for each precipitation gridded datasets.
a. General intercomparison of gridded datasets
To examine the differences between the three gridded precipitation datasets, the results from the metrics introduced in section 3b are first presented. Figure 3 shows the spatially averaged annual precipitation cycles (mm month−1) from the three datasets. The difference between NRCan and the two other datasets is visible. However, here CaPA and NRCan appear to be more synchronous with each other than with MDDELCC, mostly in the winter months (January and February) and in the months of May–October.
Figure 4 provides an overview of the spatial distribution of the 9-yr mean seasonal relative biases (%) between the three precipitation datasets. Generally, CaPA and MDDELCC present a wet bias with respect to NRCan, mostly in the center and north of the province for all seasons (Fig. 4, top and middle), except for CaPA in the wintertime showing a dry bias with respect to the NRCan dataset. This dry bias can also be seen by the more pronounced wet wintertime bias of MDDELCC with respect to CaPA (Fig. 4, bottom). Wet biases of more than 50% are obtained for CaPA and MDDELCC compared to NRCan in the central part of Quebec, more precisely, over the Saguenay–Lac-Saint-Jean watersheds for all seasons (dark blue spot in Fig. 4, top and middle). NRCan strongly underestimates the precipitation in this region compared to the other datasets (this particular aspect is addressed in the discussion section). In contrast, precipitation amounts over parts of southern Quebec are quite similar for all datasets in almost all seasons, the only difference occurring in the summertime with a dry bias for CaPA compared to NRCan, and a wet summertime bias for MDDELCC compared to CaPA. Given that there is a high density of weather stations in southern Quebec, it is expected that the biases here would be the lowest, whatever the season. A possible reason for the different behaviors found in the summertime is that summer rainfall exhibits particularly high spatiotemporal variability and the interpolation methods used to produce each dataset capture these structures differently. In addition, the weather stations used are not always the same for the three gridded datasets (as exposed in section 2b).
For every season and for each pair of gridded datasets, the ratio of variances of time series was also calculated. The results are shown in Fig. 5. Globally, in the southern part of Quebec, the ratio of variances is around unity. This is again an expected result, given the high spatial coverage of this part of Quebec by weather stations. Moving northward, differences arise in the temporal variability of precipitation. For CaPA and NRCan, the ratio of variances greatly exceeds unity (above 1.5) in much of northern Quebec, except in the winter season where the ratio of variances is below unity in the far north of Quebec. Thus, it appears that NRCan underestimates the temporal variability of precipitations compared to CaPA, except in winter when it seems to slightly overestimate the variability. Meanwhile, the variability is quite similar between MDDELCC and NRCan, again except around the Saguenay–Lac-Saint-Jean region, where NRCan strongly underestimates it. MDDELCC also seems to underestimate the variability compared to CaPA in the center and northern parts of Quebec.
Whether for the map of relative biases or the map of ratio of variances, the most considerable discrepancies between the three datasets are found in central and northern Quebec and over the Saguenay–Lac-Saint-Jean region.
b. Hydrological modeling
1) General intercomparison of simulated flows
Figure 6 shows the distributions of the differences in NSE scores for the 181 watersheds between each pair of precipitation gridded datasets used as input to the HSAMI model. The performances are generally comparable as the median and the first and third quartiles are quite similar. The differences in performance are mostly associated with the outliers. There are many positive outliers when NRCan is involved in the differences of NSE, revealing a decrease in performance when NRCan grids are used as input for HSAMI compared to CaPA and MDDELCC. The differences of NSE between CaPA and MDDELCC are small as the values lie between −0.2 and 0.2.
To better visualize the differences from the watersheds’ perspectives, performance maps were made with a color code displaying the ranges of the NSE values, shown in Fig. 7. The watersheds that perform well are common to all three datasets (e.g., the large watersheds of central and northern Quebec with NSE values greater than 0.8). The poorer performance noted in Fig. 7 for NRCan is mainly attributable to watersheds located over the Saguenay–Lac-Saint-Jean region (red area in Fig. 7, left). Here, the NSE is less than 0.5 for most of the watersheds when NRCan precipitation is used as input to the HSAMI model, which is not the case for CaPA and MDDELCC [Fig. 7, center and right (respectively)]; instead, they have NSE values that exceed 0.7 for the same watersheds. This could be related to the large differences that were observed in section 4a between the NRCan grids and the two others for the relative bias and temporal variability.
Figure 8 identifies the gridded datasets that provide the best performance in terms of NSE for each watershed. Watersheds shown in red are those for which NRCan precipitation grids lead to the best performance in terms of simulated–observed streamflow comparison. Similarly, watersheds in blue are those for which CaPA grids offer superior performances and watersheds in green are those for which MDDELCC grids are superior. Note that for some watersheds, the difference between one gridded dataset or another in terms of NSE is very small. This is the case, for instance, for the Caniapiscau and À la Baleine watersheds, for which high NSE values are achieved using any of the gridded datasets (see the two dark blue watersheds in Fig. 7, at 56°N). However, the dominance of CaPA is not statistically significant for the 181 watersheds. Indeed, the Kruskal–Wallis test was conducted on three groups of NSE. The p value obtained is 0.2458, indicating that the Kruskal–Wallis test does not reject the null hypothesis that all three groups of NSE come from the same distribution at a 5% significance level.
The predominance of CaPA is obvious, particularly for large watersheds and those located in the central and northern part of Quebec. NRCan shows the best performance for 44 watersheds representing 11% of the total watersheds area, CaPA for 70 watersheds representing 64% of the total watersheds area, and MDDELCC for 67 watersheds that represent 25% of the total watersheds area. NRCan shows the worst performance for 69 watersheds representing 40% of the total watersheds area, CaPA for 65 watersheds representing 24% of the total watersheds area, and MDDELCC for 47 watersheds representing 35% of the total watersheds area. CaPA shows the worst performance particularly for small watersheds located in the southern part of Quebec, while MDDELCC shows the best performance in the watersheds located on the south shore of the St. Lawrence River.
Figure 9 links geographical attributes of the watersheds with modeling performances. The distribution of area and mean longitudes and latitudes of the watersheds that show the best performance in terms of NSE for each gridded dataset is shown. From this figure it can be concluded that the watersheds for which NRCan outperforms the other gridded datasets are mainly small-sized watersheds located in the southern part of the province of Quebec. For CaPA, the best performances are obtained with large watersheds located in the central and northern parts of the province. Watersheds that show the best performance with MDDELCC have various sizes with a majority of small watersheds, mainly located in the southeastern part of the province, as also shown in Fig. 8. For instance, on the south shore of the St. Lawrence River, MDDELCC shows best performance in 75% of the watersheds in this area.
2) Analysis of performance improvements
All watersheds for which the NSE is below 0.6 were inventoried and reported in Tables 3 and 4, including their areas and mean annual precipitation as well as their relative bias. The watersheds are identified by the number of their hydrometric station. NSE values of less than 0.6 obtained with CaPA or MDDELCC are also systematically less than 0.6 with NRCan. Table 3 presents the watersheds that show poor performance with all three gridded datasets.
Table 4 presents the watersheds for which NSE values obtained with NRCan improved significantly with CaPA and MDDELCC. For these watersheds, the values of NSE increased from 0.40 with NRCan, to 0.81 with CaPA, and to 0.76 with MDDELCC, on average. Moreover, it can be noted that for these watersheds, NRCan shows a dry bias of at least 30% of the total precipitation compared to CaPA and MDDELCC. The poor performance of these watersheds with NRCan could be linked to NRCan’s severe underestimation of precipitation compared to other gridded data. In effect, the results in section 4a not only report the severe precipitation underestimation by NRCan, but also the temporal variability relative to CaPA and MDDELCC in the areas where most of these poor performance watersheds are found. For some watersheds presented in Tables 3 and 4, it would be necessary to question the quality of the river flow records, as the poor performance could be caused by the use of erroneous or low-quality streamflow data for calibration. For instance, this is the case for the Manouane River (gauge number 062209 in Table 4). The streamflow records for this river are known to be of poor quality because they are reconstituted from a water balance equation and subject to random fluctuations attributable to wind on the reservoir where the level is measured.
3) Bias, correlation, and relative variance of streamflows
Figure 10 shows the scatterplots that indicate the correspondence between indicators for different simulated streamflows such as RB, CC, and RV. These indicators are calculated between simulated and observed streamflows for the validation period. Regarding the RB, the scatterplots show a tendency of the values to be shifted on the side of NRCan when comparing NRCan with CaPA and MDDELCC, and on the side of MDDELCC when comparing MDDELCC with CaPA (Fig. 10, top). This implies that for the majority of watersheds, the bias for simulated streamflows is more important for NRCan than for CaPA and MDDELCC. The RB is also more important for MDDELCC than for CaPA. Regarding the CC (Fig. 10, middle), scatterplots show a tendency of the values to be shifted on the sides of CaPA and MDDELCC. This indicates that the streamflows simulated from CaPA and MDDELCC precipitation are more correlated to the observed streamflows than the ones simulated from NRCan for most of the watersheds. For the RV (Fig. 10, bottom), the scatterplots show a tendency of the values to be dispersed symmetrically around the diagonal. This means that RVs are not systematically influenced by given precipitation input datasets.
In short, this analysis reveals that the simulated streamflows from CaPA and MDDELCC input datasets show the lowest bias and are more correlated to the observed streamflows than input precipitation from NRCan. CaPA shows the lowest bias for simulated streamflows.
4) Comparison of modeling performance for extreme low and high flows
The mean absolute simulation errors in hydrological simulation based on each gridded dataset are calculated for streamflow values below the 10th percentile (i.e., MAE10) and above the 90th percentile (i.e., MAE90) of the observed streamflow series (all watersheds). The results are presented in the scatterplots of Fig. 11. For the MAE10, all values are close to the diagonal. This implies that all gridded precipitation datasets are comparable in terms of quality for dry periods. For the MAE90, the scatterplots show a slight tendency of the values to be shifted on the side of NRCan when comparing NRCan with CaPA. They are also shifted on the side of MDDELCC when comparing MDDELCC both with CaPA and with NRCan. This implies that CaPA always leads the lowest (best) MAE90 values. Since adequate simulation of high streamflow values by hydrological models are of great interest for many practical applications (flood forecasting, for instance), this is clearly an important advantage of CaPA over the two other gridded datasets. Probably because of this combination of model background with ground station observations, CaPA is better able to capture all precipitation events over entire watersheds and allow for a more accurate representation of high flows by HSAMI.
5. Discussion and conclusions
High-resolution gridded precipitation datasets have been developed by agencies throughout the world to facilitate hydrological modeling in areas where ground-based observations are scarce. These grids can be obtained from the interpolation of available ground station observations using a variety of methods (kriging, splines, and others). It is also possible to combine ground observations with alternative information, such as radar precipitation measurements or atmospheric model outputs to refine the description of the precipitation field.
Still, the choice of the most appropriate gridded dataset for a specific application can be a significant challenge, especially in areas where the spatial and temporal density of observations is low, such as in remote areas like the northern parts of Canada. This study focuses on analyzing the quality of three gridded precipitation datasets, NRCan, CaPA, and MDDELCC, for hydrological modeling purposes. Two of those datasets are obtained solely from spatial interpolation using weather stations, while the third combines information from weather stations and the GEM atmospheric model. The study area covers the entire province of Quebec. Gridded datasets were first compared with one another.
Spatial and temporal variabilities show significant differences in the northern part of Quebec. Regarding the differences in spatial variability among datasets, the results show that NRCan strongly underestimates precipitation in the central and northern parts of Quebec compared to CaPA and MDDELCC. Important differences between the three gridded datasets in the northern part of Quebec were exposed (Figs. 4, 5). More specifically, the results show that the NRCan grids are of poor quality over the Saguenay–Lac-Saint-Jean area (blue area in the center of the province; Fig. 4, middle) compared to the two other gridded datasets. NRCan grids are obtained using only ECCC weather stations, MDDELCC uses a few more stations, and CaPA incorporates a larger number of stations including stations in the Saguenay–Lac-Saint-Jean region (see section 3a).
It was also found that NRCan and MDDELCC grids underestimate the temporal and spatial variability of precipitations relative to the CaPA dataset for all seasons but winter. The difference in terms of temporal variability between CaPA and the other two datasets (Fig. 5) is especially striking. CaPA displays high variance compared to NRCan while the MDDELCC has a lower variance compared to CaPA during winter and autumn months along the eastern coast. It was already concluded from Fig. 3 that the precipitations computed from CaPA and MDDELCC are probably the better estimates of the truth.
During a second phase of the study, the three gridded precipitation datasets were used to calibrate and validate the hydrological model HSAMI on 181 watersheds displaying various physical characteristics and hydroclimatic conditions. This step made it possible to evaluate and compare the model’s response using each gridded dataset and allow for an objective comparison of interpolated precipitation datasets when no independent ground stations are available.
The assumption regarding a more detailed representation of the precipitation field in CaPA than in the other two gridded datasets (as mentioned above) was verified indirectly though hydrological modeling. Results show that the use of NRCan grids as precipitation inputs for the model leads to lower performances in terms of the NSE and for bias. This is most likely due to the underestimation of precipitation. In addition, Figs. 6 and 8 show that the use of CaPA leads to the best performances for most watersheds, mainly located in the center and north of the province of Quebec. Some important advantages were also found for CaPA, as the lowest bias with respect to observations was obtained for the simulation of streamflows in general and for the simulation of high flows in particular using CaPA as the precipitation field.
Ultimately, the results of this study illustrate how the combination of atmospheric model outputs and observations from weather stations can benefit hydrological modeling. The resulting grids are indeed able to better capture the spatial and temporal variability of the precipitation field, especially for areas where weather stations are scarce. Although CaPA is a Canadian product, the concept of combining ground-based observations with model outputs is quite generalizable. Similar datasets can be produced for almost any area around the globe. It is therefore expected that similar results could be obtained for other areas of the world with low weather station density.
In light of those results, if one was to select a single precipitation dataset, the general recommendation out of the analysis performed in this study is that gridded precipitation from CaPA should be prioritized for hydrological studies over the province of Quebec. MDDELCC comes out as the second-best overall option, but the MDDELCC product does outperform CaPA on the south shore of the St. Lawrence River. This is not surprising, since the station density is clearly higher in this region in the MDDELCC product than in the CaPA product (see Fig. 2). Hence, in this region it is recommended to use the MCCELCC product. However, CaPA is still an evolving product. Since November 2014, CaPA has included the assimilation of radar data in its interpolation scheme. According to Fortin et al. (2015), this led to further improvements in the accuracy of the resulting precipitation grids in area covered by Doppler radar, such as most of the south shore of the St. Lawrence River. This continuous improvement of CaPA is both a blessing and a curse for hydrology, since the modifications in the assimilation scheme are not retrospectively applied to previous years: consistent datasets based on the exact same sources and interpolation scheme can only cover relatively short timespans. For instance, at the time of writing this paper, there are less than 2 years of available data that incorporate radar information. This short period limits further verifications for hydrological purposes, at least for the moment, but it does not invalidate the above recommendation.
The density of the weather station networks used by NRCan and MDDELCC are quite similar for the time period under study (2002–10). Both precipitation datasets show strong similarities in the southern part of Quebec, and differences appear gradually as the density of weather stations decreases, that is to say, when moving from south to the north. Several hypotheses could justify the differences between these two precipitation databases in central and northern parts of the province of Quebec.
It should be recalled that NRCan employs thin plate smoothing splines (ANUSPLIN) while MDDELCC uses ordinary kriging with daily variogram for interpolation. Whatever the interpolation method used, the estimation of precipitation field is very unreliable when measuring stations’ network density is low (Mizukami and Smith 2012). Differences may well occur in the processing of information available in different interpolation schemes. While NRCan and MDDELCC draw from similar station networks, NRCan incorporates elevation in its interpolation scheme, which is not the case for MDDELCC, and this can lead to different results. In particular, differences are reported in this study for the Saguenay–Lac-St-Jean region. It could be hypothesized that the strong variation in the topography between the plain of Lac-Saint-Jean and the surrounding regions are the cause. Indeed, the plain of Lac-Saint-Jean whose altitude is around 200 m above sea level is surrounded by highlands and massifs that can reach altitudes greater than 1000 m above sea level.
The application of gridded datasets for hydrological modeling has been conducted in this study using only a lumped model. For large watersheds, it would be interesting to use a distributed model to further assess the influence of the three gridded precipitation datasets in the spatial representation of the hydrological cycle. In addition, following the growing body of literature advocating for ensemble rather than deterministic hydrological modeling frameworks (e.g., Krzysztofowicz 2001), multiple gridded precipitation datasets could be used jointly as successive inputs to one (or more) hydrological model(s). Each precipitation product contains useful information, and combining several datasets can lead to improved streamflow simulations, as shown, for instance, by Arsenault et al. (2016). Further studies should concentrate on these issues.
We would like to acknowledge the Data Access Integration (DAI) team for providing the data and technical support regarding CaPA. The DAI Portal is maintained through collaboration between the Canadian Centre for Climate Modelling and Analysis (CCCma), the Research Division in Numerical Weather Prediction (RDN), and the Adaptation and Climate Monitoring unit of the Meteorological Service of Canada, as part of the Atmospheric Sciences and Environmental Issues unit. We also thank the Programme de Surveillance du Climat of the MDDELCC for graciously providing the data and guidance. Funding for this research was provided by NSERC’s Collaborative Research and Development Grants program in conjunction with Hydro-Québec, Ontario Power Generation, and Rio Tinto (Grant CRDPJ435692-12). The authors thank three anonymous reviewers for their insightful comments which helped shape the paper in its current form.