A new family of copula-based concordance orderings of random pairs: Properties and nonparametric tests

The formal assessment of the stochastic dominance of a random pair with respect to another one is a question of interest in the economic analysis of populations. For example, a manager may wonder if the components of a portfolio are more associated than that of another competing portfolio, in which case the former is generally considered more at risk. In this paper, a new family of copula-based concordance orderings in the spirit of increasing convex and concave orderings of random pairs is introduced as a natural extension of the well-known concordance ordering. In addition, a complete statistical methodology to test the stochastic dominance of a random pair with respect to another one according to the new concordance orderings is developed. The proposed tests are nonparametric, consistent against all alternatives, and valid under serially dependent data satisfying the α-mixing assumption. The sampling properties of the tests are investigated with the help of Monte–Carlo simulations and their usefulness is illustrated on real multivariate data. MSC2020 subject classifications: Primary 60E15; secondary 62G10.


Introduction
Stochastic dominance is a key concept in the economic analysis of populations that allows to compare, for instance, incomes, poverty rankings and earnings. The most basic notion is that of first order stochastic dominance of a random variable Y over another variable X, which is defined as E{U (Y )} ≥ E{U (X)} for all utility functions U that are monotonically increasing; it can be shown to be equivalent to P(Y ≤ x) ≤ P(X ≤ x) for all x ∈ R. More generally, as described Such tests have been developed by McFadden (1989) for s = 1, 2 and Anderson (1996), Davidson and Duclos (2000) for s = 1, 2, 3 using a Kolmogorov-Smirnov statistic. Noting that these procedures compare distributions at a fixed number of arbitrary points, so that the tests may be inconsistent, Barrett and Donald (2003) propose a global statistic computed from the empirical distribution functions and where p-values are approximated using either a multiplier or a bootstrap method. A test for the stochastic dominance of degree s = 2 has been proposed by Eubank, Schechtman and Yitzhaki (1993) based on a necessary but not sufficient condition. Some extensions of these procedures have also been proposed. For instance, Linton, Maasoumi and Whang (2005) allow for serial dependence and residuals of linear models, while Linton, Song and Whang (2010) propose to base their decision rule on an improved bootstrap method.
Of a particular interest in this work is the stochastic dominance of a random pair (Y 1 , Y 2 ) over another pair (X 1 , X 2 ). Specifically, (Y 1 , Y 2 ) is said to stochastically dominate (X 1 , X 2 ) in the positive quadrant dependence order if the components of the former are more likely than those of the latter to take small values simultaneously; this is formally stated as P (X 1 ≤ x 1 , X 2 ≤ x 2 ) ≤ P (Y 1 ≤ x 1 , Y 2 ≤ x 2 ) for all (x 1 , x 2 ) ∈ R 2 . Under a setup of fixed marginals, Yanagimoto and Okamoto (1969) and Tchen (1980), among others, state that a pair (X 1 , X 2 ) is stochastically dominated by (Y 1 , Y 2 ) with respect to the concordance order if for all (x 1 , x 2 ) ∈ R 2 , Because of the assumption of fixed marginals, Equation (1) is equivalent to The stochastic ordering as defined in (1) is in fact closely related to the copulas underlying the joint distributions of the pairs to be compared. Specifically, according to a celebrated Theorem of Sklar (1959), there exist copulas C, D : [0, 1] 2 → [0, 1] such that for all (x 1 , x 2 ) ∈ R 2 ,

