Analysis of the Joint Link between Extreme Temperatures, Precipitation and Climate Indices in Winter in the Three Hydroclimate Regions of Southern Quebec

Abstract

We analyze the relationship between four climate variables (maximum and minimum extreme temperatures, rainfall and snowfall) measured in winter (December to March) at 17 stations from 1950 to 2000 in the three hydroclimate regions of southern Quebec, and six seasonal climate indices using canonical correlation analysis (CCA) and the copula method. This analysis yielded these major results: (1) extreme temperatures are not correlated with the amount of winter rain or snow in southern Quebec; (2) winter seasonal climate indices show better correlations with climate variables than do fall climate indices; (3) winter extreme temperatures are best correlated (positive correlation) with the Atlantic Multidecadal Oscillation (AMO) in the eastern region, but show a negative correlation with the Arctic Oscillation (AO) in the southwestern region; (4) the total amount of winter snow is best correlated (negative correlation) with the Pacific Decadal Oscillation (PDO) in the three hydroclimate regions; (5) the total amount of winter rain is best (negatively) correlated with PDO in the eastern region, but shows a positive correlation with AO in the southeast region. Finally, the copula method revealed very little change in the dependence between climate indices and climate variables in the three hydroclimate regions.

1. Introduction

A number of studies have analyzed the relationship between climate indices and winter climate variables (temperature and precipitation) in Canada and Quebec [1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21]. These studies raise five types of questions which have not been completely resolved:

• All these studies analyze temperature and precipitation (rain and snow) separately. However, it is not possible, using this approach, to determine whether the two climate variables are correlated with one another, on the one hand, and if they are correlated with the same climate indices, on the other hand. From a climate standpoint, this is important. In addition, if the two climate variables are correlated with the same climate indices, it is important to determine if this correlation changed over time dependently or independently of that of one or both climate variables. On that point, Brown [10] mentioned that, in Quebec, “evidence was founded for a shift in circulation around 1980 associated with an abrupt increase in sea level pressure (SLP) and decrease in winter precipitation, snow depth and annual maximum snow water equivalent over much of southern Québec”. This type of change could potentially change the link between climate variables and climate indices over time.
• There is some debate over which climate index is the main influence on the temporal variability of winter temperature and precipitation (snow and rain) in Quebec. Some authors argue that the El Niño/Southern Oscillation (ENSO) is the main cause of this variability (e.g., [12,13,14]), while others believe it is the North Atlantic Oscillation (NAO)/Arctic Oscillation (AO) (e.g., [6,10]), the Pacific Decadal Oscillation (PDO) (Brown, 2010) or the Pacific/North American Oscillation (PNA) [21].
• There is also some debate about the choice of seasonal climate index type that best correlates with winter temperature and precipitation. While some authors use fall climate indices, others prefer winter climate indices.
• While most of these studies have highlighted changes in stationarity (long-term trend, change in mean values) of temperature and precipitation series (e.g., [22,23,24]), none has focused on the impacts of these changes in stationarity on the evolution of the dependence between climate indices and climate variables, an aspect of prime importance in the current climate change context.

In light of the foregoing, the present study has two goals: (1) to jointly correlate winter temperature and precipitation to climate indices using canonical correlation analysis to determine whether the two types of climate variables are correlated with the same climate indices; and (2) to analyze changes in this correlation over time using the copula method. These two aspects have never been analyzed in Quebec.

2. Methods

2.1. Choice and Location of Stations

