An Assessment of Six Dissimilarity Metrics for Climate Analogs

VOLUME 52
J O U R N A L O F A P P L I E D M E T E O R O L O G Y A N D C L I M A T O L O G Y APRIL 2013
An Assessment of Six Dissimilarity Metrics for Climate Analogs
PATRICK GRENIER
Ouranos Consortium, Montreal, Quebec, Canada
ANNIE-CLAUDE PARENT
ꢀ
Civil and Water Engineering Department, Universite Laval, Quebec City, Quebec, Canada
DAVID HUARD
Ouranos Consortium, Montreal, Quebec, Canada
FRANC¸ OIS ANCTIL
ꢀ
Civil and Water Engineering Department, Universite Laval, Quebec City, Quebec, Canada
DIANE CHAUMONT
Ouranos Consortium, Montreal, Quebec, Canada
(Manuscript received 10 July 2012, in final form 2 November 2012)
ABSTRACT
Spatial analog techniques consist in identifying locations whose historical climate is similar to the anticipated future climate at a reference location. In the process of identifying analogs, one key step is the quantification of the dissimilarity between two climates separated in time and space, which involves the choice of a metric. In this study, six a priori suitable metrics are described (the standardized Euclidean distance, the Kolmogorov–Smirnov statistic, the nearest-neighbor distance, the Zech–Aslan energy statistic, the Friedman–
Rafsky runs statistic, and the Kullback–Leibler divergence) and criteria are proposed and investigated in an attempt to identify the best metric for selecting spatial analogs. The case study involves the use of numerical simulations performed with the Canadian Regional Climate Model (CRCM, version 4.2.3), from which three annual indicators (total precipitation, heating degree-days, and cooling degree-days) are calculated over 30-yr periods (1971–2000 and 2041–70). It is found that the six metrics identify comparable analog regions at a relatively large scale but that best analogs may differ substantially. For best analogs, it is shown that the uncertainty stemming from the metric choice does not generally exceed that stemming from the simulation or model choice. On the basis of the set of criteria considered in this study, the Zech–Aslan energy statistic stands out as the most recommended metric for analog studies, whereas the Friedman–Rafsky runs statistic is the least recommended.
1. Introduction
2010). Climate analogs are typically classified in temporal and spatial categories (Mearns et al. 2001). A temporal analog is a period during which the climate is similar to that during a reference period at the same location.
Such a tool may, for example, be used for seasonal forecasting by inferring that the upcoming season will be similar to the season that followed the analog to the immediate-past period (e.g., Mullan and Thompson
2006). A spatial analog is a location whose climate during a historical period is similar to the anticipated
Climate analog techniques basically involve identifying a context (location/period) the climate or weather of which is similar to that in another context. It may then be inferred that environmental variables not involved in the similarity evaluation are also analogous (Ford et al.
Corresponding author address: P. Grenier, Ouranos Consortium,
550 Sherbrooke West St., 19th floor, West Tower, Montreal, QC
H3A 1B9, Canada.
´future climate at a reference location (e.g., Ramırez-
Villegas et al. 2011). Some of the various methods
E-mail: grenier.patrick@ouranos.ca
DOI: 10.1175/JAMC-D-12-0170.1
Ó 2013 American Meteorological Society
733 734
J O U R N A L O F A P P L I E D M E T E O R O L O G Y A N D C L I M A T O L O G Y VOLUME 52
that imply analogs may be difficult to classify in this spectrum of climate simulations is not explored and because there is no treatment of the simulation biases.
The paper is divided as follows. Section 2 presents the dual scheme.
Spatial analogs are often employed to study the vulnerability of regions to projected variations in temper- experimental setup and the dataset used, including a ature and precipitation, allowing for the examination of detailed description of the six metrics under assessment. potential future conditions. Hallegatte et al. (2007) and Criteria used for comparing the metrics are also intro-
Kopf et al. (2008) used analogs to evaluate the conse- duced. Section 3 presents the results. In section 4, the quences of climate change on European urban areas. robustness of the results is discussed, and an attempt is
´
Ramırez-Villegas et al. (2011) investigated potential made to explain the results on the basis of the structural analogs for a town in Ghana as promising areas for form of the metrics. Section 5 summarizes the main comparative research on agricultural adaptation plans. conclusions. The appendix contains the list of all math-
Kalkstein and Greene (1997) sought analogs for U.S. ematical symbols used in this study. cities, discussing acclimatization to global warming. Ford et al. (2010) extend the definition of spatial analogs to interpolation methods for climate models and to vulner- 2. Materials and methods ability mapping methods, in which variables and climate/ a. Baseline setup and dataset socioeconomic relationships are interpolated or their validity is transferred by analogy. Such methods have
The data that are used represent some aspects of the been used in the context of studies on risks linked to historical (recent past) and future climates at each grid malaria (Martens et al. 1999), vulnerability to flooding cell of an array. The reference location corresponds to
(Fekete 2009), and vulnerability of the agricultural sector the grid cell for which analogs are sought, and candidate
(O’Brien et al. 2004), as well as mortality risk (Brooks locations are all grid cells contained within the investiet al. 2005). gation domain, which includes the reference location
Analog techniques often involve an evaluation of the itself. As presented in Fig. 1, the investigation domain is degree of dissimilarity between climate distributions, restricted to the following geographical limits: 358–508N, defined as ensembles of values describing relevant as- 608–1008W (excluding lake and ocean grid cells). When pects of the climate during a period and at a location. not otherwise stated, the reference location is the grid
This degree is measured by a dissimilarity metric, an cell that includes the city of Montreal (centered at apexpression that is used in this study for referring to the proximately 45.68N, 73.58W and identified by a star in method itself as well as to the scalar value resulting from Fig. 1). For some calculations, 126 reference locations the application of that method to a pair of climates. Var- are used (also identified in Fig. 1). Each location’s cliious metrics have been used and developed in the litera- mate is synthesized by a set of 90 values, accounting for ture for the calculation of spatial analogs. For example, values of three annual indicators over a 30-yr span. The Veloz et al. (2011) and Williams et al. (2007) used the historical and future periods extend over 1971–2000 and standardized Euclidean distance, whereas Kopf et al. 2041–70, respectively.
(2008) employed the Kolmogorov–Smirnov statistic.
Historical climates over the investigation domain
¨
Other metrics used in climate analog research include are relatively close with regard to the Koppen–Geiger the squared-chord distance (Gavin et al. 2003) and a set classification (Peel et al. 2007), with a mostly smooth of successive univariate criteria (Hallegatte et al. 2007). transition from one type to another. These types are
All of those dissimilarity metrics are a priori valuable for mostly Cfa (temperate with a hot summer and no dry
finding analogs. The choice is rarely justified, however, season), Dfa (cold with a hot summer and no dry seaand the sensitivity of the results to the metric is generally son), Dfb (cold with a warm summer and no dry season), left unaddressed. and Dfc (cold with a cold summer and no dry season).
This study consists of an assessment of six selected Montreal’s climate has been of Dfb type during the dissimilarity metrics, following an initial objective of twentieth century and is projected to evolve into Dfa identifying an optimal or best metric on the basis of a type during the current century (Rubel and Kottek 2010). number of criteria. The study further provides an anal- It is thus expected that Montreal’s analogs extend south ysis of relative sensitivities to the metric and climate of its location and that the progression from best to worst simulation choices. The method includes using data from analogs will appear to be relatively smooth on a map. a numerical model for ranking spatial analogs of Mon-
Most studies that characterize climates include at least treal (Quebec, Canada) within a list of candidate loca- one indicator related to temperature and one related to tions. The study does not aim at identifying the ‘‘true’’ precipitation (Kopf et al. 2008). Indicators chosen for spatial analogs for Montreal, because the complete this study are annual total precipitation (ATP), heating

APRIL 2013
G R E N I E R E T A L .
735
FIG. 1. Domain of investigation, with geographical location of the main reference location
(Montreal; star symbol) and of the 126 other reference locations (black diamonds). degree-days (HDD), and cooling degree-days (CDD).
The CRCM grid covers most of North America (the
These indicators are commonly used in climate change ‘‘AMNO’’ grid), and the chosen investigation domain impact studies (Lebassi et al. 2010; Hamlet et al. 2010), contains 3378 grid cells, from which Ncells 5 2404 have including determination of spatial analogs (Kopf et al. more than 50% of land area and are considered as an-
2008). Following the method of the latter study, base alog candidate locations. Each reference or candidate values used for degree-days calculations are set to location is identified to a single grid cell. Simulations
BHDD 5 BCDD 5 188C so that used differ only in terms of the general circulation model
(GCM) providing the CRCM its boundary conditions.
The emissions scenario is Special Report on Emissions
Scenarios (SRES) A2 in all cases except simulation 6
(SIM-06), for which SRES-A1b is used (see Nakicenovic et al. 2000). Table 1 contains relevant information on the selected simulations. The driving GCMs that are used are well known to the international climate science community. Descriptions are provided by Scinocca et al.
(2008) for the Canadian Centre for Climate Modelling and Analysis (CCCma) Coupled General Circulation
365
HDD 5 max[0, BHDD 2 Tavg(i)] and (1)
(2)
åi51
365
CDD 5 max[0, Tavg(i) 2 BCDD],
åi51 where Tavg(i) is the average temperature during yearday i (29 February is ignored). Units for ATP are millimeters per year, with snow precipitation values replaced by their liquid water mass equivalent.
´
Model, version 3 (CGCM3), by Salas-Melia et al. (2005)
´ ´ for the Centre National de Recherches Meteorologiques
(CNRM) Coupled Global Climate Model, version 3
Outputs from nine available simulations of the Canadian Regional Climate Model (CRCM, version 4.2.3) are used for the study. A regional climate model (RCM), which is based on the laws of atmospheric physics, is an appropriate tool for generating realistic climate datasets. The CRCM is a limited-area nested model that is
TABLE 1. Information on regional climate simulations used in this study. based on the fully elastic hydrostatic Euler equations, Simulation
CRCM-v4.2.3 Driver: code for simulation code center-GCM-version this study (historical/future) member which are solved using a noncentered semi-implicit semi-
Lagrangian integration scheme at a 15-min time step. The horizontal grid is defined by a polar-stereographic pro-
SIM-01 aey/afb CCCma-CGCM-v3.1 No. 1
SIM-02 aez/afc CCCma-CGCM-v3.1 No. 2
SIM-03 afa/afd CCCma-CGCM-v3.1 No. 3
SIM-04 aet/aet CCCma-CGCM-v3.1 No. 4
SIM-05 aev/aev CCCma-CGCM-v3.1 No. 5
SIM-06 agw/ahb CNRM-CM-v3 No. 1 jection with 45-km resolution (true at 608N). Vertical
resolution is variable, with 29 Gal-Chen terrain-following
levels. The model includes most of the processes typically
accounted for by RCMs (boundary layer mixing, radiation, cloud schemes, etc.). Caya and Laprise (1999) and ˆ ´ de Elia and Cote (2010) give more technical details about the CRCM.
SIM-07 agx/agx MPI-ECHAM-v5 No. 1
SIM-08 ahi/ahk MPI-ECHAM-v5 No. 2
SIM-09 ahj/ahw MPI-ECHAM-v5 No. 3 736
J O U R N A L O F A P P L I E D M E T E O R O L O G Y A N D C L I M A T O L O G Y VOLUME 52
(CM3), and by Jungclaus et al. (2005) for the German/ One concept involved in four of the metrics (DNN,
European model known as ECHAM5, implemented at DZAE, DFRR, and DKLD) is the standardized Euclidean the Max Planck Institute (MPI). distance SED(Xi, Yj) between points Xi and Yj (which should be distinguished from the metric DSED between the averages of the distributions, which will be detailed b. Dissimilarity metrics
Let one suppose a distribution cA 5 fX1, . . . , Xi, . . . , below). This quantity is
Xng, where Xi is a vector representing a point in a d-
(
)
1/2
[Xi(k) 2 Yj(k)]2 sA(k)sB(k) dimensional space and n is the number of points in the distribution, and a second distribution cB 5 fY1, . . . , d
SED(Xi, Yj) 5
,(3)
åk51
Yj, . . . , Yng, where Yj represents a point in the same space. This nomenclature is sufficient for characterizing where k is a counter for dimensions (indicators) and sA(k) [or sB(k)] is the standard deviation of cA (or cB) in dimension k. Standardization ensures that all indicators are given comparable consideration in the distance calculation. the climate distributions used here, with cA and cB associated with the climates of the reference location and of one candidate location, respectively. The setup implies n 5 30 and d 5 3, with one dimension for each of the three annual indicators (ATP, HDD, and CDD). In the context of climate analogs, a question often arising is whether cA and cB may be considered as two samples from the same parent distribution, which would support viewing them as analogs. A test for that (null) hypothesis would be considered nonparametric if it does not assume a theoretical parent distribution and multivariate whenever d . 1. Such a test often proceeds in two steps, first by calculating a dissimilarity value between samples and second by associating it with a probability
(p value) that two samples from the same distribution have a dissimilarity value at least as high as the one calculated. One classical test suitable for such analysis is the two-sample Kolmogorov–Smirnov (KS) test, in which the dissimilarity (or KS statistic) corresponds to the maximal departure between the cumulative density functions (CDF) of cA and cB and for which correspondences between the KS statistic and a p value do exist for d 5 1, 2, and 3 (e.g., Fasano and Franceschini
1987). In the approach that is presented here, it is considered that rescaling the dissimilarity measure into a p value would bring no additional relevant information.
It would be practically impossible to consider all available dissimilarity metrics, since these are numerous and are often open to a myriad of possible variants. A se-
1) STANDARDIZED EUCLIDEAN DISTANCE
Both samples are separately averaged for each of the d indicators that are used, and the SED between the two average points in the d-dimensional indicator space is next calculated. The measure of dissimilarity DSED is then
!
1/2
ꢀꢁ
2d
[cA](k) 2 [cB](k) sA(k)
DSED(cA, cB) 5
,(4)
åk51 where [](k) represents the average on all n values of a distribution (cA or cB) for the dimension k. This metric considers neither the information from individual points nor the standard deviation of the candidate distribution. The nonuse of sB(k) renders DSED mathematically nonsymmetric—that is, DSED(cA, cB)
6¼ DSED(cB, cA) in most cases, but the word ‘‘distance’’ is kept in the name of this metric. Possible values for DSED are continuous and are bounded below by zero.
2) KOLMOGOROV–SMIRNOV STATISTIC
The one-dimensional KS statistic between two samlection has then been performed, aiming at maximally ples corresponds to the largest difference between their covering the different concepts involved in these metrics CDFs and varies between 0 (for identical distributions) and at representing what had already been used in the and 1 (for nonoverlapping distributions). Although literature on climate analogs at the time of doing the there is no unique and natural extension of the concept analysis. The selection involves six metrics that may of CDF in higher dimensions, several ways of computing operate in any number of dimensions: the standardized a multidimensional KS statistic (DKS) have been pro-
Euclidean distance DSED, the Kolmogorov–Smirnov posed, notably by Friedman and Rafsky (1979), Peacock statistic DKS, the nearest-neighbor (NN) distance DNN ,
(1983), and Fasano and Franceschini (1987). For this the Zech–Aslan energy statistic DZAE, the Friedman– study, the Fasano and Franceschini variant has been
Rafsky (FR) runs statistic DFRR, and the Kullback– chosen, first because it is faster than Peacock’s algo-
Leibler divergence (KLD) DKLD. Prentice (1980) rithm, and second because the fundamental concepts describes other metrics that could have been selected as involved in Friedman and Rafsky’s procedure are alwell. ready used in another metric (DFRR). APRIL 2013 G R E N I E R E T A L .
737 nn
The three-dimensional DKS is calculated as follows.
For each of the pooled data points i 5 1, 2, . . . , 2n (a pooled distribution refers to all points from cA and cB considered as a single ensemble), the space is divided into eight octants via three orthogonal axes crossing at the point i. A dissimilarity zi(oct) is computed for each octant oct 5 1, 2, . . . , 8 as the difference between the number of points from the reference distribution and that from the candidate distribution found in this octant.
The metric DKS is then the largest zi(oct) considering all points and all octants, normalized by the maximal possible dissimilarity, which is n:
1uA 5 uB 5
R[SED(Xi, Xj)], (7b)
å å n(n 2 1) i51 j5i11 nn
1
R[SED(Yi, Yj)], and (7c)
å å n(n 2 1) i51 j5i11 nn
1n2 uAB 5
R[SED(X , Y )]. (7d) i j
å å j5i i51
The weight function is R(SED) 5 2ln(SED). Values of SED below 10212 are set to this threshold (an arbitrary very small value when compared with typically calculated SEDs). Possible values for DZAE are continuous and take negative values for close distributions
(when the pool of analog candidates is devoid of the max(fzg)
DKS
5
.(5) n
Here fzg represents the ensemble of all zi(oct) values. reference climate itself, returned values are practically all
The design of DKS implies n 1 1 possible discretized positive in this study). This metric has no upper bound. values, ranging from 0 to 1 with increments of 1/n.
5) FRIEDMAN–RAFSKY RUNS STATISTIC
3) NEAREST-NEIGHBOR DISTANCE
Friedman and Rafsky (1979) have proposed generalizations in d . 1 dimensions of the Wald–Wolfovitz
(WW) and Smirnov two-sample tests, where the ordering of the data is based on the concept of a minimal spanning tree (MST) in graph theory. In this study, the FR version of the WW test is used. In short, the FR idea is to build the MST of the pooled distribution and to assimilate the number of so-called runs to a dissimilarity metric. The algorithm starts with creating the pooled distribution, with a label identifying the original distribution for each point. Next, each pair of points is assigned an edge with length equal to the SED between them. The MST is then built, following what is known as construction A and is defined by Kruskal (1956) as follows: ‘‘Perform the following step as many times as possible: Among the edges of G [the full set of edges] not yet chosen, choose the shortest edge which does not form any loop with those edges already chosen.’’ The MST is thus the shortest network connecting all points and excluding loops. Once the MST is obtained, the number of runs NFRR (equal to 1 plus the number of edges connecting points from different original distributions) is counted. The FR runs statistic is then
The NN metric consists in counting the number of points NNN of the pooled distribution whose nearest neighbor belongs to the same original distribution (reference or candidate). For one point, the nearest neighbor is determined by calculating its SED to each of the other points [Eq. (3)] and selecting the point providing the smallest SED. The dissimilarity metric DNN is then obtained by dividing NNN by the number of points of the pooled distribution:
NNN
5
DNN .(6)
2n
The metric DNN has 2n 1 1 possible discretized values, ranging from 0 to 1 with increments of 1/(2n).
4) ZECH–ASLAN ENERGY STATISTIC
The dissimilarity measure DZAE, which is based on the work by Zech and Aslan (2003) and Aslan and Zech
(2008), stems from an analogy with electrical charge distributions. The principle is that for one positive unit and one negative unit, each divided among a number of charges tending toward infinity, the total potential energy of the combined samples is a minimum if both samples have the same distribution in space (for a oneover-distance law). The statistic used here is the variant proposed by Zech and Aslan (2003) but with point-topoint distances standardized with sAsB, as in other metrics involving an SED. The metric is given by
NFRR
DFRR 5 1 2
(8)
2n
This metric is bounded by 0 and (n 2 1)/n, with possible increments of 1/(2n). For two identical distributions, one could expect to have DFRR 5 0. A numerical recipe, however, may have difficulty in determining, for example, the shortest distance between that separating points Xn and Xi on one side and that separating points
Xn and Yj on the other side, in the case in which Xi is
DZAE 5 uA(cA) 1 uB(cB) 2 uAB(cA, cB), (7a) with 738
J O U R N A L O F A P P L I E D M E T E O R O L O G Y A N D C L I M A T O L O G Y VOLUME 52