Yuri Goegebeur

Armelle Guillou

Jing Qin

Robust estimation of the conditional stable tail dependence function

Keywords: Empirical processes, Local estimation, Multivariate extreme value statistics, Robustness, Stable tail dependence function

We propose a robust estimator of the stable tail dependence function in the case where random covariates are recorded. Under suitable assumptions, we derive the finite dimensional weak convergence of the estimator properly normalized. The performance of our estimator in terms of efficiency and robustness is illustrated through a simulation study. Our methodology is applied on a real dataset of sale prices of residential properties.

Introduction

A topic of central interest in multivariate extreme values is to measure the strength of dependence in the extremes. This can be done by using some coefficients of tail dependence or some functions, among them the Pickands dependence function or the stable tail dependence function. In the present paper, we focus on this latter function introduced by [START_REF] Huang | Statistics of bivariate extremes[END_REF], and we estimate it when the random variables of main interest are recorded along with random covariates, related to the target variables. That means that we are in the regression context where our objective is to estimate the stable tail dependence function between the response variables conditional on the covariates. This leads to the concept of conditional stable tail dependence function. Additionally, since in practice some outliers may occur in real datasets with a disturbing effect on the estimates of dependencies, we propose an estimator which is robust against observations that are atypical for the extreme dependence structure of the models under consideration. In other words, our contribution in this paper is to introduce a robust estimator of the conditional stable tail dependence function. This topic has been only partially considered in the recent literature, e.g., by [START_REF] Escobar-Bach | Bias-corrected and robust estimation of the bivariate stable tail dependence function[END_REF][START_REF] Escobar-Bach | Local estimation of the conditional stable tail dependence function[END_REF]. See also [START_REF] Gardes | Nonparametric estimation of the conditional tail copula[END_REF], de Carvalho (2016), [START_REF] Castro | Spectral Density Regression for Bivariate Extremes[END_REF], [START_REF] Castro | Time-Varying Extreme Value Dependence with Application to Leading European Stock Markets[END_REF], [START_REF] Mhalla | Regression type models for extremal dependence[END_REF], or [START_REF] Escobar-Bach | Bias correction in conditional multivariate extremes[END_REF] and [START_REF] Goegebeur | Robust nonparametric estimation of the conditional tail dependence coefficient[END_REF]. Concretely, throughout the paper, we denote by pY p1q , Y p2q q a bivariate random vector recorded along with a random covariate X P R d . Let }.} denote some norm on R d and B x prq the closed ball with respect to }.} centered at x and radius r ą 0. For j " 1, 2, we denote by F j p.|xq, the continuous conditional distribution function of Y pjq given X " x, by f X the density function of the covariate X and by x 0 a reference position such that x 0 P IntpS X q, the interior of the support S X Ă R d , which is assumed to be non-empty. Our aim in this paper is to estimate the conditional stable tail dependence function defined as lim rÑ8 rP ´1 ´F1 pY p1q |Xq ď r ´1y 1 or 1 ´F2 pY p2q |Xq ď r ´1y 2 ˇˇX " x ¯" Lpy 1 , y 2 |xq in a robust way, where we assume that the above limit exists for all x P S X . By assuming continuous marginal conditional distributions for Y p1q and Y p2q , this condition is essentially a condition on the tail behavior of the copula function underlying the joint conditional distribution of Y p1q and Y p2q given X " x. As such, the stable tail dependence function contains information about the dependence in extremes. To reach our goal, we assume that, for all x P S X , our bivariate random vector pY p1q , Y p2q q satisfies the model P ´1 ´F1 pY p1q |Xq ď y 1 , 1 ´F2 pY p2q |Xq ď y 2 ˇˇX " x ¯" y d 1 pxq 1 y d 2 pxq 2 gpy 1 , y 2 |xq p1 `δpy 1 , y 2 |xqq , (1) for any py 1 , y 2 q P r0, 1s 2 ztp0, 0qu, where d 1 pxq, d 2 pxq are positive and continuous functions such that d 1 pxq `d2 pxq " 1, gpy 1 , y 2 |xq is continuous in py 1 , y 2 , xq and homogeneous of order 0 in py 1 , y 2 q, and δpy 1 , y 2 |xq is a function of constant sign in the neighbourhood of zero, with |δp., .|xq| a bivariate regularly varying function, that is, there exists a function ξp., .|xq such that lim rÓ0 |δ|pry 1 , ry 2 |xq |δ|pr, r|xq " ξpy 1 , y 2 |xq, for all py 1 , y 2 q P r0, 8q 2 ztp0, 0qu, where the convergence is uniform in py 1 , y 2 q P p0, T s 2 and x P B x 0 pζq, for any T ą 0 and ζ ą 0. Also, ξpy 1 , y 2 |xq is assumed to be continuous in py 1 , y 2 , xq and homogeneous of order βpxq ą 0 in py 1 , y 2 q. Model (1) was introduced in a simpler context without covariates in [START_REF] Dutang | Robust and bias-corrected estimation of the coefficient of tail dependence[END_REF] and [START_REF] Escobar-Bach | Bias-corrected and robust estimation of the bivariate stable tail dependence function[END_REF], see also [START_REF] Beirlant | Bias-reduced estimators for bivariate tail modelling[END_REF], and it has its roots in [START_REF] Ledford | Modelling dependence within joint tail regions[END_REF]. Essentially (1) is a further assumption on the tail copula that underlies the joint distribution of pY p1q , Y p2q q, conditional on X " x. The approach followed in the present paper to the estimation of Lpy 1 , y 2 |xq will be non-parametric and based on local estimation of extreme value models in a neighborhood of the point of interest in the covariate space.

In real data analysis, outliers appear occasionally, and in such contexts robust methods are crucial to avoid poor performance of the usual estimators, like the maximum likelihood estimator. A huge literature exists on outlier detection and robust estimation methods, following the seminal contributions of [START_REF] Huber | Robust statistics[END_REF] and [START_REF] Hampel | Robust statistics: The approach based on influence functions[END_REF]. In the extreme value context, Dell' Aquila and Embrechts (2006) discussed some methodological aspects related to robust estimation. In particular, they showed how robust methods can improve the quality of data analysis by providing information on the atypical observations, and on the deviation from the structure of the underlying model, while guaranteeing good statistical properties of the resulting estimators computed on the complete dataset. Our aim in this paper is to estimate the conditional stable tail dependence function in a robust way, to prevent possible isolated outliers from completely disturbing the estimate. In the multivariate context, observations can be outlying with respect to the dependency structure, in the sense that they do not follow the pattern set by the majority of the data, and hence they disturb the estimation of the dependence structure. To achieve the robustness, the density power divergence criterion initially proposed by [START_REF] Basu | Robust and efficient estimation by minimizing a density power divergence[END_REF] will be used. It is defined between two density functions f and h as follows ∆ α pf, hq :"

$ ' ' ' ' & ' ' ' ' % ż R " h 1`α pyq ´ˆ1 `1 α ˙hα pyqf pyq `1 α f 1`α pyq  dy, α ą 0, ż R log f pyq hpyq f pyqdy, α " 0.
Here, f is assumed to be the true (typically unknown) density of the data, whereas h is a parametric model, depending on a vector of parameters which is estimated by minimizing the empirical version of ∆ α pf, hq. This estimator is called the minimum density power divergence (MDPD) estimator. Unlike existing methods such as minimum Hellinger distance estimation, Basu et al.'s (1998) approach avoids the use of nonparametric density estimation and the associated problem of bandwidth selection. This MDPD method only depends on a tuning parameter α which can be viewed as a trade-off between robustness and asymptotic efficiency of the estimators. When α " 0, the density power divergence is the Kullback-Leibler divergence [START_REF] Kullback | On information and sufficiency[END_REF] and the method reduces to the maximum likelihood estimation. When α " 1 it corresponds to the mean squared error or L 2 ´divergence. As such the minimum density power divergence represents a whole family of divergences, indexed by the parameter α ě 0. Thus, we introduce a nonparametric and robust estimator for the conditional stable tail dependence function when the data come from a conditional distribution whose dependence structure converges to that of a conditional extreme value distribution. Compared to related recent literature on estimation of extremal dependence the differences are as follows. [START_REF] Escobar-Bach | Bias-corrected and robust estimation of the bivariate stable tail dependence function[END_REF] consider robust estimation of the stable tail dependence function though in a context without covariates. Escobar-Bach et al. (2018a) derive a robust estimator for the Pickands dependence function in a context with covariates, where they assume that the underlying conditional distribution has a conditional extreme value copula. In [START_REF] Escobar-Bach | Local estimation of the conditional stable tail dependence function[END_REF] an estimator for the stable tail dependence function is introduced in a regression context, where it is assumed that the underlying conditional distribution has a dependence function that converges to that of a conditional extreme value distribution, though their estimator is not robust with respect to outlying observations. [START_REF] Goegebeur | Robust nonparametric estimation of the conditional tail dependence coefficient[END_REF] discuss a robust estimator for the coefficient of tail dependence in the context of random covariates. In some sense, the present paper can be viewed as a follow-up of the latter paper, although the problem considered in [START_REF] Goegebeur | Robust nonparametric estimation of the conditional tail dependence coefficient[END_REF] is simpler than the one considered in the present paper since now the aim is to estimate a dependence function rather than a single parameter. Also, in Goegebeur et al. (2020) a deterministic, i.e., non-random, intermediate threshold is used, while in the present paper we consider the more realistic situation where the intermediate threshold is taken as an intermediate conditional quantile, which complicates the asymptotic analysis considerably.

Our paper is organized as follows. In Section 2, we assume that both conditional marginal distribution functions are known, and we propose a robust estimator of the conditional stable tail dependence function for which we establish the finite dimensional weak convergence. Then, in Section 3, we consider the more realistic situation where the conditional marginal distribution functions are unknown. We estimate again in a robust way the conditional stable tail dependence function and we derive similar results as in the previous section, under some additional assumptions. The finite sample performance of our estimator in terms of efficiency and robustness is illustrated in Section 4 on a simulation experiment. Finally, in Section 5, we apply our methodology to a real dataset of sale prices of residential properties. Some concluding remarks are proposed in Section 6. The proofs of some of the main results are postponed to Section 7, whereas the others and those of some auxiliary results are given in the online Supplementary Material.

2 Estimation of Lpy 1 , y 2 |x 0 q in case of known margins

For convenience, assume that the conditional marginal distributions F 1 p.|xq and F 2 p.|xq are unit Pareto and let Z t :" mintY p1q , t 1´t Y p2q u for t :" y 1 y 1 `y2 , 0 ă t ă 1. Then, according to model (1), the conditional survival function of Z t given X " x, denoted by F Zt p.|xq, is a conditional Pareto-type model of the following form:

F Zt py|xq " G t pxqy ´1 p1 `δt py|xqq , (2)
where