This study was restricted to the period from 1950 to 2000 because very few measurements of climate variables at a large number of stations are available prior to 1950, and after 2000, climate variable measurements stopped at many stations as a result of budget cuts following program revisions by the Canadian government during the 1990s [25]. We selected every station for which temperature and precipitation measurements were taken for at least 40 years during the study period (1950–2000). Seventeen stations spread out nearly evenly throughout the southern part of the province of Quebec (Figure 1 and

Table 1. Climatic stations.

Station No.	Station	ID	Altitude (m)	Latitude (N)	Longitude (W)	Years
Southeast Hydroclimatic Region
SE1	Coaticook	7021840	259	45°09′	71°48′	51
SE2	Granby	7022802	168	45°23′	72°42′	47
SE3	Magog	7024440	274	45°16′	72°07′	51
SE4	Montréal Trudeau	7025250	36	45°28′	73°45′	51
SE5	Philipsburg	7026040	53	45°02′	73°05′	51
SE6	St Malo d′Auckland	7027520	564	45°12′	71°30′	51
East Hydroclimatic Region
E7	Bagotville	7060400	159	48°20′	71°00′	51
E8	Mont-Joli A	7055120	48	48°36′	68°13′	51
E9	Natashquan A	7045410	575	48°57′	65°31′	51
E10	Sept Îles A	7047910	310	50°13′	66°16′	51
E11	Ste Rose du Degelis	7057720	151	47°34′	68°38′	50
Southwest Hydroclimatic Region
SW12	Trois-Pistoles	7058560	58	48°09′	69°07′	49
SW13	Chelsea	7031660	112	45°31′	75°47′	51
SW14	Nominingue	7035520	305	46°23′	75°03′	51
SW15	Shawville	7038040	168	45°37′	76°28′	40
SW16	St Alban	7016800	76	46°43′	72°05′	51
SW17	St-Michel-Des-Saints	7077570	350	46°41′	73°55′	43

) were selected. These stations were subdivided into three hydroclimate regions defined for southern Quebec using principle component analysis of hydroclimate variables (e.g., [26]). These are the southeast and eastern region, primarily located on the south shore of the St. Lawrence River, and the southwestern region, located on the north shore. Temperature and precipitation (rain and snow) data were taken from [27]. These data are from Environment Canada’s Canadian Daily Rehabilitated Database [28,29,30].