2395
When the marginal distributions of (X 1 , X 2 ) and (Y 1 , Y 2 ) are continuous, C and D are unique. If in addition P(X 1 ≤ x 1 ) = P(Y 1 ≤ x 1 ) and P(X 2 ≤ x 2 ) = P(Y 2 ≤ x 2 ), then the stochastic dominance of (Y 1 , Y 2 ) over (X 1 , X 2 ) as defined in (1) reduces to the ordering of their respective copulas in the sense that C(u 1 , u 2 ) ≤ D(u 1 , u 2 ) for all (u 1 , u 2 ) ∈ [0, 1] 2 . (2) In other words, the stochastic dominance of pairs is a copula-based notion under the fixed marginals setup.
The stochastic dominance of pairs on a class F of measurable functions on R 2 has been defined by Marshall (1991), for which (Y 1 , Y 2 ) is said to dominate (X 1 , X 2 ), noted (X 1 , X 2 ) F (Y 1 , Y 2 ), if E{φ(Y 1 , Y 2 )} ≥ E{φ(X 1 , X 2 )} for all φ ∈ F for which these expectations exist. Particular cases are the sincreasing convex and concave orderings as proposed by Denuit, Lefèvre and Mesfioui (1999); the latter are related respectively, for s = (s 1 , s 2 ) ∈ N 2 , to the sets F s−ICX = {φ : φ [i1,i2] Note that the stochastic ordering in (1) corresponds to the (1, 1)-increasing convex ordering. This paper's starting point is the ascertainment that the s-increasing convex and concave orderings are not copula-based notions, even when s = (1, 1). Indeed, for (X 1 , X 2 ) (1,1)−ICX (Y 1 , Y 2 ) to imply the ordering of their copulas as stated in (2) necessitates that the two pairs belong to a Fréchet class of bivariate distributions having the same marginals. Moreover, as noted by Fernández-Ponce and Rodríguez-Griñolo (2017), the (2, 2)-convex and concave orderings of two pairs no longer entail the ordering of their corresponding copulas, even within a Fréchet class, unless additional assumptions are made on the marginals. From our point-of-view, this is a clear limitation of these orders when the prime goal is to gain an insight on the relative strength of dependence between two random couples. This paper has two main goals: (1) Introduce a family of copula-based stochastic orderings between random pairs (having possibly different marginals) and explore their properties; (2) Develop a formal statistical methodology to assess the stochastic ordering of two bivariate populations with respect to the new class, in a spirit similar to Barrett and Donald (2003) for the ordering of univariate distributions.
The paper is organized as follows. Section 2 defines the new family of copulabased orderings called s-concordance and explores their properties. Section 3 describes some interpretations of these orders, including the establishment of interesting links with conditional and unconditional versions of Spearman's rank correlation. Section 4 and Section 5 develop a whole statistical methodology to assess the s-concordance ordering of two populations when serial data are available. Section 6 investigates the size and power of the new tests under various data-generating scenarios. Section 7 shows how to adapt these tools to deal with bivariate pairs coming from the same multivariate population, and illustrates the idea on Cook & Johnson's Uranium exploration dataset; another illustration on exchange rate currencies is detailed as well. The proofs are relegated to two appendices and all the code is freely available on www.uqtr.ca/MyMatlabWebpage.

Definition of the s-concordance orders
The new family of concordance orderings that will be defined is based on the s-increasing convex ordering s−ICX as defined by Denuit, Lefèvre and Mesfioui (1999). Specifically, (X 1 , Now let (X 1 , X 2 ) be a random pair from a distribution function with continuous marginal distributions F 1 (x 1 ) = P(X 1 ≤ x 1 ) and F 2 (x 2 ) = P(X 2 ≤ x 2 ), and unique copula C such that C(u 1 , u 2 ) = P{X 1 ≤ F −1 1 (u 1 ), X 2 ≤ F −1 2 (u 2 )}. Consider another pair (Y 1 , Y 2 ) with continuous marginals G 1 , G 2 and unique copula D. LetF 1 ,F 2 ,Ḡ 1 andḠ 2 be the marginal survival functions of, respectively, X 1 , X 2 , Y 1 and Y 2 . The following definitions of orthant s-concordance orderings are based on the s-increasing convex ordering of uniformized pairs.
Similarly, the pair (Y 1 , Y 2 ) stochastically dominates (X 1 , X 2 ) according to the upper orthant s-concordance ordering, noted (X 1 , The lower orthant s-concordance ordering could as well be defined in terms of the s-increasing concave ordering s−ICV as described by Denuit, Lefèvre and Mesfioui (1999). Specifically, because of the duality relationship between the s-increasing concave and convex orders expressed in their Proposition 2.3, one has that (X 1 , X 2 ) s− o (Y 1 , Y 2 ) could have been defined similarly as (G 1 (Y 1 ), G 2 (Y 2 )) s−ICV (F 1 (X 1 ), F 2 (X 2 )).
By construction, the new orthant s-concordance orderings are marginal-free since they depend only on the copula of the pairs being compared. Indeed, for an arbitrary pair (X 1 , X 2 ) from a joint distribution with continuous marginals F 1 and F 2 , it is well known that (F 1 (X 1 ), F 2 (X 2 )) ∼ C; also, (F 1 (X 1 ),F 2 (X 2 )) ∼ C, where C(u 1 , u 2 ) = u 1 +u 2 −1+C(1−u 1 , 1−u 2 ) is the survival copula of C. This is to be contrasted with the s-increasing convex ordering that imply some ordering of the marginals. Specifically, if (X 1 , X 2 ) s−ICX (Y 1 , Y 2 ), then X j sj −ICX Y j for j = 1, 2, thus involving the marginal distributions.
From the characterization of the (1, 1)-increasing convex order, one deduces that (X 1 , , the lower orthant (1, 1)-concordance ordering corresponds to the usual concordance ordering; the exact same conclusion applies to the upper orthant (1, 1)-concordance ordering.

Basic properties
The first result concerns the duality that exists between s− o and s−uo .
As stated next, the new class of orderings is hierarchical.
Since copulas are invariant under monotone increasing transformations of the marginals, it is expected that the lower and upper orthant s-concordance orderings be invariant under such mappings. The next result establishes this basic fact, and explores situations involving monotone decreasing transformations.
, if ψ 1 is decreasing and ψ 2 is increasing.
When C = C, where C is the survival copula of C, it is said that C is a radially symmetric copula. In generic terms, radial symmetry means that the lower tail of (the density of) C has the same form as the upper tail. Knowing that, the next result which establishes the equivalence between s− o and s−uo under radially symmetric dependence structures will come as no surprise.

The fixed marginals setup and beyond
Suppose that the pairs (X 1 , X 2 ) and (Y 1 , Y 2 ) belong to the same Fréchet class of bivariate distributions with marginals F 1 and F 2 . In that case, the lower and upper orthant (1, 1)-concordance orderings are equivalent to the (1, 1)increasing convex order. Things are not as straightforward when s = (1, 1). As one can deduce from Theorem 3 of Fernández-Ponce and Rodríguez-Griñolo (2017), (X 1 , X 2 ) (2,2)−uo (Y 1 , Y 2 ) entails (X 1 , X 2 ) (2,2)−ICX (Y 1 , Y 2 ) not only if the pairs belong to the same Fréchet class, but if in addition F 1 , F 2 have decreasing densities. On the other side, a consequence of Theorem 4 of these authors is that ( if F 1 , F 2 have increasing densities. One then deduces that the only case where Proposition 2.1. Let (X 1 , X 2 ) and (Y 1 , Y 2 ) be in the same Fréchet class of bivariate distributions with continuous marginals F 1 and F 2 .
Lower orthant s-concordance One recovers Theorem 3 of Fernández-Ponce and Rodríguez-Griñolo (2017) by letting s = (2, 2) in parts (i) and (iii) of Proposition 2.1, since the 2-concavity of F −1 1 , F −1 2 means that F 1 , F 2 have decreasing densities, while their 2-convexity is equivalent to having increasing densities. For similar reasons, their Theorem 4 is a special case of parts (ii) and (iv) of Proposition 2.1.

s-concordance orderings of popular copula families
Consider the Normal, Clayton, Gumbel and Plackett copulas whose expressions are given in Table 1. These models have been extensively used for bivariate copula modeling and their properties are well-known; see Nelsen (2006) and Joe (2015), for instance. For example, the family of Normal copulas share with the Plackett dependence structures the property of radial symmetry. Table 1 The Normal, Clayton, Gumbel and Plackett copulas

Copula
Expression of the copula Parameter space a Φ θ is the cdf of the bivariate Normal with correlation θ b g θ (u 1 , u 2 ) = 1 − θ + 2θ(u 1 + u 2 ) These four copula families are parametrized such that C θ (u 1 , u 2 ) ≤ C θ (u 1 , u 2 ) for all (u 1 , u 2 ) ∈ (0, 1) 2 when θ ≤ θ . In view of Property 2.2, they are also ordered with respect to s-concordance for any s ∈ N 2 . Things become less clear, and in fact more interesting, when the goal is to stochastically compare two copulas C and D that belong to different parametric families. To this end, let (X 1 , X 2 ) and (Y 1 , Y 2 ) be random pairs with copulas C and D, respectively. Table 2 reports scenarios when (X 1 , X 2 ) and (Y 1 , Y 2 ) cannot be ordered according to the usual concordance, i.e. (X 1 , . In order to standardize the comparisons, each model has been parametrized in terms of its associated Kendall's tau, i.e. As an example of a situation where (1,1)− o and (2,2)− o , consider D being the Clayton copula with τ D = 1/3; in that case, (X 1 , X 2 ) (1,1)− o (Y 1 , Y 2 ) and Table 2 Some scenarios where the respective copulas C and D of random pairs (X 1 , X 2 ) and if the copula C of (X 1 , X 2 ) belongs to the Gumbel family with τ C ∈ (.0071, .3346]. It also happens when τ D = 2/3 and τ C ∈ (.0281, .6640], and more particularly when τ C = τ D = 1/3. Another example occurs when D is the Clayton copula with τ D = 1/3 and C belongs to the Plackett family with τ C ∈ (.1556, .3330].

Characterization of s− o and s−uo
A characterization of the s-increasing convex ordering of random pair was deduced by Denuit, Lefèvre and Mesfioui (2003) in terms of iterated distributions. Such characterizations for the s-concordance orders s− o and s−uo are developed here. The latter will prove useful later to formally test for s-concordance ordering. To this end, let ∞ ([0, 1] 2 ) be the space of bounded functions on [0, 1] 2 . For g ∈ ∞ ([0, 1] 2 ), define for each i ∈ N 2 the operator J i : As formally stated in the following result, the lower orthant s-concordance ordering of random pairs can be seen as a functional of the difference between their respectively copulas. Similarly, the upper orthant s-concordance ordering appears as a functional of their associated survival copulas.
Proposition 2.2. Let (X 1 , X 2 ) and (Y 1 , Y 2 ) be random pairs with continuous marginals and respective copulas C and D.
There is only one condition needed to establish the s-concordance ordering of two pairs when s ≤ (2, 2), since E s = {s} in that case. Proposition (2.2) then reduces to (X 1 , (4) holds.

The (2, 1) and (1, 2)-concordance orderings
As outlined by Denuit and Mesfioui (2017), the (2, 1)-increasing concave order- While condition (i) expresses the usual stochastic dominance of X 2 over Y 2 , the second condition compares the strength of the corresponding conditional shortfalls of the pairs. Specifically, (t 1 − X 1 ) + I(X 2 ≤ t 2 ) vanishes given that X 2 is larger than the threshold t 2 , so that the shortfall (t 1 − X 1 ) + with respect to the threshold t 1 becomes irrelevant. This shows some sort of compensation between the components of (t 1 − X 1 ) + I(X 2 ≤ t 2 ). Also, as shown in Proposition 3.2 of Denuit and Mesfioui (2017), the (2,1)-concave order characterizes the Rothschild-Stiglitz type of increase in risk as introduced by Guo et al. (2016).
As was noted after the statement of Definition 2.1, the lower orthant s- ). Since the marginal distributions of the pairs to be compared are, by construction, uniform on (0, 1), condition (i) becomes irrelevant and the lower orthant (2, 1)-concordance (X 1 , Clearly, the above inequality holds when the components of the pair (Y 1 , Y 2 ) are more associated than those of (X 1 , X 2 ), since then, (t 1 −F 1 (X 1 )) + I(F 2 (X 2 ) ≤ t 2 ) tends to vanish more frequently than (t 1 − G 1 (Y 1 )) + I(G 2 (Y 2 ) ≤ t 2 ). Otherwise, similar interpretations as those above can be made, but at the level of the dependence structures of the pairs that are being compared.

Consequences of s− o on Spearman's rho and other concordance measures
Several measures of dependence are concordance measures in the sense given by Scarsini (1984); see also Nelsen (2002). Generally, they can be expressed in terms of the concordance operator between two copulas as defined by where U 1 ∼ C 1 and U 2 ∼ C 2 are independent pairs. For example, the Kendall and Spearman measures of dependence of a random pair (X 1 , X 2 ) with copula C can be expressed respectively as where Π(u 1 , u 2 ) = u 1 u 2 and M (u 1 , u 2 ) = min(u 1 , u 2 ) are the copulas of independence and perfect positive dependence, respectively. In fact, concordance measures are closely linked to the concordance ordering (1,1)− o . On one side, it can be shown that if

Consequences of (2,2)− o and (2,2)−uo on conditional versions of Spearman's rho
As noted in (3), (X 1 , X 2 ) (2,2)− o (Y 1 , Y 2 ) entails the ordering of the lower orthant integrated copulas; in view of (4), (X 1 , X 2 ) (2,2)−uo (Y 1 , Y 2 ) implies a similar ordering, but with respect to upper orthant integration. In particular, In fact, this inequality holds when (X 1 , , so that the s-concordance ordering of two pairs always implies the ordering of their corresponding Spearman's rho. But still more can be said about Spearman's rho under the (2, 2)-concordance ordering. To this end, first define a version of the concordance operator in (7) constrained to the lower rectangle [0, Replacing the concordance operator Q by Q u in the definitions of Kendall and Spearman measures of association yields conditional versions of these concordance measures. Doing so for Spearman's rho, one obtains This is exactly the bivariate version of the conditional Spearman's rho as defined by Schmid and Schmidt (2007) while letting d = 2 and g := I(· ≤ u 1 , · ≤ u 2 ) in their Equation (4). As a consequence, (X 1 , Hence, based on (4), one can conclude that (

Tests of s-concordance orderings
As stated in the Introduction, this paper's second aim is to provide a nonparametric statistical methodology to formally assess the s-concordance ordering of two bivariate populations. Even in the case of the usual concordance ordering, i.e. when s = (1, 1), no procedure has been developed yet. A paper by Cebriàn, Denuit and Scaillet (2004) entitled "Testing for concordance ordering" is seemingly achieving this, but in fact the goal of these authors is to compare one bivariate population's joint distribution with a pre-specified parametric model.

Null and alternative hypotheses
For a fixed s ∈ N 2 , the goal is to test for the stochastic dominance of a random pair (Y 1 , Y 2 ) over (X 1 , X 2 ) with respect to the lower orthant s-concordance ordering. In other words, one wants to test for In view of Property 2.1, the methodology that will be developed in the sequel can easily be adapted to test for the upper orthant s-concordance ordering s−uo by considering the lower orthant dominance of (−Y 1 , −Y 2 ) over (−X 1 , −X 2 ). Now a reformulation of the null and alternative hypotheses in (9) will prove useful. To this end, let C and D be the copulas of (X 1 , X 2 ) and (Y 1 , Y 2 ), respectively. In view of Equation (6) in Proposition 2.2, the null hypothesis H 1 . This suggests basing a measure of s-concordance on some functional of Measuring the lower orthant s-concordance of a pair (X 1 , X 2 ) with respect to (Y 1 , Y 2 ) can then be based on Combination rules others than taking the maximum over i ∈ E s could be considered as well, e.g. the sum. However, the most interesting situations are those 0 holds if and only if Θ (s) κ,(C,D) = 0. The null and alternative hypotheses stated in (9) may therefore be reformulated alternatively as

Test statistics and asymptotics under α-mixing
This subsection provides an empirical version of Θ (s) κ, (C,D) and investigates its asymptotic behavior under a setup of serially dependent observations. Specifically, it will be assumed that the data at hand are realizations of strongly stationary processes that satisfy the α-mixing assumption. This notion is very general, as it is shared by many popular time series models like autoregressive and GARCH processes. Specifically, following, e.g., Bradley (2005), Carrasco and Chen (2002) or Rio (2000), consider a process (Z t ) t∈Z and define Now, let (X 11 , X 12 ), . . . , (X n1 , X n2 ) be a realization of a strongly stationary process (X t1 , X t2 ) t∈Z that is α-mixing; also assume that for all t ∈ Z, the marginal distributions of (X t1 , X t2 ) are continuous and C is its unique copula. Consider another sample (Y 11 , Y 12 ), . . ., (Y m1 , Y m2 ), independent of the first one, that is a realization of a strongly stationary α-mixing process (Y t1 , Y t2 ) t∈Z with copula D. Under these conditions, nonparametric estimators of C and D are provided by the empirical copulas, namely where n U i1 (resp. m V i1 ) is the rank of X i1 (resp. Y i1 ) among X 11 , . . . , X n1 (resp. Y 11 , . . . , Y m1 ), and similarly for n U i2 (resp. m V i2 ). An empirical plug-in version of Θ (s) κ, (C,D) defined in (11) is then As a first step, the following proposition establishes the asymptotic behavior of J i (·, ·; C n − D m ) for any i ∈ N 2 . This result can be seen somewhat as a copula version of Lemma 1 of Barrett and Donald (2003) about the iterated sample and population cdf's in the univariate case, i.e.
Before stating the result, the concept of a regular copula is reminded.

Decision rule, significance level and consistency
Based on the null and alternative hypotheses of s-concordance ordering as reformulated in (12), it is suggested to reject H  (14), and since μ κ (rg) = rμ κ (g) for r ∈ R + , Since under the null hypothesis H It follows that under H , u whereβ (s) κ is the survival function of max i∈Es μ κ (L (i) (u ). Therefore, the test whose decision rule is to reject H (s) has an asymptotic type I error of at most α. Hence, the test has a significance level equals to α as understood by Lehmann (1986) in the case of a composite null hypothesis. In the current context, it means that the test based on Θ (s) κ,(n,m) will have a rejection rate that will never exceed α for any pair of copulas C, D such that Θ (s) κ,(C,D) = 0; the asymptotic level is exactly α when C = D. A violation of H (s) 0 means that there is a set B ⊂ [0, 1] 2 of non-null Lebesgue measure such that for some i ∈ E s , J i (u The test based on Θ (s) κ,(n,m) is therefore consistent under general alternatives.

Estimation of the critical value
In order to estimate the asymptotic critical value (β (s) κ ) −1 (α), one needs to estimate the distribution function β (s) ) ≤ x}. This is not an easy task, since the limit process L (i) depends on the unknown copulas C and D under H (s) 0 . The adopted strategy will be based on the multiplier bootstrap for empirical processes as described by Kosorok (2008) and adapted to empirical copulas under α-mixing by Bücher and Ruppert (2013); the latter is a generalization to time series of the multiplier method for empirical copulas as described for instance by Rémillard and Scaillet (2009).
Definition 5.1 (Serial multipliers). A serial multiplier sample associated to sample data of size n is a realization ξ = (ξ 1 , . . . , ξ n ) of a strictly stationary process (ξ t ) t∈Z that is independent of the data process and such that (i) ξ t is independent of ξ t+h for all |h| ≥ r n , where r ∈ R is a constant and as n → ∞, n → ∞ and n /n → 0; (ii) all central moments of ξ t are bounded, E(ξ t ) = 1 and cov(ξ t , ξ t+h ) = φ(h/ n ), where φ is bounded and symmetric around zero.
In order to approximate the distribution β provides with B independent replicates of Θ (s) κ, (n,m) . The critical value (β (s) κ ) −1 (α) is finally estimated by the (1 − α)-th empirical percentile of these multiplier bootstrap replicates. The consistency of this estimator is a straightforward consequence of Proposition 5.1.

Approximation of the test statistics
From the recursive definition of J i in Equation (5), one can easily establish by induction that when g(u 1 , u 2 ) = I(a 1 ≤ u 1 , a 2 ≤ u 2 ),

2409
It follows that and similarly for J i (u 1 , u 2 ; C m ). Now the test statistic Θ (s) κ,(n,m) defined in (13), which involves the computation of μ κ {J i (·, ·; C n − D m ) for each i ∈ E s , will be based on an approximation of J i (·, ·; C n − D m ) = J i (·, ·; C n ) − J i (·, ·; D m ) on a K × K grid of [0, 1] 2 for K ∈ N taken large enough to ensure a satisfactory numerical accuracy. Specifically, consider the product intervals From the definition of μ κ in (10), an easily computable approximation is

Preliminaries
The ability of the tests of (1, 1)-concordance and (2, 2)-concordance orderings to keep their nominal 5% level, as well as their power against selected alternatives, is studied in this section. The test statistics that will be investigated are those based on the functional μ κ defined in (10) when κ ∈ {1, 2, ∞}. While μ 1 and μ 2 are global distance measures of the Cramér-von Mises type, μ ∞ (g) = sup (u1,u2)∈[0,1] 2 (g(u 1 , u 2 )) + is related to the Kolmogorov-Smirnov distance. These functionals are approximated on a grid of size K × K = 25 × 25, as described in subsections 5.2-5.3. For all the results that will be presented, the probabilities of rejection of the null hypothesis have been estimated from 1 000 replicates, each based on B = 1 000 multiplier bootstrap samples. The estimation of the partial derivativesĊ 1 ,Ċ 2 , D 1 andḊ 2 are based on a finite-difference estimator. Explicitly, the estimator ofĊ 1 is given and similarly forĊ n2 ,Ḋ m1 andḊ m2 ; in the sequel, b = 1.

Tests of (1, 1)-concordance ordering in the i.i.d. case
One first considers the size and power of the tests of concordance ordering, i.e. of H , under the i.i.d. setup where one is willing to assume serial independence in both time series. In that case, the multiplier random variables can be taken i.i.d., i.e. n = m = 1, and Exponential with mean 1. The results on the estimated probabilities of rejection of the null hypothesis are found in Table 3 for (n, m) ∈ {(100, 100), (100, 200), (200, 200)}.
The four scenarios in the upper part of Table 3 are under the null hypothesis. In the case of the first two scenarios, i.e. when C and D are normal copulas, the null hypothesis holds strictly in the sense that C = D; in that case, the three tests are rather good at holding their 5% nominal level, except for small sample sizes (n = m = 100) and a high level of dependence (τ C = τ D = 2/3). The other two scenarios are cases where H (1,1) 0 holds, but not strictly, i.e. C = D; it is therefore not surprising that the probabilities of rejection are far below the 5% nominal level of the tests.
All the other entries in Table 3 are obtained under the alternative hypothesis. Overall, the power is an increasing function of the departure from H (1,1) 0 as measured by Θ (1,1) ∞, (C,D) . However, when the value of Θ (1,1) ∞, (C,D) is small (say < .02), the power is often below the nominal level. Otherwise, the power of the three tests is very good and increases with the sample sizes, a consequence of their asymptotic consistency. It is hard to identify a statistic that is uniformly the best, but the tests based on μ 1 and μ 2 are generally preferable to that using μ ∞ . The test based on μ 2 is generally the most powerful when D is a Clayton copula, and that using μ 1 is the best when D is Plackett.

Tests of (2,2)-concordance ordering in the i.i.d. case
A study similar to that presented in subsection 6.2 has been conducted for testing H (2,2) 0 against H (2,2) 1 . Since (1,1)− o implies (2,2)− o , only the scenarios for which (1,1)− o , as identified in Table 1, have been considered . The results on the estimated probabilities of rejection of the null hypothesis are found in Table 4. First note that the four scenarios in the upper part of Table 4 concern non strict null hypotheses, i.e. situations where C = D; this explains why the probabilities of rejection are below the 5% nominal level in that case. The remaining six entries in the bottom part of Table 4 correspond to scenarios under H (2,2) 1 . As expected, the power of the tests increases with the sample sizes, as well as according to the value of Θ (2,2) ∞, (C,D) . Here, the test based on the functional μ ∞ is clearly the most powerful against all kinds of alternatives.

Size and power under Gaussian serial dependence
When there is temporal dependence between the observations, one has to rely on the serial multiplier method. In the sequel, one follows Bücher and Ruppert (2013) and let (ζ j ) j∈Z be a process of independent Gamma(q, q) ran- Table 3. Percentages of rejection, as estimated from 1 000 replicates, of the tests for the (1, 1)-concordance ordering hypothesis based on Θ   Table 4. Percentages of rejection, as estimated from 1 000 replicates, of the tests for the (2, 2)-concordance ordering hypothesis based on Θ  dom variables with q = (2 n − 1) −1 , where the bandwidth parameter is set to n = 1.1 n 1/4 . Then, for each j ∈ {1, . . . , n}, one defines ξ j as the mean of ζ j− n +1 , . . . , ζ j+ n −1 . One proceeds similarly for the second sample of size m.
The results reported in the upper part of Table 5 have been obtained for serial data generated from the lag-1 Gaussian autoregressive process where θ ∈ (−1, 1) and (ε t1 , ε t2 ) t∈Z is a process of centred independent Normal pairs with unit variances and correlation ρ ∈ (−1, 1). The middle part of Table 5 concerns the Gaussian moving-average process of order one defined by These processes are stationary and parameterized in such a way that the copula of (X t1 , X t2 ) is Normal with parameter ρ. The level of dependence of the generated time series is managed by the value of Kendall's tau via the well-known relationships ρ C = sin(πτ C /2) and ρ D = sin(πτ D /2). Only the results when τ D = 1/3 are presented, since those when τ D = 2/3 lead to similar conclusions. For comparison purposes, the results that have been already obtained under serial independence are reported here in the bottom part of Table 5.
Overall, the results are very similar to those in Table 3 in the case of i.i.d. data. In particular, if one looks at the cases when θ = 0, there is no price to pay in terms of size and power by wrongly assuming serial dependence. When τ C = 1/3, the null hypothesis H (1,1) 0 holds strictly; in that case, the three tests keep their 5% nominal level well, whatever the kind and level of serial dependence. This is an indication that the serial multiplier method is good at replicating the behavior of the test statistics under the null hypothesis.

Comparisons with a test of s-increasing convex order
Suppose a setup of fixed marginals, i.e. of pairs (X 1 , X 2 ) and (Y 1 , Y 2 ) that belong to the same Fréchet class of bivariate distributions with margins F 1 , F 2 . In that case, according to Proposition 2.1 (ii), the s-increasing convex dominance of (−Y 1 , −Y 2 ) over (−X 1 , −X 2 ) entails that (X 1 , X 2 ) s− o (Y 1 , Y 2 ) as long as F 1 is s 1 -concave and F 2 is s 2 -concave. Therefore, if in practice one is willing to assume these constraints on the marginal distributions, simplified tests for H s 0 against H s 1 could be based on the observations themselves, and not on their ranks, as is mandatory when working at the level of copulas. This procedure would avoid the estimation of the partial derivatives of copulas due to the fact that these terms are missing in the (simpler) asymptotic expression of the limit.
Specifically, suppose that the goal is to test for H s 0 : To this end, let (X 11 , X 12 ), . . . , (X n1 , X n2 ) and (Y 11 , Y 12 ), . . ., (Y m1 , Y m2 ) be realizations of the Table 5. Percentages of rejection, as estimated from 1 000 replicates, of the tests for the (1, 1)-concordance ordering hypothesis based on Θ It is indeed the case since the mapping from −X j to 1 − X j is linear increasing. What it means is that the procedure developed for testing s− o can be performed on the pairs ( X i1 , X i2 ) and ( Y i1 , Y i2 ) instead of the pairs of standardized ranks ( U i1 , U i2 ) and ( V i1 , V i2 ). However, the multiplier method has to be performed by removing the part involving the partial derivatives.
In order to evaluate how such an alternative procedure performs, some simulations have been made when s = (1, 1) and s = (2, 2) in case the marginal distribution is the Beta(1,3); the latter has a decreasing density, hence is 2concave. For simplicity, only the functional μ ∞ has been considered and the corresponding test statistics is noted Θ s ∞, (n,m) . The results are found in Table 6, where for the sake of comparison, the corresponding results for Θ s ∞, (n,m) extracted from Table 3 and Table 4, have been reproduced.
Looking at Table 6, one first notes that the test based on Θ s ∞,(n,m) holds its 5% nominal level rather well. However, somewhat surprisingly, the test based on Θ s ∞, (n,m) is much more powerful than its counterpart derived under additional assumptions on the marginals. A more detailed investigation of tests of s-increasing convex ordering would be worth of interest. Nevertheless, based on these simulation results, it seems that bringing more information about the marginal distributions do not transfer into a more powerful procedure.

Adaptation of the methodology for stochastic comparisons within the same multivariate population
The statistical methodology developed in this work can easily be adapted for the comparison of two pairs (X j , X k ) and (X j , X k ) that come as marginals of a dvariate random vector X = (X 1 , . . . , X d ) with continuous marginals F 1 , . . . , F d . If K : [0, 1] d → [0, 1] is the unique copula of X, then C(u 1 , u 2 ) = K(u (jk) ) and In that context, the copula estimators are respectively C n (u 1 , u 2 ) = K n (u (jk) ) and D n (u 1 , u 2 ) = K n (u (j k ) ), where K n is the d-dimensional empirical copula computed from X 1 , . . . , X n . According to Bücher and Volgushev (2013), as long as K is regular, i.e.K = ∂K/∂u exists and is continuous on {u ∈ [0, 1] d : 0 < u < 1} for each ∈ {1, . . . , d}, and under the same α-mixing conditions than Table 6. Percentages of rejection, as estimated from 1 000 replicates, of the tests for the (1, 1)-concordance ordering (upper panel)  those in Proposition 4.1, the empirical process K n = √ n(K n − K) converges weakly in the space ∞ ([0, 1] d ) to a limit of the form In this expression, B C is a Gaussian process on [0, 1] d with mean zero such that for U = (F 1 (X 1 ), . . . , F d (X d )), the covariance function of B C is One can then derive an adapted version of Proposition 4.1. A multiplier version of K n based on a serial multiplier sample ξ can be defined in the same line as those for C n and D m ; the counterpart of Proposition 5.1 is straightforward to obtain. From an implementation perspective, the only necessary adjustment consists in using the same multiplier sample ξ for both datasets of n pairs.

Cook & Johnson's Uranium exploration data
The Uranium exploration dataset has been first considered by Johnson (1981, 1986). It consists of concentrations of seven chemical elements measured on n = 655 water samples collected from the Montrose quadrangle of western Colorado (USA). All these samples are independent from each other. The following analyses will focus on four of these variables, namely Potassium (K), Caesium (Cs), Scandium (Sc) and Titanium (Ti).
The histograms and the scatterplots, both of the original data X 1 , . . . , X 655 and of the standardized ranks U 1 , . . . , U 655 , are found in Figure 1. Looking at the histograms of the four variables, it is clear that they are marginally quite different. Hence, if the goal is to perform stochastic comparisons among some of the pairs, it cannot reasonably be assumed that they belong to the same Fréchet class. The s-concordance orderings, which assume nothing on the marginal distributions (apart from being continuous, which is the case here), are therefore well-suited for these data.
The first analysis concerns the stochastic comparison of (Cs,Ti) with (K,Cs) using the adapted methodology of subsection 7.1 with K = 25 and B = 10, 000 i.i.d. multiplier samples; the estimation of the partial derivatives is done by letting b = 1, since the tests performed with b = 3 yielded very similar values. If one looks at Figure 1, it seems that (Cs,Ti) (1,1)− o (K,Cs) cannot hold; this is confirmed by the results of the tests that are found in Table 7 when s = (1, 1). This is also confirmed, to a certain extent, by the values of the empirical Kendall's tau, namely τ n (Cs,Ti) = .279 and τ n (K,Cs) = .200.
Nevertheless, their dependence structures can somewhat be ordered if one looks at other levels. Hence, while the ordering (2,1)− o is still rejected, the null hypotheses H   on the functionals μ 1 and μ 2 ; these hypotheses are however rejected by the test based on μ ∞ . In view of the link between the lower orthant (2,2)-concordance ordering and conditional Spearman's rho established in subsection 3.3, this sug-gests that the dependence level of (K,Cs), as measured by Spearman's rho, can be larger than that of (Cs,Ti), and vice versa, when one restricts to some lower corners of [0, 1] 2 . Replacing C with C n and developing formula (8) yields as an empirical Spearman's rho conditioned on [0, u k1 ] × [0, u k2 ], with u k = (k − 1/2)/K. The top panel of Figure 2 shows ρ Sp (Cs, Ti). It can be seen that for u 1 ≤ .7, say, Spearman's rho is larger for (K,Cs) compared to (Cs,Ti), often significantly (curve above 0); in the complementary region, it is for (Cs,Ti) that Spearman's rho is larger (curve below 0), but by a much less amount. These features could have been anticipated from the results of the tests when s = (2, 2). That the test based on μ ∞ has rejected the null hypothesis of a (2, 2)-concordance ordering may be explained by the fact that this functional can be strongly influenced by local discrepancies, while μ 1 , μ 2 are global distances. The pair (Cs,Ti) has also been stochastically compared to (Cs,Sc); the results in Table 7 are similar to those of the previous analysis. Thus, while the null hypothesis (Cs,Ti) (1,1)− o (Cs,Sc) is clearly rejected, there is nevertheless some sort of ordering at the level of conditional measures of association like Spearman's rho that can be clearly observed on the bottom panel of Figure 2.

Evolution of exchange rates
Another illustration concerns the n = 228 exchange rates of the Euro (EUR), Canada (CAN), Australia (AUS), New Zealand (NZE) and Japan (JAP) currencies as measured monthly with respect to US dollar between January 1999 and December 2017. The five time series are found at the top of Figure 3 (Japan currency has been divided by 100). The series are clearly not marginally stationary; however, since the lag-1 differentiated series are reasonably stationary (see bottom of Figure 3), the latter will be considered for the upcoming analyses. A look at the scatterplots of the lag-1 differentiated time series shown in the lower triangle of Figure 4 indicates a possible radial symmetry structure; this is confirmed by the test of radial symmetry of Bahraoui and Quessy (2017), where based on 1,000 multiplier samples, the test's p-value is estimated to 21,2%; note however that the test assumes serial independence (to date, no test exists to deal with that situation). Radial symmetry means that the orderings s− o and s−uo between two-pairs are equivalent. Because the relationship of Japan with other currencies is quite low, except maybe with the Australian currency, the former has been excluded of the following analysis; these low dependence levels can be seen from the values of Kendall's tau, i.e. τ n (EUR, JAP) = 0.341, τ n (CAN, JAP) = 0.152, τ n (AUS, JAP) = 0.226 and τ n (NZE, JAP) = 0.165. The results of the test based on Θ (s) ∞,(n,m) for each of the six possible comparisons of non-overlapping pairs of (EUR,CAN,AUS,NZE) are presented in Table 8. Here, the number of serial multiplier samples, as described in subsection 6.4, has been set to B = 1, 000 with n = 1.1 × 228 1/4 = 4. One of the conclusions is that the pair (AUS,NZE) significantly dominates (EUR,CAN) at the 5% level according to (1, 1)-concordance; this feature was expected from the respective scatterplots of standardized ranks in Figure 4. A similar conclu-sion can be made about the dominance of (CAN,NZE) over (EUR,AUS), and of (CAN,AUS) over (EUR,NZE); in these cases, however, the use of a formal test prove crucial, since the conclusion could hardly be based on looking at the scatterplots only.

Discussion
In this paper, a new family of stochastic orders that allow for marginal-free comparisons between random pairs have been introduced; these orders generalize the usual concordance ordering. The construction of this hierarchical family of orders is rooted around the concept of s-increasing convex orders computed at the level of the copula that uniquely characterizes the dependence in a random couple. It has been shown, in particular, how these orders are related to Spearman's measure of association. Also, a complete set of statistical tools has been developed to formally assess the stochastic dominance of a random pair on another pair; the proposed framework is quite general, as it allows for serially dependent data, and can also accommodate the case when the two pairs are subvectors drawn from the same multivariate population. The analysis performed on the classical Uranium exploration data is typical of the information that can be extracted in a multivariate dataset when looking from the point-of-view of s-concordance orders. Hence, while the proposed statistical methodology has clearly discarded the usual concordance ordering hypothesis for being too strong, it allowed to establish a relationship at the level of the less restrictive order s = (2, 2); in turn, this can be interpreted as the dominance of one pair on another at the level of conditional Spearman's rho when one restricts to some lower corners of [0, 1].
In a future investigation, it would be interesting to generalize the notions of positive quadrant dependence (PQD) and negative quadrant dependence (NQD) with respect to the new class of s-concordance orderings. Specifically, one could define (X 1 , X 2 ) to be s-PQD (resp. s-NQD) if (X ⊥ 1 , X ⊥ 2 ) s− o (X 1 , X 2 ) (resp. (X 1 , X 2 ) s− o (X ⊥ 1 , X ⊥ 2 )), where (X ⊥ 1 , X ⊥ 2 ) is a copy of (X 1 , X 2 ), but with independent components. The statistical tools of Section 4 and Section 5 could then be adapted in order to provide not only interesting extensions of the tests of Scaillet (2005) and Gijbels, Omelka and Sznajder (2010) to s-PQD and s-NQD, but also provide new test statistics that are valid under serial data and/or when the pairs are from the same multivariate population.
Another fruitful avenue of research would be to design an alternative bootstrap procedure in order to ensure an exact asymptotic size for the tests. A promising way would be to adapt to the current context a bootstrap procedure proposed by Linton, Song and Whang (2010) for testing univariate stochastic dominance that improves the power of the tests by Barrett and Donald (2003). Such a version for the tests developed in this work is, however, far from being straightforward. In fact, since one is working at the level of copulas, the complexity of the asymptotics is increased due to the use of ranks. In addition, the methodology would have to be adapted to serial data.

A.2. Proof of Property 2.2
Given that s−ICX is a hierarchic order, see Equation (2.21) of Denuit, Lefèvre and Mesfioui (1999), the proof is straightforward.

B.3. Proof of Proposition 4.1
Under the conditions stated, Bücher and Volgushev (2013) obtained the weak convergence in ∞ ([0, 1] 2 ) of C n = √ n(C n −C) to C; similarly, D m = √ m(D m − D) converges weakly to D. Strictly speaking, because the two samples are independent, these two convergences are simultaneous, i.e. the pair of processes (C n , D m ) converges to (C, D). From there, where ω n,m = n/(n + m). Since the operator J i is continuous, one can then conclude that L (i) n,m L (i) = J i (·, ·; √ 1 − ω C − √ ω D). That this convergence happens jointly for any i ∈ N 2 is obvious from the definition of L