G t pxq :" ˆt 1 ´t ˙d2 pxq g ˆ1, t 1 ´t ˇˇx ˙, δ t py|xq :" δ ˆ1 y , t 1 ´t 1 y ˇˇx ˙.
Note that δ t p.|xq is regularly varying at infinity with index ´βpxq, i.e., δ t puy|xq{ δ t pu|xq Ñ y ´βpxq as u Ñ 8 for all y ą 0. Additionally, we assume in the sequel the following classical condition. Assumption pD t q For all x P S X , the conditional survival function of Z t given by ( 2) is such that |δ t p.|xq| is normalized regularly varying with index ´βpxq, i.e.,

δ t py|xq " C t pxq exp ˆż y 1 ε t pu|xq u du ˙,
with C t pxq P R and ε t py|xq Ñ ´βpxq as y Ñ 8. Moreover, we assume y Þ Ñ ε t py|xq to be a continuous function.

Under Assumption pD t q we have that δ t py|xq is differentiable and hence F Zt p.|xq has a density function. This condition is a restriction of the Karamata representation of regularly varying functions (see, e.g., Corollary 2.1 in [START_REF] Resnick | Heavy-tailed phenomena[END_REF]. We now turn to the estimation of Lpy 1 , y 2 |x 0 q. The above implies that

Lpy 1 , y 2 |x 0 q " lim ∆Ó0 1 ∆ ! P ´1 ´F1 pY p1q |x 0 q ď ∆ y 1 ˇˇX " x 0 ¯`P ´1 ´F2 pY p2q |x 0 q ď ∆ y 2 ˇˇX " x 0 P ´1 ´F1 pY p1q |x 0 q ď ∆ y 1 , 1 ´F2 pY p2q |x 0 q ď ∆ y 2 ˇˇX " x 0 ¯) " y 1 `y2 ´lim ∆Ó0 1 ∆ P ˆZ1´t ě 1 ∆ y 1 ˇˇX " x 0 " y 1 `y2 ´y1 G 1´t px 0 q.
Thus, estimating the conditional stable tail dependence function requires the estimation of G 1´t px 0 q. To reach this goal, note that from (2), we deduce that

G 1´t px 0 q " k n U Z 1´t pn{k|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q , (3)
where U Z 1´t p.|x 0 q is the conditional tail quantile function defined as U Z 1´t p.|x 0 q :" infty : F Z 1´t py|x 0 q ě 1 ´1{.u, and k is an intermediate sequence such that k Ñ 8 and k{n Ñ 0. If p U Z 1´t p.|x 0 q and p δ 1´t p.|x 0 q are estimators for U Z 1´t p.|x 0 q and δ 1´t p.|x 0 q, respectively, then by the plug-in method we derive the following estimator for G 1´t px 0 q:

p G 1´t,k px 0 q :" k n p U Z 1´t pn{k|x 0 q 1 `p δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q , ( 4
)
which yields a simple estimator for the conditional stable tail dependence function:

p L k py 1 , y 2 |x 0 q :" y 1 `y2 ´y1 p G 1´t,k px 0 q. ( 5
)
Recall that we want to propose a robust estimator. To this aim, we will adjust the MDPD criterion to the local estimation context. Remark that F Z 1´t belongs to the class of distribution functions of [START_REF] Beirlant | Second-order refined peaks-over-threshold modelling for heavy-tailed distributions[END_REF]. Thus, the distribution of the relative excesses Z 1´t {u n given Z 1´t ą u n can, for u n large, be approximated by an extended Pareto distribution (EPD) function given by Hpy; δ 1´t pu n |x 0 q, βpx 0 qq :"

# 1 ´y´1 " 1 `δ1´t pu n |x 0 q ´1 ´y´βpx 0 q ¯ı´1 , y ą 1, 0, y ď 1,
where δ 1´t pu n |x 0 q ą maxt´1, ´1{βpx 0 qu. Moreover, using Proposition 2.3 in [START_REF] Beirlant | Second-order refined peaks-over-threshold modelling for heavy-tailed distributions[END_REF] the approximation error is uniformly opδ 1´t pu n |x 0 qq for u n Ñ 8. Using this property, one can estimate δ 1´t pu n |x 0 q with the MDPD approach as follows.

Starting from pY p1q i , Y p2q i , X i q, i " 1, . . . , n, independent copies of pY p1q , Y p2q , Xq, we obtain pZ 1´t,i , X i q, i " 1, . . . , n, independent copies of pZ 1´t , Xq, and fit the density function h associated with H and defined for y ą 1 as h py; δ 1´t pu n |x 0 q, βpx 0 qq :" y ´2 " 1 `δ1´t pu n |x 0 q ´1 ´y´βpx 0 q ¯ı´2

ˆ"1 `δ1´t pu n |x 0 q ´1 ´p1 ´βpx 0 qqy ´βpx 0 q ¯ı , locally to the relative excesses Z 1´t,i { p U Z 1´t pn{k|x 0 q, i " 1, . . . , n, given that Z 1´t,i ą p U Z 1´t pn{k|x 0 q. Here, p U Z 1´t pn{k|x 0 q is the natural estimator for U Z 1´t pn{k|x 0 q defined as p U Z 1´t pn{k|x 0 q :" infty : p F Z 1´t py|x 0 q ě 1 ´k{nu where, for p F Z 1´t py|x 0 q, we use the kernel-type estimator

p F Z 1´t py|x 0 q :" 1 n ř n i"1 K hn px 0 ´Xi q1l tZ 1´t,i ďyu 1 n ř n i"1 K hn px 0 ´Xi q , ( 6
)
with K hn p.q :" Kp.{h n q{h d n , K a joint density on R d and h n a positive non-random sequence satisfying h n Ñ 0 as n Ñ 8. This leads to the minimum density power divergence estimator, p δ n,1´t :" p δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q, for δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q, and defined as the point minimizing the empirical density power divergence, that is, for α ą 0

p ∆ α,1´t pδ 1´t |x 0 q :" 1 k n ÿ i"1 K hn px 0 ´Xi q # ż 8 1 h 1`α py; δ 1´t , βqdy ´ˆ1 `1 α ˙hα ˜Z1´t,i p U Z1´t pn{k|x 0 q ; δ 1´t , β ¸+ 1l tZ1´t,ią p U Z 1´t pn{k|x0qu .
In our proposed procedure, we only estimate δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q with the density power divergence criterion, while the second order rate parameter will be fixed at a value, denoted β, which can be either the true value or a mis-specified one. Fixing this second order rate parameter β at some value is quite common when fitting second order models like (2) to data, see, e.g., [START_REF] Feuerverger | Estimating a tail exponent by modelling departure from a Pareto distribution[END_REF], [START_REF] Gomes | Bias-reduction and explicit semi-parametric estimation of the tail index[END_REF], [START_REF] Dutang | Robust and bias-corrected estimation of the coefficient of tail dependence[END_REF][START_REF] Escobar-Bach | Bias-corrected and robust estimation of the bivariate stable tail dependence function[END_REF].

To study the asymptotic behaviour of our estimator defined in (5), we need to assume some classical conditions due to the regression context, which are nowadays well-known in the conditional extreme value framework.

First, the density f X and the functions appearing in F Z 1´t py|xq need to satisfy the following Hölder conditions.

Assumption pH 1´t q There exist positive constants

M f X , M G 1´t , M C 1´t , M ε 1´t , η f X , η G 1´t , η C 1´t
and η ε 1´t , such that for all x, z P S X :

|f X pxq ´fX pzq| ď M f X }x ´z} η f X , |G 1´t pxq ´G1´t pzq| ď M G 1´t }x ´z} η G 1´t , |C 1´t pxq ´C1´t pzq| ď M C 1´t }x ´z} η C 1´t , sup yě1 |ε 1´t py|xq ´ε1´t py|zq| ď M ε 1´t }x ´z} ηε 1´t .
Then, we have also to impose a condition on the kernel function K, which is a standard condition in local estimation.

Assumption pK 1 q K is a bounded density function on R d , with support S K included in the unit ball in R d .

We have now all the ingredients to state the main result of this section, namely the joint weak convergence of the estimators p L k py 1,j , y 2,j |x 0 q, j " 1, . . . , J, after proper normalization. In the sequel, weak convergence is denoted by .

Theorem 2.1 Assume pD 1´t j q and pH 1´t j q for j " 1, . . . , J, pD 0.5 q, pH 0.5 q, pK 1 q, x 0 P IntpS X q with f X px 0 q ą 0, and y Þ Ñ F Z 1´t j py|x 0 q, j " 1, . . . , J, are strictly increasing. Consider sequences k Ñ 8 and h n Ñ 0 as n Ñ 8 such that k{n Ñ 0, kh d n Ñ 8, h

ηε 1´t 1 ^¨¨¨^ηε 1´t J ^ηε 0.5 n log n k Ñ 0, a kh d n h η f X ^ηG 1´t 1 ^¨¨¨^η G 1´t J n Ñ 0, a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q| Ñ 0, j " 1, . . . , J. Then, for n Ñ 8, we have, b kh d n ¨p L k py 1,1 , y 2,1 |x 0 q ´Lpy 1,1 , y 2,1 |x 0 q . . . p L k py 1,J , y 2,J |x 0 q ´Lpy 1,J , y 2,J |x 0 q ‹ ' ¨L1 . . . L J ‹ ',
where, for j " 1, . . . , J,

L j :" ´y1,j G 1´t j px 0 q W 1´t j p1q f X px 0 q `y1,j G 1´t j px 0 qc " 2α ż 1 0 " W 1´t j pzq z ´W1´t j p1q  z 2α dz ´p1 `βqp2α `βq ż 1 0 " W 1´t j pzq z ´W1´t j p1q  z 2α`β dz * , c :" p1 `2αqp1 `2α `βqp1 `2α `2βq β 2 p1 `β `4α 2 `2αβq 1 f X px 0 q ,
and W 1´t j pyq, j " 1, . . . , J, are zero centered Gaussian processes with

EpW 1´t j pyqW 1´t j 1 pyqq " }K} 2 2 f X px 0 q ˜max ˜G1´t j px 0 q y , G 1´t j 1 px 0 q y ¸¸´d 1 px 0 q ˆ˜max ˜tj 1 ´tj G 1´t j px 0 q y , t j 1 1 ´tj 1 G 1´t j 1 px 0 q y ¸¸´d 2 px 0 q ˆg ¨1 max ˆG1´t j px 0 q y , G 1´t j 1 px 0 q y ˙, 1 max ˆtj 1´t j G 1´t j px 0 q y , t j 1 1´t j 1 G 1´t j 1 px 0 q y ˙ˇˇˇˇˇˇˇx 0 ‹ ‹ ' .
In practice, the conditional marginal distribution functions F 1 p.|x 0 q and F 2 p.|x 0 q are unknown. The aim of the next section is to extend our results to this new framework. The proof of the new theorem will be given in the online Supplementary Material.

3 Estimation of Lpy 1 , y 2 |x 0 q in case of unknown margins

We consider the general framework where F 1 p.|xq and F 2 p.|xq are unknown conditional margins. We want to mimic what has been done in the previous section by transforming the margins into approximate unit Pareto distributions. To this aim, we define

q Z 1´t :" min # 1 1 ´p F n,1 pY p1q |Xq , 1 ´t t 1 1 ´p F n,2 pY p2q |Xq + ,
where the estimators p F n,j , j " 1, 2, are defined as p F n,j py|x 0 q :"

1 n ř n i"1 K cn px 0 ´Xi q1l tY pjq i ďyu 1 n ř n i"1 K cn px 0 ´Xi q , ( 7
)
with c n a positive non-random sequence satisfying c n Ñ 0 as n Ñ 8. Note that this estimator has the same form as (6) but with a bandwidth c n which needs to be different from h n and the kernel used here is, for simplicity, the same as the one used in the MDPD method.

A similar estimator as the one defined in (5) can be proposed for the robust estimation of the conditional stable tail dependence function in case of unknown margins:

q L k py 1 , y 2 |x 0 q :" y 1 `y2 ´y1 q G 1´t,k px 0 q, where q G 1´t,k px 0 q :" k n p U q Z 1´t pn{k|x 0 q 1 `q δ 1´t p p U q Z 1´t pn{k|x 0 q|x 0 q
, with q δ 1´t p p U q Z 1´t pn{k|x 0 q|x 0 q the MDPD estimator based on q Z 1´t .

To establish the joint weak convergence of the estimators q L k py 1,j , y 2,j |x 0 q, j " 1, . . . , J, after proper normalization, we need to impose again some assumptions, in particular a Hölder-type condition on each marginal conditional distribution function F j , j " 1, 2.

Assumption pF m q. There exist M F j ą 0 and η F j ą 0 such that |F j py|xq ´Fj py|zq| ď M F j }x ´z} η F j , for all y P R, all px, zq P S X ˆSX and j " 1, 2.

Concerning the kernel K a stronger assumption than pK 1 q is needed.

Assumption pK 2 q. K satisfies Assumption pK 1 q, there exists δ, m ą 0 such that B 0 pδq Ă S K and Kpuq ě m for all u P B 0 pδq, and K belongs to the linear span (the set of finite linear combinations) of functions k ě 0 satisfying the following property: the subgraph of k, tps, uq : kpsq ě uu, can be represented as a finite number of Boolean operations among sets of the form tps, uq : qps, uq ě ϕpuqu, where q is a polynomial on R d ˆR and ϕ is an arbitrary real function.

The latter assumption has already been used in [START_REF] Giné | Rates of strong uniform consistency for multivariate kernel density estimators[END_REF] or [START_REF] Giné | Weighted uniform consistency of kernel density estimators[END_REF].

As stated in these contributions, it is satisfied by Kpxq " φtapxqu, a being a polynomial and φ a bounded real function of bounded variation (see, e.g., [START_REF] Nolan | U-processes: rates of convergence[END_REF]. This is also the case, e.g., if the graph of K is a pyramid (truncated or not), or if K " 1l r´1,1s d , etc.

The main result of the paper is given in the below theorem.