Six climate indices were selected that have been shown to influence the spatial and temporal variability of temperature and precipitation in Quebec in particular, and North America in general (e.g., [5,6,10,12]. These are the Atlantic Multidecadal Oscillation (AMO), the Arctic Oscillation (AO), the North Atlantic Oscillation (NAO), the Pacific Decadal Oscillation (PDO), the Pacific/North American Oscillation (PNA) and the Southern Oscillation (SOI). Data for these indices were taken from [31].

2.2. Statistical Analysis

2.2.1. Variable (Temperature and Precipitation) and Climate Index Series

Two temperature series were assembled, namely, daily extreme maximum and minimum temperatures. These two series consist respectively of mean values of daily maximum (Tmax) and minimum (Tmin) temperatures measured from December to March during the period from 1950 to 2000. Two precipitation series were also assembled: a rainfall series and a snowfall series. Both series are comprised of the total daily sum, respectively, of rainfall (TRN) and snowfall (TSN) measured from December to March during the period from 1950 to 2000. It should be mentioned that, in southern Quebec, depending on the region, winter starts in December and ends in March.

Fours series were assembled for each climate index, including two for fall and two for winter. For fall, the first series consists of monthly means for the months of October to December (OND), and the second series, of monthly means for the months of September to December (SOND). For winter, the first series consists of monthly means for the months of January to March (JFM), and the second series, of monthly means for the months of December to March (DJFM). We used the two types of climate index series in order to compare our results with those of previous works which mostly used quarterly means, on the one hand, and because of the fact that the precipitation totals and mean values of temperature that we analyzed were calculated over four-month periods, on the other hand.

2.2.2. Canonical Analysis

Canonical correlation analysis (CCA) is widely used in climatology and hydrology to explore the links between two sets of variables (e.g., [32]). This method, which is extensively described in the literature (e.g., [33]), is generally used for large data matrices (>100). In our study, the maximum number of years for which climate variables have been measured is 51, which statistically does not lend itself to canonical multivariate analysis. To get around this technical problem, we grouped all stations from a given hydroclimate region into a single matrix. Thus, for each hydroclimate region, canonical analysis was applied to a matrix for which the number of rows is the product of the number of stations in the region and the number of years of climate variable measurement at each station, and the number of columns is the sum of the five climate variables and six climate indices (11 columns in total). The main advantage of this approach is that it allows simultaneous correlation of multiple climate indices with multiple climate variables measured at different stations. In addition, the analysis of stations with differing numbers of years of observations is also possible. This is what makes this study novel. As a result, there is no loss of information that would otherwise result from grouping stations in homogeneous climate regions, a commonly used approach in climatology to analyze the temporal variability of climate variables. In addition, with the temporal approach, it is possible to include stations for which data for some years are missing, thus optimizing the use of available measurement data for the different stations. Finally, unlike the regionalization approach, the temporal approach can be used to analyze simultaneously multiple climate variables. It should be recalled that canonical analysis is used to relate a group of dependent variables Y (in the present case, five climate variables) and a group of independent variables X (in the present case, six climate indices). These two groups will be replaced, respectively, by new canonical variables V and W. Correlations will be derived between the new canonical variables on one hand, and the original variables and new canonical variables of each group, on the other hand. Data analysis was carried out using the Statistical Analysis System (SAS, version 9) software.

2.2.3. Copula Analysis

The copula analysis is used to constrain the evolution of the dependence between two variables. Some variables may show changes in mean or variance over time (change in stationarity). This type of change may have a significant effect on the dependence between two variables over time. However, studies in the literature based on the correlation between climate variables and climate indices never take into account the possible impact of this type of change on the dependence between variables. In the present study, the copula method was used to analyze potential changes that may have occurred over time in the dependence between climate indices and climate variables. The use of this method in climatology remains very limited [22].

The dependence in a random vector (X, Y) is contained in its corresponding copula function C. Specifically, the celebrated theorem of Sklar ensures that there exists a unique $C:{\left[0,1\right]}^2\stackrel{}{\to }\left[0,1\right]$ such that

$$P\left(X\le x,Y\le y\right)=C\left\{P\left(X\le x\right),P\left(Y\le y\right)\right\}.$$

(1)

Quessy et al. [34] developed a testing procedure to identify a change in the copula (i.e., dependence structure) of a bivariate series $\left(X_1,Y_1\right),...,\left(X_n,Y_n\right)$. In our study, X represents one of the six climate indices, and Y, one of the four climate variables significantly correlated with climate index X. The idea is based on Kendall’s tau, which is a nonparametric measure of dependence. Let ${\widehat{T}}_{1:T}$ be Kendall’s tau measured for the first T observations and ${\widehat{T}}_{T+1:n}$ be Kendall’s tau for the remaining $n-T$ observations. The proposed test statistic is

$$M_n={\max }_{1<T<n}\frac{T\left(n-t\right)}{n\sqrt{n}}\left|{\widehat{T}}_{1:T}-{\widehat{T}}_{T+1:n}\right|$$

(2)

i.e., a maximum weighted difference between the Kendall’s tau. Since $M_n$ depends on the unknown distribution of the observations, the so-called multiplier re-sampling method is used for the computation of p-values. Specifically, for $n$ sufficiently large (n > 50), this method yields independent copies $M_n^{(1)},...M_n^{(N)}$ of $M_n$. Then, a valid p-value for the test is given by the proportion of $M_n^{(i)}$’s larger than $M_n$. For more details, see Quessy et al. [34]. Usually, one can expect that the series X₁,…, X_n and Y₁,…, Y_n are subject to changes in the mean and/or variance following, e.g., the smooth-change model [35]. If such changes are detected, the series must be stabilized in order to have (approximately) constant means and variances. Finally, a change in the degree of dependence between two series is statistically significant when $M_n$ > Vc, where Vc is the critical value derived from observational data. This method was applied to constrain changes in dependence between those climate indices best correlated with temperature and winter precipitation at each station in the three hydroclimate regions. Data analysis was carried out using a statistical program developed in Matlab (version R2014) by the Université du Québec à Trois-Rivières statistics laboratory.

3. Results

3.1. Relationship between Variables and Indices Climatic

Canonical coefficient values derived from data matrices using fall and winter climate indices are compared in

Table 2. Comparison of canonical correlation coefficients (r) as a function of seasonal climate indices.

Canonical Axes	Winter Climate Indices (DJFM)			Fall Climate Indices (SOND)
	r	F	P > F	r	F	P > F
Eastern Region
CC1	0.524	5.92	<0.0001	0.485	4.83	<0.0001
CC2	0.324	2.66	0.0006	0.261	2.13	0.0074
CC3	0.134	0.80	0.5985	0.176	1.33	0.2258
CC4	0.060	0.36	0.7830	0.065	0.42	0.7416
Southeast Region
CC1	0.515	7.03	<0.0001	0.445	3.67	<0.0001
CC2	0.389	4.62	<0.0001	0.183	1.35	0.1662
CC3	0.197	2.35	0.0171	0.144	1.26	0.2601
CC4	0.152	2.34	0.0739	0.114	1.31	0.2723
Southwest Region
CC1	0.539	6.25	<0.0001	0.574	4.71	<0.0001
CC2	0.399	4.21	<0.0001	0.196	1.04	0.4076
CC3	0.235	2.68	0.0069	0.136	0.83	0.5754
CC4	0.187	2.75	0.0433	0.102	0.80	0.4929

r = correlation between canonical axes (CC); statistically significant r at 5% probability level are shown in bold; DJFM = December to March; SOND = September to December.

. Results for quarterly fall indices (SON) are similar to those for seasonal fall indices (SOND, calculated over four months) presented in Table 2, and the same is true for winter climate indices. For this reason, only results for the SOND and DJFM indices are presented. It can be seen from Table 2 that canonical coefficient of correlation values (r) calculated with winter indices are generally higher than those calculated with fall indices in the three hydroclimate regions, implying that winter indices are better correlated with the four climate variables than fall indices. In addition, in the three hydroclimate regions, the number of statistically significant canonical axes is higher for winter than for fall. From a statistical as well as a climate standpoint (data interpretation), results obtained with winter climate indices appear more interesting than those obtained with fall climate indices. Finally, the number of statistically significant canonical axes varies between regions for winter climate indices. In the eastern hydroclimate region, only the first two axes are statistically significant at the 5% level, while three and four axes are significant, respectively, for the southeast and southwestern regions. For the southeast region, however, the last axis is statistically significant at the 10% level.

The two types of seasonal indices lead to the same results in terms of the relationships between climate variables and climate indices. Due to these similar results, we only interpret results derived with winter climate indices in order to limit the number of tables to analyze and any redundancy in the results. The structure coefficients derived using fall and winter climate indices are shown in

Table 3. Structure coefficients derived using seasonal climate indices (DJFM) in three hydroclimatic regions. Climatic variables.

Variables	East Region				Southeast Region				Southwest Region
	V1	V2	V3	V4	V1	V2	V3	V4	V1	V2	V3	V4
Tmax	0.950	0.147	0.211	−0.180	0.956	0.041	0.146	0.250	0.857	0.276	−0.077	0.429
Tmin	0.750	0.048	0.372	0.547	0.641	0.272	−0.566	0.441	0.383	0.434	−0.354	0.734
TRN	−0.084	0.636	0.721	−0.262	0.431	−0.800	−0.319	0.268	0.283	0.216	−0.930	−0.080
TSN	−0.009	0.701	−0.580	0.415	−0.141	0.094	−0.260	−0.950	0.021	−0.968	−0.032	0.249
EV (%)	36.80	22.99	25.98	14.31	38.26	18.11	12.77	30.78	24.04	31.20	24.92	19.81

EV = Explained variance; the highest structure coefficient values that are statistically significant at the 5% level are shown in bold. TRN = total rainfall; TSN = total snowfall.

(climatic variables) and 4 (climatic indices) for the three hydroclimate regions.

For the eastern region (Table 3 and

Table 4. Structure coefficients derived using winter seasonal (DJFM) climate indices in the three climatic regions. Climatic indices.

Variables	East Region				Southeast Region				Southwest Region
	W1	W2	W3	W4	W1	W2	W3	W4	W1	W2	W3	W4
AMO	0.837	−0.013	−0.217	−0.203	0.460	0.423	0.112	0.145	0.582	−0.051	0.573	0.311
AO	−0.739	0.222	−0.137	−0.342	0.341	−0.722	0.507	−0.078	0.306	−0.089	−0.102	−0.910
NAO	−0.512	0.279	0.127	−0.483	0.520	−0.414	0.266	0.351	0.475	0.300	−0.083	−0.550
PDO	−0.005	−0.618	0.739	−0.244	−0.307	−0.120	0.022	0.842	−0.304	0.903	0.248	−0.057
PNA	0.100	−0.098	0.948	−0.021	0.164	0.000	0.063	0.570	0.153	0.779	0.063	0.031
SOI	−0.106	−0.030	−0.657	0.587	−0.160	0.493	0.469	−0.438	−0.147	−0.711	0.364	0.091
EV (%)	25.50	8.66	32.64	13.26	12.41	18.81	9.41	22.93	13.25	33.80	9.05	20.66

EV = Explained variance; the highest structure coefficient values that are statistically significant at the 5% level are shown in bold. AMO = Atlantic Multidecadal Oscillation; AO = Arctic Oscillation; NAO = North Atlantic Oscillation; PDO = Pacific Decadal Oscillation; PNA = Pacific/North American Oscillation; SOI = Southern Oscillation.

), only the first two canonical axes are statistically significant (see Table 2). It follows that V1 and V2 on one hand, and W1 and W2 on the other hand, are statistically significant canonical variables. V1 is positively correlated with maximum and minimum temperatures and V2 is correlated, albeit moderately, to precipitation (snow and rain). For the canonical variables extracted from the climate index data matrix, W1 is positively correlated with AMO, but negatively correlated with AO. W2 is negatively correlated with PDO. Since V1 is correlated with W1 and V2 is correlated with W2, it follows that maximum and minimum temperatures are positively correlated with AMO and, to a lesser extent, negatively correlated with AO. Precipitation is positively correlated with PDO.

In the southeast region (Table 3 and Table 4), the first three canonical axes are statistically significant (see Table 2). It can be seen from Table 4 that V1 is positively correlated with maximum and minimum temperatures and V2 is negatively correlated with rainfall. The amount of snow is correlated with V4, the canonical axis of which is significant at the 10% level. While the W1 canonical variable is not significantly correlated with any climate index, W2 is negatively correlated with AMO. The last canonical variable, W3, is not correlated with any of the climate indices. W4 is positively correlated with PDO. Unlike the previous hydroclimate region, temperatures in the southeast region are not correlated with any climate index.

In the southwest region (Table 3 and Table 4), V1 is positively correlated with maximum temperature, V2 and V3 are negatively correlated, respectively, with the amounts of snow and rainfall, and V4 is positively correlated with minimum temperature. W1 is not significantly correlated with any climate index. In contrast, W2 is positively correlated with PDO and PNA, but negatively correlated with SOI. W1 and W3 are not correlated with any climate index. W4 is negatively correlated with AO. Therefore, the total amount of winter snow is strongly correlated with PDO, and minimum temperature is correlated with AO. The other two climate variables are not correlated with any climate index.

Table 5. Comparison of the link between climate indices and climate variables in the three southern Quebec hydroclimate regions (1950–2000). Summary of canonical correlation analysis results.

Climate Variables	Eastern Region	Southeast Region	Southwest Region
Tmax	AMO (+)	-	-
Tmin	AMO (+)	-	AO (−)
TRN	PDO (−)	AO (+)	-
TSN	PDO (−)	PDO (−)	PDO (−)

These are the climate indices that show the best correlation with climate variables. (+) = positive correlation; (−) = negative correlation.

presents a summary of canonical correlation analysis results.

This analysis reveals that PDO is the only climate index that influences the amount of winter snow in the three hydroclimate regions. Temperatures and amount of rain are correlated with different climate indices depending on the region.

3.2. Analysis of the Evolution of the Dependence between Climate Indices and Climate Variables

Copula method results are shown in

Table 6. Results of the copula method used to analyze changes in dependence between climate variables and winter climate indices (DJFM) (1950–2000) in the eastern region.

	AMO-Tmax			AMO-Tmin			TRN-PDO			PDO-TSN
	Mn	VC	p-Value	Mn	VC	p-Value	Mn	VC	p-Value	Mn	VC	p-Value
Bagotville	0.658	0.728	0.084	0.640	0.839	0.187	0.384	0.766	0.678	0.951	0.845	0.029
Mont-Joli	0.570	0.774	0.200	0.704	0.694	0.045	0.610	0.768	0.193	0.809	0.905	0.090
Natashqu	0.555	0.782	0.2670	0.7481	0.8262	0.0940	0.6216	0.8053	0.2170	0.6977	0.8745	0.1730
Sept Îles	0.6098	0.7502	0.1500	0.6851	0.7956	0.1260	0.6077	0.8419	0.2360	0.6468	0.8078	0.1780
Ste Rose	0.5474	0.7212	0.2300	0.7197	0.7397	0.0620	0.5275	0.8467	0.4430	0.6244	0.7771	0.1730
Trois-P	0.7355	0.8290	0.0920	0.7225	0.7289	0.0550	0.6008	0.8331	0.2580	0.9202	0.8437	0.0270

Significant values of M_n are shown in bold. VC = critical value.

,

Table 7. Results of the copula method used to analyze changes in dependence between climate variables and winter climate indices (DJFM) (1950–2000) in the southeast region.

Stations	AO-TRN			PDO-TSN
	Mn	VC	p-Value	Mn	VC	p-Value
Coaticook	0.7505	0.7728	0.0660	0.6200	0.8329	0.2470
Granby	0.6945	0.9181	0.1960	0.6365	0.8934	0.2800
Magog	0.8295	0.8492	0.0560	0.7133	0.8265	0.1150
Montréal	0.7503	0.9421	0.1840	0.4785	0.7370	0.4040
Philipsburg	0.6651	0.8833	0.2570	0.5385	0.8427	0.3920
St Malo d’Auckland	0.9994	0.9091	0.0230	0.4434	0.8386	0.6070

Significant values of M_n are shown in bold.

,

Table 8. Results of the copula method used to analyze changes in dependence between climate variables and winter climate indices (DJFM) (1950–2000) in the southwest region.

Stations	AO-Tmin			PDO-TSN
	Mn	VC	p-Value	Mn	VC	p-Value
Chelsea	0.3960	0.9349	0.7890	0.4641	0.8345	0.5470
Nominingue	0.6679	0.8877	0.1980	0.7742	0.9119	0.1250
St-Alban	0.4463	0.9276	0.6600	0.4220	0.8826	0.7040
Saint-Michel-Des-Saints	0.5854	0.9912	0.4360	0.5657	0.7874	0.250
Shawville	0.6945	0.9636	0.2240	0.5152	0.9396	0.5860

and

Table 9. Comparison of Kendall’s tau and coefficient of linear correlation values between climate indices and climate variables before and after shifts.

Stations	Variables	Shift Years	Kendall’s Tau		Linear Coefficient of Correlation
			Before Shift Year	After Shift Year	Before Shift Year	After Shift Year
Bagothville	PDO-TSN	1970	0.316	−0.297	0.443	−0.470
Mont-Joli	AMO-Tmin	1967	0.224	0.529	0.306	0.729
St-Malo	AO-TRN	1992	0.191	0.287	0.193	0.347
Trois-Pistoles	PDO-TRN	1978	0.070	−0.357	0.190	−0.415

TSN = Total snowfall; TRN = Total rainfall.

as well as in Figure 2, Figure 3 and Figure 4. The correlation between climate variables and climate indices changed very little over time at the scale of Quebec during the period from 1950 to 2000. As far as maximum temperatures are concerned, no change is observed over time in the relationship between this climate variable and climate indices for any of the three hydroclimate regions (Table 6, Table 7 and Table 8). For minimum temperatures, a shift is observed only at the Mont-Joli station located in the eastern hydroclimate region (Table 6). The maximum value of M_n exceeded that of V_c in 1967 (Figure 2). In fact, after that year, the correlation between minimum temperatures and the AMO index increased significantly at that station (Table 9). In the other two hydroclimate regions, no change in the relationship between minimum temperatures and climate indices is observed at any of the stations (Table 7 and Table 8). As far as the total amount of rainfall is concerned, a shift in the relationship between it and the PDO index is observed at the Trois-Pistoles station (Table 6), and between it and the AO index, at the Saint-Malo station (Table 7). These shifts occurred in 1978 and 1992, respectively, at the Trois-Pistoles (Table 9) and Saint-Malo d’Auckland (Figure 3 and Table 9) stations, and they also resulted in an increase in coefficient of correlation values (absolute values) at both stations. As far as the total amount of snow is concerned, a shift in the relationship between this climate variable and the PDO index is only observed at the Bagotville station (Table 6). The correlation between these two variables was positive before 1970, but became negative after that year (Figure 4 and Table 9). Finally, in the southwest region, no change is observed in the dependence between climate indices and climate variables.

4. Discussion and Conclusions

Canonical correlation analysis (CCA) was used to analyze simultaneously the correlation between winter temperature and precipitation (climate variables) on one hand, and their correlation with six climate indices, on the other. The first original contribution of this study is the use of the CCA method to simultaneously analyze the correlation between temperature, precipitation and climate indices measured at multiple stations over all the years of measurement. In this way, no information is lost by first grouping the stations (through principal component analysis or hierarchical classification, for instance) into homogeneous climate regions. The second original contribution of the study is the use of the copula method to constrain changes in dependence between climate indices and climate variables over time. Results obtained using these two methods contribute to a better understanding of controversial issues relating to winter climate patterns in southern Quebec. The study yielded the following notable results.

As regards the selection of seasonal climate indices that are best correlated with climate variables, canonical analysis showed that winter climate indices are better correlated with winter temperature and precipitation in the three hydroclimate regions in southern Quebec. Thus, winter indices yielded generally higher canonical correlation coefficients and a larger number of statistically significant canonical axes then fall climate indices. Other authors have already shown that winter NAO indices show a better correlation with winter temperatures and precipitation in Quebec, among other places (e.g., [6,10,12]). In a climate perspective, winter climate conditions are therefore primarily influenced by winter climate indices.

As regards the relationship between the various climate variables (temperature and precipitation), the study shows that there is no statistically significant correlation between winter temperature and precipitation in Quebec. Thus, maximum and minimum temperatures correlated with one another are not correlated with the amount of rain or snow. A commonly accepted assumption is that increasing winter temperature induces a decrease in the amount of snow and a concomitant increase in the amount of rain in winter. However, the link between temperature and amount of precipitation is not as simple as commonly suggested. From a climate standpoint, the amount of snow or rain in winter does not depend exclusively on air temperature variability. Three other factors are also important: the amount of water vapor available in the atmosphere, the frequency and magnitude of convective movements, and the frequency and persistence of cold and warm air masses. In Quebec, winter rainfalls are associated with the presence of warm and humid air masses coming from the United States and/or the Atlantic Ocean. However, the amount of rain produced by these warm and humid air masses during a given season depends primarily on their frequency in the region and the amount of water vapor they contain. Hence, although it has been shown that winter temperatures in Quebec are on the rise, this warming does not necessarily lead to higher rainfall in winter, because this temperature increase affects minimum nighttime temperatures much more strongly than it does maximum daytime temperatures [23]. In winter, although minimum temperatures are increasing, they still most commonly fall below 0 °C, which precludes an increase in the amount or frequency of precipitation as rain during winter in Quebec. Finally, the impact of this nighttime warming on the increased frequency of warm and humid air masses that produce rain has yet to be demonstrated in Quebec. The amount of snow in winter in Quebec primarily depends, for its part, on the frequency of cold fronts associated with cyclogenesis (polar front). However, this frequency is independent of the nighttime warming observed in Quebec as such, because these fronts are part of larger regional and global scale air mass patterns. This may account for the absence of any significant relationship between temperature and amount of snow in Quebec.

As regards the link between climate variables (temperature and precipitation) and climate indices, the study reveals that temperature and precipitation during winter in Quebec are not correlated with the same climate indices. In other words, they are not affected by the same climate factors, which may in part account for the absence of correlation observed between the two climate variables. For temperature, CCA showed that maximum and minimum temperatures are best positively correlated with AMO only in the eastern region, which has a maritime-like climate. In the southeast region, temperatures are not correlated with any climate index, while in the southwest region, only minimum temperatures show a better negative correlation with AO. It should be recalled that the impacts of these two climate indices on temperature and precipitation have been extensively studied in North America (e.g., [36,37,38,39]). AMO describes the temporal variability of Atlantic Ocean surface temperatures in the Northern Hemisphere. This positive correlation suggests that, when Atlantic surface water temperatures are above normal (positive anomaly), maximum and minimum temperatures during winter in Quebec increase, likely as a result of energy transfer from the ocean to the North American continent. As far as AO is concerned, its influence on temperature, precipitation and streamflow in rivers in the northeastern part of North America has already been described by [6]. According to these authors, this influence in eastern Canada is due to the fact that when this index is in a positive phase in winter, temperatures are below normal (negative correlation) due to stronger polar air in the region. However, the influence of this climate index does not seem to affect the whole south shore (eastern and southeast regions) of the St. Lawrence River. As regards total winter rain, they show a better negative correlation with PDO in the eastern region, but a positive correlation with AO in the southeast region. This latter positive correlation goes against the explanation proposed by Kingston et al. [6]. Finally, PDO is the only climate index that shows a better negative correlation with the total amount of snow in all three hydroclimate regions. Brown [10] observed a strong negative correlation between PDO and the total duration of snow cover in the fall in the western part of southern Quebec. The climate impacts of PDO have been analyzed at the global scale by Mantua and Hare [40], according to whom “many of the climate anomalies associated with PDO are broadly similar to those connected with ENSO variation”. Thus, in the North American Great Lakes region (including Quebec), the positive (warm) phase of PDO coincides with negative precipitation anomalies, which accounts for the negative correlation observed between amount of snow and PDO in the three Quebec hydroclimate regions. This negative correlation is also observed between El Niño events and winter precipitation in Quebec, among other places (e.g., [11,12,13,14,17]). According to Shabbar [12], “during El Niño winters, the polar jet stream over northeastern Canada (including Quebec) keeps the Arctic air mass over the high Artic and northeastern Arctic, resulting in colder-than-normal and drier-than-normal winters”.

As far as changes in the correlation between climate indices and climate variables are concerned, the copula method showed that this correlation did not change significantly over time despite a change in air mass circulation in Quebec observed during the 1980s [10]. However, for four of the stations analyzed, there is a significant change in the correlation between climate indices and climate variables. This change, however, which is generally characterized by an increase in coefficient of correlation values, predates the 1980s decade, and its causes have not been elucidated because human factors (different measurement devices, operators, measurement sites, etc.) cannot be called upon to explain them.

In summary, the study highlights the fact that PDO is the only climate index that influences the amount of snow in winter throughout Quebec, although the climate mechanisms that account for this influence remain unclear. At the scale of the province, temperatures and the total amount of rain in winter are not correlated with the same climate indices.