Theorem 3.1 Assume that there exists b ą 0 such that f X pxq ě b, @x P S X Ă R d , f X is bounded, pD 1´t j q, pH 1´t j q for j " 1, . . . , J, pD 0.5 q, pH 0.5 q, pK 2 q, pF m q hold, and that y Þ Ñ F Z 1´t j py|x 0 q, j " 1, . . . , J, are strictly increasing at x 0 P IntpS X q non-empty. Consider sequences k Ñ 8, h n Ñ 0 and c n Ñ 0 as n Ñ 8, such that k{n Ñ 0, kh d n Ñ 8, h

ηε 1´t 1 ^¨¨¨^ηε 1´t J ^ηε 0.5 n log n k Ñ 0, a kh d n h η f X ^ηG 1´t 1 ^¨¨¨^η G 1´t J n Ñ 0, a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q| Ñ 0, j " 1, .
. . , J. Assume also that there exists an ε ą 0 such that for n sufficiently large inf xPS X λ ptu P B 0 p1q : x ´cn u P S X uq ą ε,

where λ denotes the Lebesgue measure, and for some q ą 1 and

0 ă η ă minpη F 1 , η F 2 q n d h d n k max ˜d | log c n | q nc d n , c η n ¸ÝÑ 0, as n Ñ 8. (9)
Then, we have

b kh d n ¨q L k py 1,1 , y 2,1 |x 0 q ´Lpy 1,1 , y 2,1 |x 0 q . . . q L k py 1,J , y 2,J |x 0 q ´Lpy 1,J , y 2,J |x 0 q ‹ ' ¨L1 . . . L J ‹ ',
where L j , j " 1, . . . , J, are defined as in Theorem 2.1.

Note that the conditions ( 8) and ( 9) are needed to measure the discrepancy between the conditional distribution function F j py|xq and its empirical kernel version p F n,j py|xq, j " 1, 2, uniformly in px, yq, see, e.g., Lemma 3.1 in Escobar-Bach et al. (2018a).

Simulation study

The aim of this section is to illustrate the performance of our robust estimator q L k py 1 , y 2 |x 0 q with a simulation study. The two following models will be considered. Model 1. The logistic copula model Cpu 1 , u 2 |xq " e ´rp´ln u 1 q x `p´ln u 2 q x s 1{x , u 1 , u 2 P r0, 1s, x ě 2.

We take X " U r2, 10s, and combine this copula model with Fréchet distributions for Y p1q and Y p2q : F j pyq " e ´y´1{γ j , y ą 0, j " 1, 2.

We set γ 1 " 0.25 and γ 2 " 0.5. This model corresponds to Lpy 1 , y 2 |xq " py x 1 `yx 2 q 1{x . Model 2. The conditional distribution of pY p1q , Y p2q q given X " x is that of

p|Z 1 | γ 1 pxq , |Z 2 | γ 2 pxq q,
where pZ 1 , Z 2 q follow a bivariate standard Cauchy distribution with density function

f pz 1 , z 2 q " 1 2π p1 `z2 1 `z2 2 q ´3{2 , pz 1 , z 2 q P R 2 .
We take X " U r0, 1s, and set γ 1 pxq " 0.25 `0.125 sinp2πxq,

γ 2 pxq " 0.1 `0.1x.
This model corresponds to Lpy 1 , y 2 |xq " a y 2 1 `y2 2 . Contamination will be introduced by adding observations that follow a different dependency structure, namely contamination according to the following mixture model

F ε py 1 , y 2 |xq " p1 ´εq F py 1 , y 2 |xq `ε F c py 1 , y 2 |xq,
where ε denotes the fraction of contamination, F is the distribution function of Model 1 or Model 2 described above, and F c is the contamination distribution function. Given X " x, the distribution function F c used with Model 1 is

F c py 1 , y 2 |xq " 1 2 te ´y´1 1 `e´y ´1 2 u1l ty 1 ě0,y 2 ě0u ,
which corresponds to a contamination on the axes, whereas for Model 2, where the dependence is not very strong, we propose to use a diagonal contamination to highlight the effect of contamination. Concretely, that means a distribution function F c of the following form

F c py 1 , y 2 |xq " e ´tminpy 1 ,y 2 qu ´1 , y 1 , y 2 ą 0,
corresponding to the distribution function of completely dependent unit Fréchet random variables. We want to estimate the extreme dependence structure of F py 1 , y 2 |xq in presence of contamination coming from F c py 1 , y 2 |xq.

Note that the logistic and Cauchy models have already been considered in Escobar-Bach et al. ( 2017) with a similar model (1) as ours, but in a framework without covariates. This model, naturally extended to the regression context in the present paper, is also satisfied for these two conditional models, with βpx 0 q " 1 for Model 1 and βpx 0 q " 2 for Model 2. Additionally, we can also check that the Hölder-type conditions pH 1´t q are satisfied. Concerning the conditional marginal distribution functions in the two models, they are standard heavy-tailed distributions (see, e.g., [START_REF] Beirlant | Second-order refined peaks-over-threshold modelling for heavy-tailed distributions[END_REF], and satisfy our Assumption pF m q.

To compute our estimates q L k py 1 , y 2 |x 0 q, first, we have to transform the margins into approximate unit Pareto distributions using the kernel-type empirical distribution functions given in (7). To this aim, we need to choose a kernel K and to select the bandwidths c n for each of the margins. Since the kernel has almost no impact on the results, we use in ( 7) and also in our MDPD procedure, the same bi-quadratic function Kpxq :" 15 16 p1 ´x2 q 2 1l txPr´1,1su , which satisfies our Assumption pK 2 q. Concerning the bandwidth c n , a cross validation criterion, already used in an extreme value context by [START_REF] Daouia | Kernel estimators of extreme level curves[END_REF], is performed, where c n,j :" arg min

r c j PC n ÿ i"1 n ÿ k"1 " 1l ! Y pjq i ďY pjq k ) ´r F n,´i,j pY pjq k |X i q  2 , j " 1, 2,
where C is a grid of values of r c j and r F n,´i,j py|xq :"

ř n k"1,k‰i K r c j px ´Xk q1l tY pjq k ďyu ř n k"1,k‰i K r c j px ´Xk q .
The bandwidth h n is taken as h n " minpc n,1 , c n,2 q `k n ˘1{d | logrminpc n,1 , c n,2 qs| ´ξ , where ξd ą q, in order to satisfy the condition:

n c h d n k d | log c n | q nc d n Ñ 0,
coming from (9) in our main theorem.

In the minimization of the empirical density power divergence, we fix β at the value 1, i.e., the true value for Model 1 and a mis-specified value for Model 2. The parameter k which determines the threshold p U q Z 1´t pn{k|x 0 q is selected by an automated procedure based on minimizing the standard deviation of the estimates q L k py 1 , y 2 |x 0 q computed in a moving window over the range for k, see, e.g., [START_REF] Goegebeur | Bias-corrected estimation for conditional Paretotype distributions with random right censoring[END_REF]. In all the settings, C " R X ˆt0.05, 0.075, . . . , 0.3u, where R X is the range of the covariate X, and ξ " 1.1. Figures 1, 2 and 3 illustrate the boxplots of the estimates q L k py 1 , y 2 |x 0 q based on 500 samples of size n " 1 000 for py 1 , y 2 q P tpl{10, 1 ´l{10q, l " 1, . . . , 9u and for three values of the covariate: x 0 " 3, 5 and 9, in case of Model 1. The columns of the figures represent the two fractions of contamination: ε " 0% (left) and ε " 10% (right), and the rows the three values of α, namely, from the top to the bottom, α " 0.1, 0.5 and 1. Figures 4, 5 and 6 are constructed similarly for Model 2 and the three covariate values: x 0 " 0.2, 0.5 and 0.8. Each time, the true function Lpy 1 , y 2 |x 0 q is computed at the same positions tpl{10, 1 ´l{10q, l " 1, . . . , 9u and connected with a blue line. Based on these simulations, we can draw the following conclusions:

• Overall, our robust estimator performs quite well, but of course, the results depend on the model, the covariate position and the fraction of contamination. In Model 1, Lpy 1 , y 2 |x 0 q depends on the covariate, but the marginal distributions do not. On the contrary, for Model 2, Lpy 1 , y 2 |x 0 q does not depend on the value of x 0 but the marginal distributions do;

• For all models, when ε " 0, the best results are obtained when α " 0.1. This result was expected since this value is close to 0, the value which leads to the maximum likelihood estimator, which is efficient (but not robust). On the contrary, in case of contamination, increasing α is crucial to get more robustness, the central box remaining closer to the true value for large values of α compared to α " 0.1;

• For Model 1, the contamination on the axes pulls slightly the estimates up, whereas, on the contrary, for Model 2, the diagonal contamination pulls the estimates a bit down, as expected. The estimation results are good for all covariate positions but exhibit more variability at py 1 , y 2 q close to p1{2, 1{2q.

• We have also considered data with 20% contamination but for such a high percentage of contamination the estimation procedure did not perform well anymore.

INSERT FIGURES 1-6 HERE

Real data analysis

In this section we illustrate the robust estimator for Lpy 1 , y 2 |x 0 q on the Ames housing dataset [START_REF] De Cock | Ames, Iowa: Alternative to the Boston Housing Data as an End of Semester Regression Project[END_REF], which is publicly available at https://www.kaggle.com/c/house-prices-advanced-regression-techniques. This dataset contains information on the sale of individual residential property in Ames, Iowa, from 2006 to 2010. The dataset has n " 2930 observations on a large number of variables (23 nominal, 23 ordinal, 14 discrete, and 20 continuous) involved in assessing home values.

We estimate the conditional stable tail dependence function of the variables sale price (Y p1q ) and above grade living area in square feet (Y p2q ) conditional on the original construction year of the property (X). When estimating a residential property's market value, living area is an important element to consider, since a bigger property will positively impact its valuation. Indeed, many buyers look at the sales price divided by the square footage of a property, which is a usual indicator of the value of a property. This thus motivates the study of the measure of dependence between the sale price of a residential property and the above grade living area in square feet. In Figure 7 we show the scatterplot of sale price versus above grade living area.

The scatterplot shows overall a positive association between the two variables, though there are also some observations that are atypical for the dependence structure and hence may disturb the estimation of the extremal dependence. One knows that properties that are newer often appraise at a higher value. Indeed, the fact that some parts of the property, like the plumbing, the electrical installations, and the roof are newer can generate savings for a buyer. For example, if a roof has a 20-year warranty, then that is money an owner will save over the next two decades, compared to an older home that may need a roof replaced in just a few years. This is illustrated in Figure 9 where the sale prices tend to be larger in case of an original construction year in 2004 (right panel) compared to 1946 (left panel) for a given living area. Thus, the original construction year is an important covariate which should be taken into account when estimating the measure of dependence between Y p1q and Y p2q .

INSERT FIGURE 7 HERE

We estimate Lpy 1 , y 2 |x 0 q in a robust way with the proposed local minimum density power divergence method. The estimation is implemented with the bi-quadratic kernel function and the same cross validation criterion for c 1 and c 2 as described in the simulation section. Also the bandwidth h n is here determined by cross validation. In Figure 8 we show the estimates q L k py 1 , y 2 |x 0 q with py 1 , y 2 q P pl{20, 1 ´l{20q, l " 1, . . . , 19, for α " 0.1 (blue), 0.5 (black) and 1 (green), for the years of original construction 1946 (left) and 2004 (right). For the construction year 1946 we see that the non robust estimates, corresponding to α " 0.1, are in line with the robust estimates obtained with α " 0.5 and α " 1, which indicates that the data used for the estimation did not contain disturbing observations. For the construction year 2004, the non-robust estimates with α " 0.1 are somehow different from those obtained with α " 0.5 and α " 1, where the latter two are similar, which indicates potential outliers in the data used for the estimation. This is confirmed by the scatterplots of the data used in the local estimation, given in Figure 9. Indeed, for construction year 2004, there are two observations with an above grade living area greater than 4000 square feet (which is very high) for a corresponding sale price lower than 200 000 dollars (which is not credible). These two observations are clearly outliers, and since they are far away from the main cloud, that are atypical for the dependence structure, while the scatterplot of the data used for estimation at x 0 " 1946 does not indicate outliers. For x 0 " 2004, these two outlying observations were removed and the estimates for Lpy 1 , y 2 |x 0 q were calculated again. The result of this is shown in Figure 10, where we see that the estimates obtained with the three values of α are now closer together, as expected.

INSERT FIGURES 8-10 HERE

Concluding remarks

In this paper we introduced a robust nonparametric estimator for the stable tail dependence function when next to the variables of main interest, Y p1q and Y p2q , there is also a random covariate, X. The work proposed here provides a series of interesting open questions which will lead to further investigations, among them:

• Outlier detection. To reach the goal of robustness, we adapted the idea of MDPD estimation to our context. As illustrated in the original [START_REF] Basu | Robust and efficient estimation by minimizing a density power divergence[END_REF] paper, in the density power divergence criterion, the estimating equations consist generally of likelihood score functions with a relative-to-the model down-weighting for outlying observations. Thus, if an observation is unusual relative to the proposed model then its contribution to the estimating equations gets less weight and as such its influence on the estimation results becomes dampened. However, although the methodology avoids that possible isolated outliers can completely disturb the estimation, it does not allow to identify which observations are the outliers in the dataset. In other words, the takeaway message of our paper is to compute our estimator for several values of α, among them some small α like α " 0.1. If there is almost no difference in the estimates, that means that there are no outliers in the dataset and in that case maximum likelihood method should be used since it is efficient. On the contrary, if there are differences, that means that we have to take care because of the presence of outliers. If we want to know which observation is an outlier, a heuristic approach would be to draw the scatterplot of the data and to try to visualize the observations which seem to be far away from the main cloud. Then, we could remove them one after one, and compute again our estimator for different values of α. If the estimates become this time close to each other, that means that the observations removed were indeed outliers. This is the strategy used in our real dataset for the year 2004. As an alternative to this heuristic approach we could investigate in future research outlier detection on the basis of the empirical influence function of a robust estimator for tail dependence (like the one introduced in the present paper). This was pursued in [START_REF] Hubert | Detecting influential data points for the Hill estimator in Pareto-type distributions[END_REF] in the context of identifying influential observations for the Hill estimator in univariate extreme value statistics.

• Change-points. In some applications, we are faced with events that can cause structural changes in the underlying model. The statistical analysis for detecting such changes is referred to as change-point analysis. It has been recently considered in the multivariate extreme value framework, see, e.g., [START_REF] De Carvalho | Tracking change-points in multivariate extreme[END_REF] or [START_REF] Drees | Statistical inference on a changing extreme value dependence structure[END_REF]. Since the traditional methods for identifying change-points can struggle with the presence of outliers, robust methods should be developed in that context, based on, e.g., the MDPD method used in this paper. A starting point for this new topic of research might be the recent paper of [START_REF] Song | Sequential change point test in the presence of outliers: the density power divergence based approach[END_REF] to be adapted to the context of extreme values.

• Parametric models. In the present paper we have applied the MDPD estimation method locally in order to obtain a nonparametric estimator for tail dependence. The MDPD method can also be used for fitting completely parametric extreme value models to data.

In the nonparametric approach one lets 'the data speak for themselves' and as such we get a preview of the extreme dependence structure. This nonparametric estimate could also be useful to evaluate the fit of parametric models.

• Other types of divergences. In the paper we used the density power divergence method of [START_REF] Basu | Robust and efficient estimation by minimizing a density power divergence[END_REF] to obtain a robust estimate. The basic idea of the density power divergence is to introduce a density power weight in the estimation procedure. This idea is also at the basis of other robust methods like those based on the γ divergence [START_REF] Fujisawa | Robust parameter estimation with a small bias against heavy contamination[END_REF] and the β divergence [START_REF] Mihoko | Robust blind source separation by beta divergence[END_REF]. The development of estimation procedures for extreme value problems based on the latter types of divergences, and a comparison of their performance with that of estimators based on the density power divergence is a topic of future research.

Appendix

The minimization of the empirical density power divergence p ∆ α,1´t pδ 1´t |x 0 q is based on its derivative. Direct computations show that all the terms appearing in this derivative have the following form S n,1´t ps|x 0 q :" 1 k

n ÿ i"1 K hn px 0 ´Xi q ˜Z1´t,i p U Z 1´t pn{k|x 0 q ¸s 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu
for s ă 0.

Assuming F Z 1´t py|x 0 q is strictly increasing in y, we can rewrite this main statistic as follows:

S n,1´t ps|x 0 q " 1 k n ÿ i"1 K hn px 0 ´Xi q # 1 `ż Z 1´t,i p U Z 1´t pn{k|x 0 q s u s´1 p U s Z 1´t pn{k|x 0 q du + 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu " 1 k n ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu `1 k n ÿ i"1 K hn px 0 ´Xi q # ż Z 1´t,i p U Z 1´t pn{k|x 0 q s u s´1 p U s Z 1´t pn{k|x 0 q du + 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu " 1 k n ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu `ż 8 p U Z 1´t pn{k|x 0 q 1 k n ÿ i"1 K hn px 0 ´Xi q s u s´1 p U s Z 1´t pn{k|x 0 q 1l tuăZ 1´t,i u du " 1 k n ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu `ż 8 p U Z 1´t pn{k|x 0 q 1 k n ÿ i"1 K hn px 0 ´Xi q s u s´1 p U s Z 1´t pn{k|x 0 q 1l tF Z 1´t pZ 1´t,i |x 0 qă k n n k F Z 1´t pu|x 0 qu du " 1 k n ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu `ż 1 0 1 k n ÿ i"1 K hn px 0 ´Xi qsz ´1´s 1l tF Z 1´t pZ 1´t,i |x 0 qă k n n k F Z 1´t pz ´1 p U Z 1´t pn{k|x 0 q|x 0 qu dz " T n,1´t ps n,1´t p1|x 0 q|x 0 q `ż 1 0
T n,1´t ps n,1´t pz|x 0 q|x 0 q s z ´1´s dz,

where

T n,1´t py|x 0 q :" 1 k n ÿ i"1 K hn px 0 ´Xi q1l tF Z 1´t pZ 1´t,i |x 0 qă k n yu , y P p0, T s, s n,1´t pz|x 0 q :" n k F Z 1´t ´z´1 p U Z 1´t pn{k|x 0 q ˇˇx 0 ¯.
Thus we start this appendix with some auxiliary results allowing us to study the statistic T n,1´t py|x 0 q and subsequently in Section 7.2 we establish the weak convergence of S n,1´t ps|x 0 q. Finally, in Section 7.3, Theorem 2.1 will be established. The proof of Theorem 3.1 from Section 3 is deferred to the online Supplementary Material.

Auxiliary results in case of known margins

First we establish the joint weak convergence of processes W n,1´t j :" t a kh d n rT n,1´t j py|x 0 q ýf X px 0 qs; y P p0, T su, j " 1, . . . , J.

Lemma 7.1 Assume pD 1´t j q and pH 1´t j q for j " 1, . . . , J, pD 0.5 q, pH 0.5 q, pK 1 q, x 0 P IntpS X q with f X px 0 q ą 0, and y Þ Ñ F Z 1´t j py|x 0 q, j " 1, . . . , J, are strictly increasing. Consider sequences k Ñ 8 and h n Ñ 0 as n Ñ 8 such that k{n Ñ 0, kh d n Ñ 8, h

ηε 1´t 1 ^¨¨¨^ηε 1´t J ^ηε 0.5 n log n k Ñ 0, a kh d n h η f X ^ηG 1´t 1 ^¨¨¨^η G 1´t J n
Ñ 0, and for j " 1, . . . , J,

a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q|h η C 1´t j n Ñ 0 and a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q|h ηε 1´t j n log n k Ñ 0.
Then, for n Ñ 8, we have pW n,1´t 1 , . . . , W n,1´t J q pW 1´t 1 , . . . , W 1´t J q, in J pp0, T sq, for any T ą 0.

Lemma 7.2 Under the assumptions of Lemma 7.1, for any sequence u

pjq n satisfying b kh d n ¨F Z 1´t j pU Z 1´t j pn{k|x 0 q|x 0 q F Z 1´t j pu pjq n |x 0 q ´1' Ñ c j P R,
as n Ñ 8, j " 1, . . . , J, we have

¨bnh d n F Z 1´t 1 pu p1q n |x 0 q ˜p F Z 1´t 1 pu p1q n |x 0 q F Z 1´t 1 pu p1q n |x 0 q ´1. . . b nh d n F Z 1´t J pu pJq n |x 0 q ˜p F Z 1´t J pu pJq n |x 0 q F Z 1´t J pu pJq n |x 0 q ´1¸‹ ‹ ‹ ‹ ‹ ‹ ‹ ' 1 f X px 0 q ¨W1´t 1 p1q . . . W 1´t J p1q ‹ '.
Lemma 7.3 Assume pD 1´t j q and pH 1´t j q for j " 1, . . . , J, pD 0.5 q, pH 0.5 q, pK 1 q, x 0 P IntpS X q with f X px 0 q ą 0, and y Þ Ñ F Z 1´t j py|x 0 q, j " 1, . . . , J, are strictly increasing. Consider sequences k Ñ 8 and h n Ñ 0 as n Ñ 8 such that k{n Ñ 0, kh d n Ñ 8, h

ηε 1´t 1 ^¨¨¨^ηε 1´t J ^ηε 0.5 n log n k Ñ 0, a kh d n h η f X ^ηG 1´t 1 ^¨¨¨^η G 1´t J n Ñ 0, a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q| Ñ 0, j " 1, . . . , J. Then, we have b kh d n ¨p U Z 1´t 1 pn{k|x 0 q U Z 1´t 1 pn{k|x 0 q ´1 . . . p U Z 1´t J pn{k|x 0 q U Z 1´t J pn{k|x 0 q ´1 ‹ ‹ ‹ ‹ ' 1 f X px 0 q ¨W1´t 1 p1q . . . W 1´t J p1q ‹ '.
7.2 Joint weak convergence of S n,1´t j ps j |x 0 q, j " 1, ..., M

We have now all the ingredients to state the joint weak convergence of S n,1´t j ps j |x 0 q, j " 1, . . . , M . Note that we allow for the possibility that t j " t j 1 for j ‰ j 1 , but of course the statistics S n,1´t j ps j |x 0 q, j " 1, . . . , M , must be different. This is due to the fact that, for a given value of t, the study of the MDPD estimator p δ n,1´t requires the joint convergence in distribution of several statistics S n,1´t ps|x 0 q, with different values of s.

Theorem 7.1 Under the conditions of Theorem 2.1, we have, for s 1 , . . . , s M ă 0,

b kh d n ¨Sn,1´t 1 ps 1 |x 0 q ´1 1´s 1 f X px 0 q . . . S n,1´t M ps M |x 0 q ´1 1´s M f X px 0 q ‹ ' ¨s1 ş 1 0 " W 1´t 1 pzq z ´W1´t 1 p1q ı z ´s1 dz . . . s M ş 1 0 " W 1´t M pzq z ´W1´t M p1q ı z ´sM dz ‹ ‹ ‹ ' .
To prove this Theorem 7.1, we start to establish the weak convergence of an individual statistic S n,1´t ps|x 0 q, properly normalized. We have the following decomposition

b kh d n ˆSn,1´t ps|x 0 q ´1 1 ´s f X px 0 q " ż 1 0 rW 1´t pzq ´W1´t p1qs s z ´1´s dz `"b kh d n rT n,1´t ps n,1´t p1|x 0 q|x 0 q ´sn,1´t p1|x 0 qf X px 0 qs ´W1´t ps n,1´t p1|x 0 qq * `tW 1´t ps n,1´t p1|x 0 qq ´W1´t p1qu `bkh d n ps n,1´t p1|x 0 q ´1q f X px 0 q `ż 1 0 " b kh d n rT n,1´t
ps n,1´t pz|x 0 q|x 0 q ´sn,1´t pz|x 0 qf X px 0 qs ´W1´t ps n,1´t pz|x 0 qq * s z ´1´s dz `ż 1 0 rW 1´t ps n,1´t pz|x 0 qq ´W1´t pzqs s z ´1´s dz (11)

`fX px 0 q b kh d n ż 1 0
rs n,1´t pz|x 0 q ´zs s z ´1´s dz ":

ż 1 0 rW 1´t pzq ´W1´t p1qs s z ´1´s dz `6 ÿ i"1 T i,k . (12)
We study the terms separately. Clearly, using Lemma 5.2 from [START_REF] Goegebeur | A Weissman-type estimator of the conditional marginal expected shortfall[END_REF] we have that for n large, with arbitrary large probability,

|T 1,k | ď sup yPp0,2s ˇˇˇb kh d n rT n,1´t py|x 0 q ´yf X px 0 qs ´W1´t pyq ˇˇˇ, (13)
and

|T 4,k | ď sup yPp0,2s
ˇˇˇb kh d n rT n,1´t py|x 0 q ´yf X px 0 qs ´W1´t pyq ˇˇˇˇˇˇż

1 0 s z ´1´s dz ˇˇˇ, (14)
and hence, by Lemma 7.1 combined with the Skorohod construction we obtain T 1,k " o P p1q and T 4,k " o P p1q.

Using again Lemma 5.2 in [START_REF] Goegebeur | A Weissman-type estimator of the conditional marginal expected shortfall[END_REF] with continuity, we have

|T 2,k | " o P p1q. ( 15
)
Concerning T 3,k , we can use the following decomposition:

T 3,k " b kh d n « F Z 1´t p p U Z 1´t pn{k|x 0 q|x 0 q F Z 1´t pU Z 1´t pn{k|x 0 q|x 0 q ´1ff f X px 0 q " b kh d n » - ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 ´1fi fl 1 `δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q f X px 0 q `bkh d n « 1 `δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q ´1ff f X px 0 q " b kh d n » - ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 ´1fi fl 1 `δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q f X px 0 q `bkh d n δ 1´t pU Z 1´t pn{k|x 0 q|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q $ & % » - δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q δ 1´t pU Z 1´t pn{k|x 0 q|x 0 q ´˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q fi fl `» - ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q ´1fi fl , .
-

.

": ´bkh d n « p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ´1ff f X px 0 qp1 `oP p1qq `bkh d n δ 1´t pU Z 1´t pn{k|x 0 q|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q T p1q 3,k .
By Proposition B.1.10 in de [START_REF] De Haan | Extreme value theory[END_REF], for n large, with arbitrary large probability, we have for ε, ξ ą 0

|T p1q 3,k | ď ε ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q˘ξ `˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q `1. (16)
In the above, the notation a ˘' means a ' if a ě 1 and a ´' if a ă 1. This implies by Lemma 7.3 and our conditions that

T 3,k ´W1´t p1q. (17)
Concerning now T 5,k , we have for any δ P p0, 1q small

|T 5,k | ď ż δ 0 |W 1´t ps n,
Finally, concerning T 6,k , we have

T 6,k " f X px 0 q b kh d n ż 1 0 « F Z 1´t pz ´1 p U Z 1´t pn{k|x 0 q|x 0 q F Z 1´t pU Z 1´t pn{k|x 0 q|x 0 q ´zff s z ´1´s dz " f X px 0 q b kh d n $ & % ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 ´1, . - s ż 1 0 z ´sdz `fX px 0 q b kh d n ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 ż 1 0 ˜1 `δ1´t pz ´1 p U Z 1´t pn{k|x 0 q|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q ´1¸s z ´sdz ": ´fX px 0 q s 1 ´s b kh d n ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ´1¸p 1 `oP p1qq `fX px 0 q b kh d n ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 δ 1´t pU Z 1´t pn{k|x 0 q|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q T p1q 6,k with |T p1q 6,k | ď |s| ż 1 0 ˇˇˇˇˇδ 1´t pz ´1 p U Z 1´t pn{k|x 0 q|x 0 q δ 1´t pU Z 1´t pn{k|x 0 q|x 0 q ´˜z ´1 p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q ˇˇˇˇˇz ´sdz `|s| ż 1 0 ˇˇˇˇˇ˜z ´1 p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q ´1ˇˇˇˇˇˇz
´sdz " O P p1q, using arguments similar to those for T p1q 3,k . Consequently, using again Lemma 7.3, we deduce that

T 6,k ´s 1 ´s W 1´t p1q. (19)
Combining decomposition ( 12) with ( 13)-( 19), the proof of the marginal weak convergence of S n,1´t ps|x 0 q, properly normalized, is achieved. The joint weak convergence of p a kh d n rS n,1´t j ps j |x 0 q ´fX px 0 q{p1 ´sj qs, j " 1, . . . , M q follows from Lemmas 7.1 and 7.3, respectively.

Proof of Theorem 2.1

Again we first consider the case of a single estimator p by (16). This implies that b kh d n ´p L k py 1 , y 2 |x 0 q ´Lpy 1 , y 2 |x 0 q

b kh d n ´p L k py 1 , y 2 |x 0 q ´Lpy 1 , y 2 |x 0 q ¯" ´y1 b kh d n ´p G 1´t,k px 0 q ´G1´t px 0 q " ´y1 b kh d n ˜k n p U Z 1´t pn{k|x 0 q 1 `p δ n,1´t ´G1´t px 0 q " ´y1 G 1´t px 0 q b kh d n ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q 1 `δ1´t pU Z 1´t pn{k|x 0 q|x 0 q 1 `p δ n,1´t ´1" ´y1 G 1´t px 0 q b kh d n ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ´1y 1 G 1´t px 0 q b kh d n ´p δ n,1´t ´δ1´t pU Z 1´t pn{k|x 0 q|x 0 q ¯1 1 `p δ n,1´t `y1 G 1´t px 0 q p δ n,1´t ´δ1´t pU Z 1´t pn{k|x 0 q|x 0 q 1 `p δ n,1´t b kh d n ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ´1¸. Now remark that b kh d n ˇˇδ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q ´δ1´t pU Z 1´t pn{k|x 0 q|x 0 q ˇ" b kh d n ˇˇδ 1´t pU Z 1´t pn{k|x 0 q|x 0 q ˇˇˇˇˇˇˇδ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q δ 1´t pU Z 1´t pn{k|x 0 q|x 0 q ´1ˇˇˇˇ" o P p1q,
" ´y1 G 1´t px 0 q b kh d n ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ´1y 1 G 1´t px 0 q b kh d n ´p δ n,1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ¯`o P p1q.
Using the fact that

b kh d n ¨p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ´1 p δ n,1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ‹ ‹ ' ¨W1´tp1q f X px0q c ´2α ş 1 0 " W1´tpzq z ´W1´t p1q ı z 2α dz ´p1 `βqp2α `βq ş 1 0 " W1´tpzq z ´W1´t p1q ı z 2α`β dz ¯‹ b kh d n ´p L k py 1 , y 2 |x 0 q ´Lpy 1 , y 2 |x 0 q ¯ ´y1 G 1´t px 0 q W 1´t p1q f X px 0 q `y1 G 1´t px 0 qc " 2α ż 1 0 " W 1´t pzq z ´W1´t p1q  z 2α dz ´p1 `βqp2α `βq ż 1 0 " W 1´t pzq z ´W1´t p1q  z 2α`β dz * .
Now, concerning the finite dimensional convergence, it follows from Lemma 7.3 combined with the following theorem which states the joint behavior of the MDPD estimator p δ n,1´t j , j " 1, ..., J, and whose proof is deferred to the online Supplementary Material: Theorem 7.2 Under the conditions of Theorem 2.1, with probability tending to one, there exists sequences of solutions p p δ n,1´t j q ně1 , j " 1, . . . , J, to the MDPD estimating equations such that

¨p δ n,1´t 1 ´δ1´t 1 p p U Z 1´t 1 pn{k|x 0 q|x 0 q . . . p δ n,1´t J ´δ1´t J p p U Z 1´t J pn{k|x 0 q|x 0 q ‹ ' P Ñ 0.
Moreover, for the consistent solution sequences one has that

b kh d n ¨p δ n,1´t1 ´δ1´t1 p p U Z1´t 1 pn{k|x 0 q|x 0 q . . . p δ n,1´t J ´δ1´t J p p U Z1´t J pn{k|x 0 q|x 0 q ‹ ' c ¨2α ş 1 0 " W1´t 1 pzq z ´W1´t1 p1q ı z 2α dz ´p1 `βqp2α `βq ş 1 0 " W1´t 1 pzq z ´W1´t1 p1q ı z 2α`β dz . . . 2α ş 1 0 " W1´t J pzq z ´W1´t J p1q ı z 2α dz ´p1 `βqp2α `βq ş 1 0 " W1´t J pzq z ´W1´t J p1q ı z 2α`β dz ‹ ‹ ‹ '
, where c is defined in Theorem 2.1. Figure 10: Ames housing dataset: estimates of Lpy 1 , y 2 |x 0 q with x 0 " 2004 for α " 0.1 (blue), 0.5 (black) and 1 (green) after removal of the two outlying observations.

SUPPLEMENTARY MATERIAL

8 Proof of Lemma 7.1

From Corollary 5.1 in [START_REF] Goegebeur | A Weissman-type estimator of the conditional marginal expected shortfall[END_REF], with a slight adjustment for the left continuity of T n,1´t j py|x 0 q, we have that W n,1´t j W 1´t j in pp0, T sq, for j " 1, . . . , J. To obtain the joint convergence of pW n,1´t 1 , . . . , W n,1´t J q we need to verify the finite dimensional convergence and joint tightness. We start by proving the finite dimensional convergence. Let

T n :" » - - T n,1´t 1 py 1 |x 0 q . . .
T n,1´t J py J |x 0 q fi ffi fl and T :"

» - - y 1 f X px 0 q . . . y J f X px 0 q fi ffi fl . Since b kh d n pT n ´T q " b kh d n pT n ´EpT n qq `bkh d n pEpT n q ´T q,
and a kh d n pEpT n q ´T q Ñ 0, by Proposition 5.1 in [START_REF] Goegebeur | A Weissman-type estimator of the conditional marginal expected shortfall[END_REF] and our assumptions, it is sufficient to study the weak convergence of a kh d n pT n ´EpT n qq. To this aim we use the Cramér-Wold device, and show that

Λ n :" ψ T b kh d n pT n ´EpT n qq N p0, ψ T Σψq
for all ψ :" pψ 1 , . . . , ψ J q T P R J , where the elements of Σ are as in the statement of Lemma 7.1. By a straightforward rearrangement of terms we have

Λ n " n ÿ i"1 c h d n k J ÿ j"1 ψ j " K hn px 0 ´Xi q1l tF Z 1´t j pZ1´t j ,i|x0qă k n yj u ´E ˆKhn px 0 ´Xq1l tF Z 1´t j pZ1´t j |x0qă k n yj u ˙ ": n ÿ i"1 V i,n .
Since V 1,n , . . . , V n,n are independent and identically distributed random variables, we have VarpΛ n q " nVarpV 1,n q, and hence

VarpΛ n q " J ÿ j"1 J ÿ j 1 "1 ψ j ψ j 1 nh d n k C j,j 1 , (20)
where

C j,j 1 :" Cov ˆKhn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j u , K hn px 0 ´Xq1l tF Z 1´t j 1 pZ 1´t j 1 |x 0 qă k n y j 1 u ˙.
We have

E ˆK2 hn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j ,F Z 1´t j 1 pZ 1´t j 1 |x 0 qă k n y j 1 u " E ˆK2 hn px 0 ´XqP ˆF Z 1´t j pZ 1´t j |x 0 q ă k n y j , F Z 1´t j 1 pZ 1´t j 1 |x 0 q ă k n y j 1 ˇˇˇX ˙˙. Now let m p1q n :" max ˆUZ 1´t j ˆn ky j ˇˇˇx 0 ˙, U Z 1´t j 1 ˆn ky j 1 ˇˇˇx 0 ˙˙, m p2q n :" max ˆtj 1 ´tj U Z 1´t j ˆn ky j ˇˇˇx 0 ˙, t j 1 1 ´tj 1 U Z 1´t j 1 ˆn ky j 1 ˇˇˇx 0 ˙˙.
Then, for x P S X ,

P ˆF Z 1´t j pZ 1´t j |x 0 q ă k n y j , F Z 1´t j 1 pZ 1´t j 1 |x 0 q ă k n y j 1 ˇˇˇX " x ˙" P ´Y p1q ą m p1q n , Y p2q ą m p2q n ˇˇX " x ¯.
We have

U Z 1´t py|x 0 q " G 1´t px 0 qyr1 `a1´t py|x 0 qs,
where a 1´t p.|x 0 q is regularly varying with index ´βpx 0 q, and hence

m p1q n " n k max ˜G1´t j px 0 q y j " 1 `a1´t j ˆn ky j ˇˇˇx 0 ˙ , G 1´t j 1 px 0 q y j 1 " 1 `a1´t j 1 ˆn ky j 1 ˇˇˇx 0 ˙" : n k r m p1q n , m p2q n " n k max ˜tj 1 ´tj G 1´t j px 0 q y j " 1 `a1´t j ˆn ky j ˇˇˇx 0 ˙ , t j 1 1 ´tj 1 G 1´t j 1 px 0 q y j 1 " 1 `a1´t j 1 ˆn ky j 1 ˇˇˇx 0 ˙" : n k r m p2q n .
Using our model (1) for pY p1q , Y p2q q, gives then

P ´Y p1q ą m p1q n , Y p2q ą m p2q n ˇˇX " x " k n p r m p1q n q ´d1 pxq p r m p2q n q ´d2 pxq g ˜1 r m p1q n , 1 r m p2q n ˇˇˇˇx ¸«1 `δ ˜k n 1 r m p1q n , k n 1 r m p2q n ˇˇˇˇx ¸ff ,
and thus

E ˆK2 hn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j ,F Z 1´t j 1 pZ 1´t j 1 |x 0 qă k n y j 1 u " k nh d n ż S K K 2 pvqp r m p1q n q ´d1 px 0 ´hnvq p r m p2q n q ´d2 px 0 ´hnvq g ˜1 r m p1q n , 1 r m p2q n ˇˇˇˇx 0 ´hn v «1 `δ ˜k n 1 r m p1q n , k n 1 r m p2q n ˇˇˇˇx 0 ´hn v ¸ff f X px 0 ´hn vqdv. Write δ ˜k n 1 r m p1q n , k n 1 r m p2q n ˇˇˇˇx 0 ´hn v ¸" δ ˆk n , k n ˇˇˇx 0 ˙δ `k n , k n ˇˇx 0 ´hn v δ `k n , k n ˇˇx 0 ˘«ξ ˜1 r m p1q n , 1 r m p2q n ˇˇˇˇx 0 ´hn v δ ´k n 1 r m p1q n , k n 1 r m p2q n ˇˇx 0 ´hn v δ `k n , k n ˇˇx 0 ´hn v ˘´ξ ˜1 r m p1q n , 1 r m p2q n ˇˇˇˇx 0 ´hn v ¸fi ffi
fl .

By our model (1), the assumptions of Lemma 7.1, and the fact that δ 0.5 p.|xq satisfies pD 0.5 q and pH 0.5 q, one has that

δ ˜k n 1 r m p1q n , k n 1 r m p2q n ˇˇˇˇx 0 ´hn v ¸" δ ˆk n , k n ˇˇˇx 0 ˙«ξ ˜1 r m p1q n , 1 r m p2q n ˇˇˇˇx 0 ´hn v ¸`op1q ff ,
where the op1q term is uniform in v P S K . Hence, by Lebesgue's dominated convergence theorem

E ˆK2 hn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j ,F Z 1´t j 1 pZ 1´t j 1 |x 0 qă k n y j 1 u " k nh d n $ & % }K} 2 2 f X px 0 q ˜max ˜G1´t j px 0 q y j , G 1´t j 1 px 0 q y j 1 ¸¸´d 1 px 0 q ˆ˜max ˜tj 1 ´tj G 1´t j px 0 q y j , t j 1 1 ´tj 1 G 1´t j 1 px 0 q y j 1 ¸¸´d 2 px 0 q ˆg ¨1 max ˆG1´t j px 0 q y j , G 1´t j 1 px 0 q y j 1 ˙, 1 max ˆtj 1´t j G 1´t j px 0 q y j , t j 1 1´t j 1 G 1´t j 1 px 0 q y j 1 ˙ˇˇˇˇˇˇˇx 0 ‹ ‹ ' `op1q , / / . / / -
.

By taking into account that

E ˆKhn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j u ˙" k n E `Tn,1´t j py j |x 0 q ˘" k n y j f X px 0 qp1 `op1qq,
where the last step follows from Proposition 5.1 in [START_REF] Goegebeur | A Weissman-type estimator of the conditional marginal expected shortfall[END_REF], we have that

nh d n k C j,j 1 Ñ }K} 2 2 f X px 0 q ˜max ˜G1´t j px 0 q y j , G 1´t j 1 px 0 q y j 1 ¸¸´d 1 px 0 q ˆ˜max ˜tj 1 ´tj G 1´t j px 0 q y j , t j 1 1 ´tj 1 G 1´t j 1 px 0 q y j 1 ¸¸´d 2 px 0 q ˆg ¨1 max ˆG1´t j px 0 q y j , G 1´t j 1 px 0 q y j 1 ˙, 1 max ˆtj 1´t j G 1´t j px 0 q y j , t j 1 1´t j 1 G 1´t j 1 px 0 q y j 1 ˙ˇˇˇˇˇˇˇx 0 ‹ ‹ ' ,
which shows the convergence of the variance in (20).

To ensure the convergence in distribution of Λ n to a normal random variable, we have to verify the Lyapounov condition for triangular arrays of random variables (Billingsley, 1995, p. 362).

In the present context this simplifies to showing that nE|V 1,n | 3 Ñ 0. We have

E|V 1,n | 3 ď ˆhd n k ˙3{2 $ & % E » - ˜J ÿ j"1 |ψ j |K hn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j u ¸3fi fl `3E » - ˜J ÿ j"1 |ψ j |K hn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j u ¸2fi fl ˆE ˜J ÿ j"1 |ψ j |K hn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j u 4 « E ˜J ÿ j"1 |ψ j |K hn px 0 ´Xq1l tF Z 1´t j pZ 1´t j |x 0 qă k n y j u ¸ff3
, .

-.

With arguments similar to the ones used above when studying the covariance terms, one finds that E|V 1,n | 3 " Op1{pn a kh d n qq, and hence nE|V 1,n | 3 Ñ 0. This establishes the finite dimensional convergence. The joint tightness follows from the individual tightness (similarly to Lemma 1 in [START_REF] Bai | Multivariate limit theorems in the context of long-range dependence[END_REF].

9 Proof of Lemma 7.2 Note that p F Z 1´t j pu pjq n |x 0 q F Z 1´t j pu pjq n |x 0 q " 1 1 n ř n i"1 K hn px 0 ´Xi q 1 n ř n i"1 K hn px 0 ´Xi q1l tZ 1´t j ,i ąu pjq n u F Z 1´t j pu pjq n |x 0 q " 1 1 n ř n i"1 K hn px 0 ´Xi q F Z 1´t j ´UZ 1´t j pn{k|x 0 q|x 0 F Z 1´t j pu pjq n |x 0 q T n,1´t j ´n k F Z 1´t j ´upjq n |x 0 ¯ˇˇx 0 ¯.
Then, by using Cramér-Wold device, we have for all ψ :" pψ 1 , . . . , ψ J q T P R J J ÿ j"1

ψ j c nh d n F Z 1´t j pu pjq n |x 0 q » - p F Z 1´t J pu pjq n |x 0 q F Z 1´t J pu pjq n |x 0 q ´1fi fl " J ÿ j"1 ψ j c nh d n F Z 1´t j pu pjq n |x 0 q ˆ» - 1 1 n ř n i"1 K hn px 0 ´Xi q F Z 1´t j ´UZ 1´t j pn{k|x 0 q|x 0 F Z 1´t j pu pjq n |x 0 q T n,1´t j ´n k F Z 1´t j ´upjq n |x 0 ¯ˇˇx 0 ¯´1 fi fl " J ÿ j"1 ψ j a kh d n " T n,1´t j ´n k F Z 1´t j ´upjq n |x 0 ¯ˇˇx 0 ¯´n k F Z 1´t j pu pjq n |x 0 qf X px 0 q ı ´W1´t j ´n k F Z 1´t j ´upjq n |x 0 ¯1 n ř n i"1 K hn px 0 ´Xi q `J ÿ j"1 ψ j b kh d n ¨g f f f e F Z 1´t j ´UZ 1´t j pn{k|x 0 q|x 0 F Z 1´t j pu pjq n |x 0 q ´1‹ ' ˆTn,1´t j ´n k F Z 1´t j ´upjq n |x 0 ¯ˇˇx 0 ¯´n k F Z 1´t j pu pjq n |x 0 qf X px 0 q 1 n ř n i"1 K hn px 0 ´Xi q `J ÿ j"1 ψ j W 1´t j ´n k F Z 1´t j ´upjq n |x 0 ¯¯´W 1´t j p1q 1 n ř n i"1 K hn px 0 ´Xi q `J ÿ j"1 ψ j W 1´t j p1q 1 n ř n i"1 K hn px 0 ´Xi q ´c k n a nh d n " 1 n ř n i"1 K hn px 0 ´Xi q ´fX px 0 q ‰ 1 n ř n i"1 K hn px 0 ´Xi q J ÿ j"1 ψ j g f f f e F Z 1´t j pu pjq n |x 0 q F Z 1´t j ´UZ 1´t j pn{k|x 0 q|x 0 ¯ J ÿ j"1 ψ j W 1´t j p1q f X px 0 q ,
under our assumptions combined with Lemma 7.1 and the Skorohod construction.

10 Proof of Lemma 7.3

First, remark that

P ˜bkh d n ˜p U Z 1´t j pn{k|x 0 q
U Z 1´t j pn{k|x 0 q ´1¸ď z j , @j " 1, . . . , J "

P ¨cnh d n F Z 1´t j pa pjq n |x 0 q ¨p F Z 1´t j pa pjq n |x 0 q F Z 1´t j pa pjq n |x 0 q ´1' ď c nh d n F Z 1´t j pa pjq n |x 0 q ¨F Z 1´t j pU Z 1´t j p n k |x 0 q|x 0 q F Z 1´t j pa pjq n |x 0 q ´1' , @j " 1, . . . , J ',
where a pjq n :" U Z 1´t j pn{k|x 0 qp1 `zj { a kh d n q. Then, using Lemma 7.2, it remains to show that, for all j " 1, . . . , J, we have

b kh d n ¨F Z 1´t j pU Z 1´t j pn{k|x 0 q|x 0 q F Z 1´t j pa pjq n |x 0 q ´1' ÝÑ c j P R.
This convergence holds since, under our assumptions, we have

b kh d n ¨F Z 1´t j pU Z 1´t j pn{k|x 0 q|x 0 q F Z 1´t j pa pjq n |x 0 q ´1' " z j 1 `δ1´t j pU Z 1´t j pn{k|x 0 q|x 0 q 1 `δ1´t j pa pjq n |x 0 q ´bkh d n δ 1´t j pU Z 1´t j pn{k|x 0 q|x 0 q 1 `δ1´t j pa pjq n |x 0 q ˜δ1´t j pa pjq n |x 0 q δ 1´t j pU Z 1´t j pn{k|x 0 q|x 0 q ´1Ý Ñ z j . Now, since c nh d n F Z 1´t j pa pjq n |x 0 q ¨F Z 1´t j pU Z 1´t j pn{k|x 0 q|x 0 q F Z 1´t j pa pjq n |x 0 q ´1' ÝÑ z j ,
by Lemma 7.2 and continuity we deduce that

P ˜bkh d n ˜p U Z 1´t j pn{k|x 0 q
U Z 1´t j pn{k|x 0 q ´1¸ď z j , @j " 1, . . . , J ¸ÝÑ P ˆW1´t j p1q f X px 0 q ď z j , @j " 1, . . . , J

˙.

This achieves the proof of Lemma 7.3.

11 MDPD calculations and proof of Theorem 7.2

Derivatives

We need to compute the two first derivatives of the empirical divergence. Direct computations yield

d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t " p1 `αq # ż 8 1 h α py; δ 1´t , βq d dδ 1´t hpy; δ 1´t , βqdy 1 k n ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu ´1 k n ÿ i"1 K hn px 0 ´Xi qh α´1 ˜Z1´t,i p U Z 1´t pn{k|x 0 q ; δ 1´t , β ḑ dδ 1´t h ˜Z1´t,i p U Z 1´t pn{k|x 0 q ; δ 1´t , β ¸1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 q u + " p1 `αq # 2αβ p1 `2αqp1 `2α `βq 1 k n ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu `1 k n ÿ i"1 K hn px 0 ´Xi q ˜Z1´t,i p U Z 1´t pn{k|x 0 q ¸´2α 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu ´p1 `βq 1 k n ÿ i"1 K hn px 0 ´Xi q ˜Z1´t,i p U Z 1´t pn{k|x 0 q ¸´2α´β 1l tZ 1´t,i ą p U Z
ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu ´α ´1 k n ÿ i"1 K hn px 0 ´Xi qh α´2 ˜Z1´t,i p U Z 1´t pn{k|x 0 q ; δ 1´t , β ḑ dδ 1´t h ˜Z1´t,i p U Z 1´t pn{k|x 0 q ; δ 1´t , β ¸¸2 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 q u ´1 k n ÿ i"1 K hn px 0 ´Xi qh α´1 ˜Z1´t,i p U Z 1´t pn{k|x 0 q ; δ 1´t , β ḑ2 dδ 2 1´t h ˜Z1´t,i p U Z 1´t pn{k|x 0 q ; δ 1´t , β ¸1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 q u + " p1 `αq # αβ 2 p´7 `β `4α 2 `2αβq p1 `2αqp1 `2α `βqp1 `2α `2βq 1 k n ÿ i"1 K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu ´1 `α k n ÿ i"1 K hn px 0 ´Xi q ˜Z1´t,i p U Z 1´t pn{k|x 0 q ¸´2α 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu `2p1 `αqp1 `βq 1 k n ÿ i"1 K hn px 0 ´Xi q ˜Z1´t,i p U Z 1´t pn{k|x 0 q ¸´2α´β 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu ´r1 `α `2β `2αβ `αβ 2 ´β2 s ˆ1 k n ÿ i"1 K hn px 0 ´Xi q ˜Z1´t,i p U Z 1´t pn{k|x 0 q ¸´2α´2β 1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu
`OP pδ 1´t qu .

Asymptotic behavior of the derivatives

In this section we work under the assumptions of Theorem 2.1.

b kh d n d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t | δ 1´t "δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q " p1 `αq " 2αβ p1 `2αqp1 `2α `βq ˆbkh d
n rs n,1´t p1|x 0 q ´1s f X px 0 q `W1´t p1q 2αβ p1 `2αqp1 `2α `βq ˆbkh d n rT n,1´t ps n,1´t p1|x 0 q|x 0 q ´sn,1´t p1|x 0 qf X px 0 qs ´W1´t ps n,1´t p1|x 0 qq 2αβ p1 `2αqp1 `2α `βq rW 1´t ps n,1´t p1|x 0 qq ´W1´t p1qs `bkh

d n " S n,1´t p´2α|x 0 q ´1 1 `2α f X px 0 q  ´p1 `βq b kh d n " S n,1´t p´2α ´β|x 0 q ´1 1 `2α `β f X px 0 q  `OP ˆbkh d n δ 1´t ´p U Z 1´t pn{k|x 0 q ˇˇx 0 ¯˙* .
Now, remark that by Proposition B.1.10 in de [START_REF] De Haan | Extreme value theory[END_REF], for n large, with arbitrary large probability, we have for ε, ξ ą 0

b kh d n ˇˇδ 1´t ´p U Z 1´t pn{k|x 0 q ˇˇx 0 ¯ˇ" b kh d n ˇˇδ 1´t ´UZ 1´t pn{k|x 0 q ˇˇx 0 ¯ˇˇˇˇˇˇˇˇδ 1´t ´p U Z 1´t pn{k|x 0 q ˇˇx 0 δ1´t ´UZ 1´t pn{k|x 0 q ˇˇx 0 ¯ˇˇˇˇď b kh d n ˇˇδ 1´t ´UZ 1´t pn{k|x 0 q ˇˇx 0 ¯ˇˇ» -ε ˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q˘ξ `˜p U Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q fi fl " o P p1q.
Using Lemma 5.2 from Goegebeur et al. (2021), Lemma 7.1 and (16) in the proof of Theorem 7.1, both combined with the Skorohod construction, and Theorem 7.1, yields

b kh d n d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t | δ 1´t "δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q ´2αp1 `αq ż 1 0 " W 1´t pzq z ´W1´t p1q  z 2α dz `p1 `βqp2α `βqp1 `αq ż 1 0 " W 1´t pzq z ´W1´t p1q  z 2α`β dz. (21)
Also

d 2 p ∆ α,1´t pδ 1´t |x 0 q dδ 2 1´t | δ 1´t "δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q " β 2 p1 `αqp1 `β `4α 2 `2αβq p1 `2αqp1 `2α `βqp1 `2α `2βq f X px 0 q
`αβ 2 p1 `αqp´7 `β `4α 2 `2αβq p1 `2αqp1 `2α `βqp1 `2α `2βq rT n,1´t ps n,1´t p1|x 0 q|x 0 q ´sn,1´t p1|x 0 qf X px 0 qs `αβ 2 p1 `αqp´7 `β `4α 2 `2αβq p1 `2αqp1 `2α `βqp1 `2α `2βq rs n,1´t p1|x 0 q ´1s f X px 0 q ´p1 `αq 2 " S n,1´t p´2α|x 0 q ´fX px 0 q 1 `2α  `2p1 `αq 2 p1 `βq

" S n,1´t p´2α ´β|x 0 q ´fX px 0 q 1 `2α `β  ´p1 `αq " 1 `α `2β `2αβ `αβ 2 ´β2 ‰ ˆ"S n,1´t p´2α ´2β|x 0 q ´fX px 0 q 1 `2α `2β  `OP ´δ1´t ´p U Z 1´t pn{k|x 0 q ˇˇx 0 ¯" β 2 p1 `αq " 1 `β `4α 2 `2αβ ‰ p1 `2αqp1 `2α `βqp1 `2α `2βq f X px 0 q `oP p1q. ( 22
)
11.3 Proof of Theorem 7.2

We limit the proof to deriving the asymptotic properties of a single MDPD estimator p δ n,1´t . The joint asymptotic behavior of p p δ n,1´t 1 , . . . , p δ n,1´t J q follows then from Theorem 7.1.

(i) Existence and consistency

The idea is to adjust the arguments used to prove existence and consistency of solutions of the likelihood estimating equations, see, for instance, Theorems 3.7 and 5.1 in Chapter 6 of [START_REF] Lehmann | Theory of Point Estimation[END_REF], to the MDPD framework. To this aim, we start to prove that, for any r ą 0 sufficiently small, we have

P ´p ∆ α,1´t ´δ1´t ´p U Z 1´t pn{k|x 0 q ˇˇx 0 ¯ˇˇx 0 ¯ă p ∆ α,1´t pδ 1´t |x 0 q for δ 1´t " δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q ˘rÝ Ñ 1, (23)
as n Ñ 8. By applying a Taylor series expansion, we have

p ∆ α,1´t pδ 1´t |x 0 q ´p ∆ α,1´t pδ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q|x 0 q " ´δ1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ¯d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t ˇˇˇˇδ 1´t "δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q `1 2 ´δ1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ¯2 d 2 p ∆ α,1´t pδ 1´t |x 0 q dδ 2 1´t ˇˇˇˇδ 1´t "δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q `1 6 ´δ1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ¯3 d 3 p ∆ α,1´t pδ 1´t |x 0 q dδ 3 1´t ˇˇˇˇδ 1´t " r δ 1´t ": T 1,n `T2,n `T3,n ,
where r δ 1´t is an intermediate value between δ 1´t and δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q. According to (21), we have |T 1,n | ď r 3 , whereas according to ( 22), there exists c ą 0 such that T 2,n ą cr 2 with probability tending to 1. Additionally, tedious computations allow us to show that sup

δ 1´t Prδ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q´r;δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q`rs ˇˇˇˇd 3 p ∆ α,1´t pδ 1´t |x 0 q dδ 3 1´t ˇˇˇˇă M,
with arbitrary large probability, from which we can deduce that |T 3,n | ď M r 3 {6 with probability tending to 1. Combining all these bounds, we deduce that with probability tending to 1, we have T 1,n `T2,n `T3,n ą cr 2 ´p1 `M {6qr 3 , which yields (23).

To complete the proof, we adjust the line of argumentation of Theorem 3.7 in Chapter 6 of [START_REF] Lehmann | Theory of Point Estimation[END_REF]. To this aim, for any r ą 0 sufficiently small, we define

S n prq :" ! p ∆ α,1´t ´δ1´t ´p U Z 1´t pn{k|x 0 q ˇˇx 0 ¯ˇˇx 0 ¯ă p ∆ α,1´t pδ 1´t |x 0 q for δ 1´t " δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q ˘r) .
For v P S n prq, since p ∆ α,1´t pδ 1´t |x 0 q is differentiable with respect to δ 1´t , there exists r δ n,1´t P ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ´r; δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q `rw here p ∆ α,1´t pδ 1´t |x 0 q achieves a local minimum. According to (23), we have PpS n prqq Ñ 1 for any small enough r, and hence there exists a sequence r n Ó 0 such that PpS n pr n qq Ñ 1 as n Ñ 8. Now, let p δ n,1´t :" r δ n,1´t if v P S n pr n q and arbitrary otherwise. Since v P S n pr n q implies

d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t ˇˇˇˇδ 1´t " p δ n,1´t
" 0, we have

P ¨d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t ˇˇˇˇδ 1´t " p δ n,1´t
" 0 'ě PpS n pr n qq ÝÑ 1, as n Ñ 8, which establishes the existence part.

Concerning now the consistency part, note that for any r ą 0 and n large enough such that r n ď r, we have

P ´ˇˇp δ n,1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ˇˇă r ¯ě P ´ˇˇp δ n,1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ˇˇă r n ě PpS n pr n qq ÝÑ 1,
as n Ñ 8, whence the consistency of the estimator sequence.

(ii) Asymptotic normality By definition we have

d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t ˇˇδ 1´t " p δ n,1´t
" 0. Thus a Taylor series expansion around δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q combined with the boundedness of the third derivative of p ∆ α,1´t pδ 1´t |x 0 q with respect to δ 1´t leads to 0

" d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t ˇˇˇˇδ 1´t "δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q `d2 p ∆ α,1´t pδ 1´t |x 0 q dδ 2 1´t ˇˇˇˇδ 1´t "δ 1´t p p U Z 1´t pn{k|x 0 q|x 0 q ´p δ n,1´t ´δ1´t p p U Z 1´t pn{k|x 0 q|x 0 q ¯p1 `oP p1qq from which we deduce that b kh d n ´p δ n,1´t ´δ1´t p p U Z1´t pn{k|x 0 q|x 0 q " ´¨d 2 p ∆ α,1´t pδ 1´t |x 0 q dδ 2 1´t ˇˇˇˇδ 1´t"δ1´tp p U Z 1´t pn{k|x0q|x0q '´1 ˆbkh d n d p ∆ α,1´t pδ 1´t |x 0 q dδ 1´t ˇˇˇˇδ 1´t"δ1´tp p U Z 1´t pn{k|x0q|x0q p1 `oP p1qq c " 2α ż 1 0 " W 1´t pzq z ´W1´t p1q  z 2α dz ´p1 `βqp2α `βq ż 1 0 " W 1´t pzq z ´W1´t p1q  z 2α`β dz * .
Proposition 12.1 Assume that there exists b ą 0 such that f X pxq ě b, @x P S X Ă R d , f X is bounded, pD 1´t j q, pH 1´t j q for j " 1, . . . , J, pD 0.5 q, pH 0.5 q, pK 2 q, pF m q hold, and that y Þ Ñ F Z 1´t j py|x 0 q, j " 1, . . . , J, are strictly increasing at x 0 P IntpS X q non-empty. Consider sequences k Ñ 8, h n Ñ 0 and c n Ñ 0 as n Ñ 8, such that k{n Ñ 0, kh d n Ñ 8, h

ηε 1´t 1 ^¨¨¨^ηε 1´t J ^ηε 0.5 n log n k Ñ 0, a kh d n h η f X ^ηG 1´t 1 ^¨¨¨^η G 1´t J n Ñ 0, a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q|h η C 1´t j n Ñ 0 and a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q|h ηε 1´t j n
log n k Ñ 0 for j " 1, . . . , J. Under conditions ( 8) and ( 9), we have

´| W n,1´t 1 , . . . , | W n,1´t J ¯ pW 1´t 1 , . . . , W 1´t J q ,
in J pp0, T sq, for any T ą 0.

Then, the next ingredients required in order to proof the weak convergence of q S n,1´t ps|x 0 q are lemmas similar to Lemma 7.3 from the present paper and Lemma 5.2 in [START_REF] Goegebeur | A Weissman-type estimator of the conditional marginal expected shortfall[END_REF], but this time for the variable q Z 1´t instead of Z 1´t . This is the aim of Lemmas 12.2 and 12.3 below for which we need first to show the weak convergence of p F q Z 1´t p.|x 0 q correctly normalized.

Lemma 12.1 Under the assumptions of Proposition 12.1, for any sequence u

pjq n satisfying b kh d n ¨F Z 1´t j pU Z 1´t j pn{k|x 0 q|x 0 q F Z 1´t j pu pjq n |x 0 q ´1' Ñ c j P R,
as n Ñ 8, j " 1, . . . , J, we have

¨bnh d n F Z 1´t 1 pu p1q n |x 0 q ˜p F q Z 1´t 1 pu p1q n |x 0 q F Z 1´t 1 pu p1q n |x 0 q ´1. . . b nh d n F Z 1´t J pu pJq n |x 0 q ˜p F q Z 1´t J pu pJq n |x 0 q F Z 1´t J pu pJq n |x 0 q ´1¸‹ ‹ ‹ ‹ ‹ ‹ ‹ ' 1 f X px 0 q ¨W1´t 1 p1q . . . W 1´t J p1q ‹ '.
Lemma 12.2 Assume that there exists b ą 0 such that f X pxq ě b, @x P S X Ă R d , f X is bounded, pD 1´t j q, pH 1´t j q, for j " 1, . . . , J, pD 0.5 q, pH 0.5 q, pK 2 q, pF m q hold, and that y Þ Ñ F Z 1´t j py|x 0 q, j " 1, . . . , J, are strictly increasing at x 0 P IntpS X q non-empty. Consider sequences k Ñ 8, 8) and ( 9), we have

h n Ñ 0 and c n Ñ 0 as n Ñ 8, such that k{n Ñ 0, kh d n Ñ 8, h ηε 1´t 1 ^¨¨¨^ηε 1´t J ^ηε 0.5 n log n k Ñ 0, a kh d n h η f X ^ηG 1´t 1 ^¨¨¨^η G 1´t J n Ñ 0, a kh d n |δ 1´t j pU Z 1´t j p n k |x 0 q|x 0 q| Ñ 0, j " 1, . . . , J. Under conditions (
b kh d n ¨p U q Z 1´t 1 pn{k|x 0 q U Z 1´t 1 pn{k|x 0 q ´1 . . . p U q Z 1´t J pn{k|x 0 q U Z 1´t J pn{k|x 0 q ´1 ‹ ‹ ‹ ‹ ‹ ' 1 f X px 0 q ¨W1´t 1 p1q . . . W 1´t J p1q ‹ '.
From Lemma 12.2, we can show now the uniform convergence in probability of q s n,1´t pz|x 0 q towards z for any z P p0, T s.

Lemma 12.3 Under the assumptions of Lemma 12.2, for any T ą 0, we have sup zPp0,T s |q s n,1´t pz|x 0 q ´z| " o P p1q.

12.1 Proof of the auxiliary results in case of unknown margins Proof of Proposition 12.1. Firstly, we consider the weak convergence of a single process t a kh d n p q T n,1´t py|x 0 q ´yf X px 0 qq; y P p0, T su, where for simplicity of notation we have ignored the index j from t. We use the decomposition b kh d n ´q T n,1´t py|x 0 q ´yf X px 0 q " b kh d n pT n,1´t py|x 0 q ´yf X px 0 qq `bkh d n ´q T n,1´t py|x 0 q ´Tn,1´t py|x 0 q ´E " q T n,1´t py|x 0 q ´Tn,1´t py|x 0 q ıb kh d n E " q T n,1´t py|x 0 q ´Tn,1´t py|x 0 q ı ":

3 ÿ i"1 Q i,1´t py|x 0 q. (24)
According to Lemma 7.1, we have

Q 1,1´t py|x 0 q W 1´t pyq, (25)
To this aim, we will make use of empirical process theory with changing function classes, see for instance [START_REF] Van Der Vaart | Weak convergence and empirical processes, with applications to statistics[END_REF]. We start by introducing some notation. Let P be the distribution measure of pY p1q , Y p2q , Xq, and denote the expected value under P as P f :" ş f dP for any real-valued measurable function f : R 2 ˆRd Ñ R. For a function class F, define now the covering number N pF, L 2 pQq, τ q as the minimal number of L 2 pQq-balls of radius τ needed to cover the class of functions F and the uniform entropy integral as Jpδ, F, L 2 q :"

ż δ 0 c log sup QPQ N pF, L 2 pQq, τ }F } Q,2 q dτ,
where Q is the set of all probability measures Q for which 0 ă }F } 2 Q,2 :" ş F 2 dQ ă 8 and F is an envelope function for the class F.

Let

I n :" tg y,δ,n : y P p0, T s, δ P Hu where H :" δ " pδ 1 , δ 2 q; δ : R ˆR ˆSX Ñ R 2 ( , and

g y,δ,n pv 1 , v 2 , uq :" c nh d n k K hn px 0 ´uqq y,δ,n pv 1 , v 2 , uq, q y,δ,n pv 1 , v 2 , uq :" 1l tF Z 1´t pZ δ pv 1 ,v 2 ,uq|x 0 qă k n yu , with Z δ pv 1 , v 2 , uq :" min ˆ1 |1 ´δ1 pv 1 , v 2 , uq| , 1 ´t t 1 |1 ´δ2 pv 1 , v 2 , uq| ˙.
For convenience, denote δ n :" ´p F n,1 , p F n,2 ¯and δ 0 :" pF 1 , F 2 q. According to Lemma 3.1 in Escobar-Bach et al. (2018a), if r n :" maxp a | log c n | q {nc d n , c η n q, we have r ´1 n |δ n ´δ0 | converges in probability towards the null function H 0 :" t0u in H, endowed with the norm }δ} H :" }δ 1 } 8 `}δ 2 } 8 for any δ P H. We want to apply Theorem 2.3 in van der [START_REF] Van Der Vaart | Empirical processes indexed by estimated functions. Asymptotics: Particles, Processes and Inverse Problems[END_REF]. To this aim, we consider the class E n py, bq :" tg y,δ 0 `rnδ,n ´gy,δ 0 ,n : δ P H, }δ} H ď bu , K hn px 0 ´uq1l tUZ 1´t pn{pkyq|x 0 qPrminpZ δ 0 `rnδ pv 1 ,v 2 ,uq,Z δ 0 pv 1 ,v 2 ,uqq,maxpZ δ 0 `rnδ pv 1 ,v 2 ,uq,Z δ 0 pv 1 ,v ď G n py, bqpv 1 , v 2 , uq.

We need to show Assertion 1: sup yPp0,T s ? nP G n py, b n q ÝÑ 0 for every b n Ñ 0 and Assertion 2: Let G n be an envelope for tG n py, bq : y P p0, T su, we have A) P G Thus @ε ą 0, we have

P G 2 n tG n ě ε ? nu ď 1 n α{2 ε α E ´G2`α n ď C ε α n k 1 pkh d n q α{2
ˆ"E Combining (25), ( 26) and ( 27) we have shown that b kh d n ´q T n,1´t py|x 0 q ´yf X px 0 q ¯ W 1´t pyq, in pp0, T sq.

Secondly, to obtain the joint weak convergence of p | W n,1´t 1 , . . . , | W n,1´t J q we need to establish the finite dimensional weak convergence together with tightness. The finite dimensional weak convergence can be shown by the Cramér-Wold device, as in the proof of Lemma 7.1 along with decomposition (24). The joint tightness follows from the individual tightness (similarly to Lemma 1 in [START_REF] Bai | Multivariate limit theorems in the context of long-range dependence[END_REF].

Proof of Lemma 12.1. This proof is omitted since it is similar to the proof of Lemma 7.2, by using Proposition 12.1 instead of Lemma 7.1 combining with the Skorohod construction (keeping the same notation).

Proof of Lemma 12.2. Again, this proof is omitted since it is similar to the proof of Lemma 7.3, by using Lemma 12.1 instead of Lemma 7.2.

Proof of Lemma 12.3. By definition, we have |q s n,1´t pz|x 0 q ´z| " ˇˇˇˇˇF Z 1´t ´z´1 p U q Z 1´t pn{k|x 0 q ˇˇx 0 F Z 1´t `UZ 1´t pn{k|x 0 q|x 0 ˘´z ˇˇˇˇ" ˇˇˇˇˇz ˜p U q Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 1 `δ1´t ´z´1 p U q Z 1´t pn{k|x 0 q ˇˇx 0 1 `δ1´t `UZ 1´t pn{k|x 0 q|x 0 ˘´z ˇˇˇˇď

z ˇˇˇˇˇ˜p U q Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 ´1ˇˇˇˇˇˇˇˇˇˇˇˇ1 `δ1´t ´z´1 p U q Z 1´t pn{k|x 0 q ˇˇx 0 1 `δ1´t `UZ 1´t pn{k|x 0 q|x 0 ˘ˇˇˇˇž ˇˇˇˇˇ1 `δ1´t ´z´1 p U q Z 1´t pn{k|x 0 q ˇˇx 0 1 `δ1´t `UZ 1´t pn{k|x 0 q|x 0 ˘´1 ˇˇˇˇď z ˇˇˇˇˇ˜p U q Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´1 ´1ˇˇˇˇˇˇˇˇˇˇˇˇ1 `δ1´t ´z´1 p U q Z 1´t pn{k|x 0 q|x 0 1 `δ1´t `UZ 1´t pn{k|x 0 q|x 0 ˘ˇˇˇˇž ˇˇˇˇδ1´t `UZ 1´t pn{k|x 0 q|x 0 1 `δ1´t `UZ 1´t pn{k|x 0 q|x 0 ˘ˇˇˇˇ$ & % ˇˇˇˇˇ˜z ´1 p U q Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q ´1ˇˇˇˇˇδ 1´t ´z´1 p U q Z 1´t pn{k|x 0 q|x 0 δ1´t
`UZ 1´t pn{k|x 0 q|x 0 ˘´˜z ´1 p U q Z 1´t pn{k|x 0 q U Z 1´t pn{k|x 0 q ¸´βpx 0 q ˇˇˇˇˇ, .

-.

By Proposition B.1.10 in de [START_REF] De Haan | Extreme value theory[END_REF], for n large, with arbitrary large proba-

Figure 1 :

 1 Figure1: Logistic model, x 0 " 3: no contamination (left) and 10% axis contamination (right), and α " 0.1 (first row), α " 0.5 (second row) and α " 1 (third row).

Figure 2 :

 2 Figure2: Logistic model, x 0 " 5: no contamination (left) and 10% axis contamination (right), and α " 0.1 (first row), α " 0.5 (second row) and α " 1 (third row).

Figure 3 :

 3 Figure3: Logistic model, x 0 " 9: no contamination (left) and 10% axis contamination (right), and α " 0.1 (first row), α " 0.5 (second row) and α " 1 (third row).

Figure 4 :

 4 Figure4: Cauchy model, x 0 " 0.2: no contamination (left) and 10% diagonal contamination (right), and α " 0.1 (first row), α " 0.5 (second row) and α " 1 (third row).

Figure 5 :

 5 Figure5: Cauchy model, x 0 " 0.5: no contamination (left) and 10% diagonal contamination (right), and α " 0.1 (first row), α " 0.5 (second row) and α " 1 (third row).

Figure 6 :

 6 Figure6: Cauchy model, x 0 " 0.8: no contamination (left) and 10% diagonal contamination (right), and α " 0.1 (first row), α " 0.5 (second row) and α " 1 (third row).

Figure 7 :

 7 Figure 7: Ames housing dataset: scatterplot of sale price versus above grade living area.

Figure 8 :

 8 Figure8: Ames housing dataset: estimates of Lpy 1 , y 2 |x 0 q with x 0 " 1946 (left) and x 0 " 2004 (right) for α " 0.1 (blue), 0.5 (black) and 1 (green).

Figure 9 :

 9 Figure 9: Ames housing dataset: scatterplot of data used for local estimation at x 0 " 1946 (left) and x 0 " 2004 (right).

  in pp0, T sq. The next step consists in showing that sup yPp0,T s |Q 2,1´t py|x 0 q| " o P p1q.

  with envelope function given by G n py, bqpv 1 , v 2 , uq :" Pr1´F 1 pv 1 |uq´rnb;1´F 1 pv 1 |uq`rnbs * Pr t 1´t p1´F 2 pv 2 |uq´rnbq; t 1´t p1´F 2 pv 2 |uq`rnbqs * since |g y,δ 0 `rnδ,n pv 1 , v 2 , uq ´gy,δ 0 ,n pv 1 , v 2 , uq|" px 0 ´uq|q y,δ 0 `rnδ,n pv 1 , v 2 , uq ´qy,δ 0 ,n pv 1 , v 2 , uq| " px 0 ´uq ˇˇ1l tZ δ 0 `rnδ pv 1 ,v 2 ,uqąU Z 1´t pn{pkyq|x 0 qu ´1l tZ δ 0 pv 1 ,v 2 ,uqąU Z 1´t pn{pkyq|x 0 qu ˇď c nh d n k

n

  tG n ě ε ? nu Ñ 0, @ε ą 0; C) sup yPp0,T s P G 2 n py, bq Ñ 0; D) J pd n , tG n py, bq : y P p0, T su, L 2 q Ñ 0, @d n OE 0.Proof of Assertion 1. Remark that ? nP G n py, Pr1´F 1 pY p1q |Xq´rnbn,1´F 1 pY p1q |Xq`rnbns Pr t 1´t p1´F2pY p2q |Xq´rnbnq, t 1´t p1´F2pY p2q |Xq`rnbnqs n E rK hn px 0 ´Xqs .Thus Assertion 1 is satisfied as soon as r n n

  pn{pkT q|x 0 q `rn b  which tends again to 0 under our assumptions.Proof of Assertion 2 C). Following the lines of proof of Assertion 1Proof of Assertion 2 D). We follow the lines of proof of Theorem 2.2 in[START_REF] Escobar-Bach | Local estimation of the conditional stable tail dependence function[END_REF]. The class of functions on r0, 1su Ñ 1l ty 1 ď1´uďy 2 u , y 1 ă y 2 ( , is a V C-class.This allows us to prove that there exist positive constants C and V such that sup QPQ N ptG n py, bq :y P p0, T su , L 2 pQq, τ }G n } Q,2 q ď C ˆ1 τ ˙V ,from which Assertion 2 D) follows. This achieves the proof of (26px 0 ´Xq ´1l t1´F Z 1´t p qZ 1´t |x 0 qă k n yu ´1l t1´F Z 1´t pZ 1´t |x 0 qă k px 0 ´Xq ˇˇ1l t1´F Z 1´t p q Z 1´t |x 0 qă k n yu ´1l t1´F Z 1´t pZ 1´t |x 0 qă k n yu ˇˇı px 0 ´Xq ˇˇ1l t1´F Z 1´t pZ δn |x 0 qă k n yu ´1l t1´F Z 1´t pZ δ 0 |x 0 qă k n yu ˇˇı ď ? nP G n py, bq,for n large enough. This implies that sup py, bq ÝÑ 0 by Assertion 1 since it is clear that b n Ñ 0 can be replaced by any fixed value b as soon as

  1´t pz|x 0 qq ´W1´t pzq| |s| z ´1´s dz

	`ż 1	|W 1´t ps n,1´t pz|x 0 qq ´W1´t pzq| |s| z ´1´s dz
	δ
	#				+	ż δ
	ď |s|	sup	|W 1´t ps n,1´t pz|x 0 qq| `sup	|W 1´t pzq|	z ´1´s dz
	zPp0,δs	zPp0,δs			0
				ż 1
	`|s| sup	|W 1´t ps n,1´t pz|x 0 qq ´W1´t pzq|	z ´1´s dz
	zPpδ,1s		δ
	" o P p1q.

  1´t pn{k|x 0 qu

				`OP pδ 1´t qu
	and
	d 2	p ∆ α,1´t pδ 1´t |x 0 q
		dδ 2 1´t
	" p1 `αq	# α	ż 8 1	h α´1 py; δ 1´t , βq	ˆd dδ 1´t	hpy; δ 1´t , βq ˙2 dy	1 k	n ÿ i"1	K hn px 0 ´Xi q1l tZ 1´t,i ą p U Z 1´t pn{k|x 0 qu
			`ż 8 1	h α py; δ 1´t , βq	d 2 dδ 2 1´t	hpy; δ 1´t , βqdy	k 1	n

  Prminp|1´F 1 pv 1 |uq´rnδ 1 |,1´F 1 pv 1 |uqq,maxp|1´F 1 pv 1 |uq´rnδ 1 |,1´F 1 pv 1 |uqqs Pr t 1´t minp|1´F 2 pv 2 |uq´rnδ 2 |,1´F 2 pv 2 |uqq, t 1´t maxp|1´F 2 pv 2 |uq´rnδ 2 |,1´F 2 pv 2 |uqqs

	ď	c	nh d n k	K hn px 0 ´uq
		1 pn{pkyq|x 0 q `1l " ˆ"1l t U Z 1´t 1
				U Z 1´t	pn{pkyq|x 0 q
					2 ,uqqsu

*

, .

Acknowledgements

The authors sincerely thank the editor, associate editor and the referees for their helpful comments and suggestions that led to substantial improvement of the paper. The research of Armelle Guillou was supported by the French National Research Agency under the grant ANR-19-CE40-0013-01/ExtremReg project and an International Emerging Action (IEA-00179). Computation/simulation for the work described in this paper was supported by the DeIC National HPC Centre, SDU.

Auxiliary results in case of unknown margins

We introduce the same key statistic as S n,1´t ps|x 0 q but, this time, defined with Z 1´t replaced by q Z 1´t , that is q S n,1´t ps|x 0 q :" 1 k

Similarly as in (10) of the present paper, assuming that F Z 1´t py|x 0 q is strictly increasing in y, we can show that q S n,1´t ps|x 0 q " q T n,1´t pq s n,1´t p1|x 0 q|x 0 q `ż 1 0 q T n,1´t pq s n,1´t pz|x 0 q|x 0 q s z ´1´s dz, where q T n,1´t py|x 0 q :" 1 k

The aim of the next theorem is to show the joint weak convergence of q S n,1´t j ps j |x 0 q, j " 1, . . . , M .

Theorem 12.1 Under the same conditions as Theorem 3.1, we have for s

, where the processes W 1´t j , j " 1, . . . , M , are as in Theorem 2.1.

Note that this limiting process is the same as the one obtained in Theorem 7.1, with additional conditions needed to measure the discrepancy between the conditional distribution function F j py|xq and its empirical kernel version p F n,j py|xq, j " 1, 2, uniformly in px, yq.

To prove Theorem 12.1, we need, as a preliminary result, the weak convergence of statistics q T n,1´t py|x 0 q. This is the aim of the next proposition. Let

T n,1´t j py|x 0 q ´yf X px 0 q ¯; y P p0, T s * , for j " 1, . . . , J.

bility, we have for ε, ξ ą 0 and z ď T

Using Lemma 12.2, Lemma 12.3 follows.

12.2 Proof of Theorem 12.1.

This proof is omitted since it follows exactly the same lines of proof as those for Theorem 7.1 but with the auxiliary results in case of unknown margins.

12.3 Proof of Theorem 3.1.

This proof is omitted since it follows exactly the same lines of proof as those for Theorem 2.1 with the auxiliary results in case of unknown margins.