# Seasonal bias-correction of daily precipitation over France using a stitch model designed for robust extremes representation

Philippe Ear, Elena Di Bernardino, Thomas Laloë, Magali Troin, Adrien Lambert

<figure><img src="image_1.png" /></figure>

To cite this version:

Philippe Ear, Elena Di Bernardino, Thomas Laloë, Magali Troin, Adrien Lambert. Seasonal bias-correction of daily precipitation over France using a stitch model designed for robust extremes representation. 2025. hal-05004515

HAL Id: hal-05004515

https://hal.science/hal-05004515v1

Preprint submitted on 26 Mar 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Seasonal bias-correction of daily precipitation over
France using a stitch model designed for robust
extremes representation

Ear Philippe¹,²*†, Di Bernardino Elena¹†, Laloë Thomas¹†, Troin Magali²†, Lambert Adrien²†

¹Laboratoire Jean-Alexandre Dieudonné, Université Côte d'Azur, Nice, 06108, France.

²Hydroclimat, Aubagne, 13400, France.

*Corresponding author(s). E-mail(s): philippe.ear@hydroclimat.com;

†These authors contributed equally to this work.

Abstract

Highly resoluted and accurate daily precipitation data are required for impact
models to perform adequately and to correctly measure high-risk events' impact.
In order to produce such data, bias-correction is often needed. Most of those
statistical methods correct the probability distributions of daily precipitation
by modeling them using either empirical or parametric distributions. A recent
semi-parametric model based on a penalized Berk-Jones (BJ) statistical test
which allows for an automatic and personalized splicing of parametric and non-
parametric has been developed. This method, called Stitch-BJ model, was found
to be able to model daily precipitation correctly and showed interesting potential
in a bias-correction setting. In the present study, we will consolidate these results
by taking into account the seasonal properties of daily precipitation in an out-of-
sample context, and by considering dry days probabilities in our methodology. We
evaluate the performance of the Stitch-BJ method in this seasonal bias-correction
setting against more classical models such as the Gamma, Exponentiated Weibull
(ExpW), Extended Generalized Pareto (EGP) or empirical distributions. The
Stitch-BJ distribution was able to consistently perform as good or better than
all the other models over the validation set, including the empirical distribution,
which is often used due to its robustness.

**Keywords:** Bias correction; Extreme value theory; Goodness-of-Fit; Parametric distribution; Precipitation modeling
# 1 Introduction

Highly resoluted and robust daily precipitation data are key for hydrological models, impact models or to produce climate indicators (Casanueva et al. 2016; Shayeghi et al. 2024). Insufficiently resoluted models are not able to produce realistic extreme events due to the smoothing effect of models and them not being convection-permitting (Fosser et al. 2024). Such extreme events can often lead to dramatic catastrophes with important economic and human impacts if not correctly assessed (Alfieri and Thielen 2015; Sangati and Borga 2009). However, highly resoluted and robust data are quite scarce and often suffer from either short time series or small available areas. Most common datasets are issued from Global Circulation Models (GCM) which provide global coverage over a long period but at a coarse resolution (typically from 50km to 200km), which is not resoluted enough for impact studies or to produce realistic extremes (Henckes et al. 2018; Prein et al. 2016). Regional Climate Models (RCM) are able to bridge this resolution gap to some extent, with some RCM producing data with a resolution as high as 2x2km. However, data from GCM and RCM suffer from bias from multiple sources (Wang et al. 2014) which need to be corrected (Şan et al. 2023; Xu et al. 2015), and increasing a model resolution may not increase its skill in modeling extreme events compared to its lower resolution counterpart (Bador et al. 2020).

Multiple bias correction methods exist in the literature (Déqué 2007; Michelangeli et al. 2009) but the most used is often referred as the Quantile Mapping (QM) bias correction method, which is paired to empirical distributions and referred to as Empirical Quantile Mapping (EQM). However, other studies have used parametric distributions instead of the empirical one (Li et al. 2021; Naveau et al. 2016; Mamalakis et al. 2017) with promising results. The use of the empirical distribution is often justified by being distribution-free and its robustness. Variation to the EQM method while being fully empirical have also been explored (Velasquez et al. 2020; Byun and Hamlet 2024). Moreover, classical parametric distributions often fail to correctly model daily precipitation over large and highly resoluted study areas. This is especially true for study areas with a variety of climate zones, such as metropolitain France, containing mountainous regions, plains, as well as regions with a mediterranean climate (Joly et al. 2010; Strohmenger et al. 2024). Empirical models cannot however extrapolate values outside the sample's minimum and maximum, which prevent from producing new extremes that are likely to happen, especially in a climate change context (Boé et al. 2007; Déqué 2007; Byun and Hamlet 2024).

Let us remind the quantile mapping method :

$$ \hat{x}_{obs,fut} = F_{obs,ref}^{-1}(F_{mod,ref}(x_{mod,fut})). \quad (1) $$

where $F_{data,period}$ refers to the cumulative distribution function (CDF) of the considered dataset data for the period period. Here the datasets considered are the CERRA-Land reanalysis for obs and the ERA5-Land reanalaysis for mod, while the ref period is from the 01/01/1985 to 31/12/2009 and the fut period is from 01/01/2010 to 31/12/2020. More details on the datasets can be found in Section 2.1. In the past few years, new parametric and semi-parametric models have been studied to bridge
the gap between classical parametric distributions with their abilities to extrapolate extreme values, and the robustness and adaptability of empirical distributions. Such models may be fully parametric, with the presence of a stitch (merging of two distributional models at a given cutting point) or not (Naveau et al. 2016; Mamalakis et al. 2017; Derdour et al. 2022), or semi-parametric (Langousis et al. 2016; Holthuijzen et al. 2022; Trentini et al. 2023; Ear et al. 2025). In order to correctly evaluate the performance of a model, an out-of-sample validation is required, with a separation of training and validation data. Additionally, data is also often separated into months or seasons to remove intra-annual seasonality (Katiraie-Boroujerdy et al. 2020) and increase confidence on the data's stationarity. Finally, most studies only focus on correcting the wet days intensity, but only a few are also considering dry days probabilities (Vrac et al. 2016). The study from Ear et al. (2025) focused on introducing a semi-parametric model referred as to the Stitch-BJ model, but was not applied in a more applicative bias-correction context. In the present paper, we will study the bias-correction performance of the Stitch-BJ model, by comparing it to multiple parametric and non-parametric distributions and applying multiple modifications, while taking into account the following operational constraints:

1. The CERRA-Land and ERA5-Land datasets were separated into a training and validation period (01/01/1985 - 31/12/2009 and 01/01/2010 - 31/12/2020) as discussed in Section 2.1. This separation makes it possible to include the empirical distribution in the bias-correction performance comparison. As already remarked, in this study, CERRA-Land is used in Equation (1) as *obs* data and ERA5-Land as *mod* data;

2. A separation using meteorological seasons DJF (i.e. December-January-February), MAM, JJA and SON was used in order to take daily precipitation's seasonality into account and increase the time series' stationarity (see Section 2.1 for details);

3. Correction of dry days probability is included in the bias correction using the Singularity Stochastic Removal from Vrac et al. (2016) as described in Section 2.2.

We will then study the bias-correction performance using the quantile mapping in (1) paired with the Gamma, Exponentiated Weibull, Extended Generalized Pareto, empirical and Stitch-BJ models, all described in Section 2.3.

*Structure of the paper.* In Section 2.1, we present the datasets for daily precipitation over France used in this study and explain how and why they were divided for the bias-correction evaluation. Sections 2.2 and 2.3 describe the methods used for respectively the dry days' proportion correction and the distributions for wet days modeling. In Sections 3.1, 3.2 and in Appendix A, the seasonal bias correction results are presented where we used the CERRA-Land dataset to correct the ERA5-Land reanalysis with respect to several metrics. Section 3.3 compares the differences in performance depending on the season. An analysis of a selected location is presented in Section 3.4. Finally, the conclusion and discussion on the performances and limitations are included in Section 4.
# 2 Materials and Methods

## 2.1 Daily precipitation datasets

In this study, we consider the whole of metropolitan France, covering over 550000km².
Two reanalyses were used to evaluate the bias-correction performance of our models:

* The ERA5-Land (ERA5L) reanalysis (Muñoz-Sabater et al. 2021) with a resolution of 0.1° × 0.1°, spanning from 01/01/1950 to 31/12/2021;

* The Copernicus Regional Reanalysis for Europe (CERRA-Land, CERRAL) reanalysis (Verrelle et al. 2021) with a resolution of 5.5 × 5.5km and spanning from 01/01/1984 to 31/07/2021.

Both reanalysis were produced by the European Center for Medium-Range Weather Forecast (ECMWF) and are freely available through their Climate Data Store (see the *Data Availability Statement* section at the end of the present paper).

The ERA5L reanalysis is used as the model to be corrected (*mod* data as in (1)), while CERRAL is used as a gridded observational dataset (*obs* data as in (1)). The aim is to use bias-correction methods to match as closely ERA5L data to CERRAL data. The datasets used here are not free of bias (Pelosi 2023; Guo et al. 2024), they have been chosen in this study for illustration purposes but can be swapped for any other datasets that may have a different resolution or bias. Datasets were chosen due to their similar time-period availability as well as their similar spatial coverage for bias-correction performance evaluation.

The CERRAL data was interpolated into a regular longitude/latitude grid (approximately 0.05° × 0.05°) using the linear scattered interpolant (*scatteredInterpolant*) from MATLAB. The ERA5L data was interpolated to the same grid using the CloughTocher2D interpolator from Python's library SciPy.

We selected a common period for both reanalyses, starting from 01/01/1985 to 31/12/2020. Both reanalyses were separated using the meteorological seasons : DJF, MAM, JJA, SON as in Katiraie-Boroujerdy et al. (2020) and Gutiérrez et al. (2019). Discussion on the impact of subsampling on quantile mapping performance can be found in Reiter et al. (2018). The seasonality must be taken into account because of the intra-annual non-stationarity of daily precipitation. In our study area (Metropolitain France), summers see in general the least amount of precipitation, while autumn sees the most precipitation (Chaouche et al. 2010), as seen in Figure 1.

We performed two-sample t-tests between the means of each season for the period 01-12-1985 to 30/11/2009 to try to detect significative differences in the inter-seasonal mean. Our results in Table 1 show an important rejection rate of the test at a significance level α = 5% for all pairings, with a minimum of 35% of locations rejected (JJA-MAM on ERA5L data), and a maximum of 80% (SON-MAM on ERA5L data) of rejected locations. Such differences in means contradict the constant mean hypothesis for weak-stationarity in time series.

However, when testing the intra-seasonal means with a monthly division, we manage to greatly reduce the rejection rate for our Student's t-tests as shown in Table 2. In the table, each season is decomposed into its 3 months component e.g; DJF is decomposed into December, January and February, and are then referred to as respectively
<table><caption>**Table 1:** Proportion of rejected locations of the Student's t-test for seasonal means differences. Note that values in **bold** correspond to ERA5L results, while values in *italic* correspond to CERRAL results.</caption><thead><tr><th></th><th>DJF</th><th>MAM</th><th>JJA</th><th>SON</th></tr></thead><tbody><tr><th>DJF</th><td>x</td><td>0.54</td><td>0.52</td><td>0.50</td></tr><tr><th>MAM</th><td>0.50</td><td>x</td><td>0.35</td><td>0.80</td></tr><tr><th>JJA</th><td>0.60</td><td>0.49</td><td>x</td><td>0.55</td></tr><tr><th>SON</th><td>0.51</td><td>0.79</td><td>0.49</td><td>x</td></tr></tbody></table>

1, 2 and 3. The two samples Student's t-test is then applied on every pair of months for the given season. In the example, to get the proportion of rejection for December against February, one would need to look at the DJF table, on either the first column and last line for CERRAL data, or the last column and first line for ERA5L data. The maximum rejection proportion is at 33% (for March against April in the MAM subtable on ERA5L data) which is lower than the lowest rejection rate from Table 1. While a constant mean is not sufficient for stationarity, separation into seasons seems to be an acceptable precautionary step to take to improve the homogeneity of our data.

<table><caption>**Table 2:** Proportion of rejected locations of the Student's t-test for monthly means differences. Note that values in **bold** correspond to ERA5L results, while values in *italic* correspond to CERRAL results.</caption><thead><tr><th>DJF</th><th>1</th><th>2</th><th>3</th><th>MAM</th><th>1</th><th>2</th><th>3</th></tr></thead><tbody><tr><th>1</th><td>x</td><td>0.05</td><td>0.32</td><th>1</th><td>x</td><td>0.33</td><td>0.26</td></tr><tr><th>2</th><td>0.10</td><td>x</td><td>0.17</td><th>2</th><td>0.19</td><td>x</td><td>0.12</td></tr><tr><th>3</th><td>0.25</td><td>0.06</td><td>x</td><th>3</th><td>0.24</td><td>0.18</td><td>x</td></tr></tbody></table><table><thead><tr><th>JJA</th><th>1</th><th>2</th><th>3</th><th>SON</th><th>1</th><th>2</th><th>3</th></tr></thead><tbody><tr><th>1</th><td>x</td><td>0.11</td><td>0.11</td><th>1</th><td>x</td><td>0.14</td><td>0.30</td></tr><tr><th>2</th><td>0.17</td><td>x</td><td>0.03</td><th>2</th><td>0.13</td><td>x</td><td>0.15</td></tr><tr><th>3</th><td>0.10</td><td>0.08</td><td>x</td><th>3</th><td>0.25</td><td>0.18</td><td>x</td></tr></tbody></table>

Moreover, a separation into a training and validation set was done with the period from 01/01/1985 to 31/12/2009 (25 years) used for training and the period 01/01/2010 to 31/12/2020 (11 years) is used for validation. Figure 1 not only shows the intra-annual variability of daily precipitation, justifying the seasonal separation of our datasets but also the wet bias (in terms on total cumulated intensity here) the ERA5L dataset suffers compared to the CERRAL data. This can be noticed in panel (b) where for all seasons, the CERRAL dataset produced a median annual cumulated precipitation lower than ERA5L. More details on the wet or dry bias different datasets may suffer and how to correct them are discussed in Section 2.2.

The JJA season produced much less total precipitation (15% to 30% less) than the other seasons for both training and validation period. This further shows the need
<figure><img src="image_3.png" /><figcaption>**Fig. 1:** Yearly mean of seasonal cumulated precipitation for ERA5L and CERRAL over the training (a) and validation (b) periods.</figcaption></figure>

to seasonally divide our data as stationary time-series should have a constant mean throughout.

This separation was chosen to keep a sufficient amount of wet days in the dry season for the training period since the fitting ability of parametric distributions is greatly dependent on the number of available samples (Braunstein 1992). Using 25 years for the training period, some locations can receive as few as 104 (resp. 126) rainy days during the JJA season for respectively the ERA5L and CERRAL datasets.

## 2.2 On the correction of number of dry days and rain's probability

Most of the daily precipitation data coming from climate models (GCM, RCM or reanalysis) suffer from the drizzle effect (Chen et al. 2021; Gutowski et al. 2003; Argüeso et al. 2013) with too many low precipitations occurrences. The drizzle effect, which happens with numerical models producing a high number of rainy days of very small intensity, can have an important impact on the total precipitation amount over an extended period. The spatial and temporal distribution of precipitation for a given location is also affected, with daily precipitation-induced indicators such as the Cumulative Dry Days (CDD) being greatly impacted by this drizzle effect (Maraun 2013).

While both considered datasets share similarities in their model and forcing data, the differences in number of wet days and overall daily rain probability can be important. These differences can be due to both intrinsic differences from the model, but most likely from the difference in spatial resolution.

In Figure 2, we have the number of wet days (days with more than 1mm of rain) for each season, for both datasets in the training (a) and validation (b) period. There is a clear wet bias from the ERA5L model, with not only a higher median of the number of wet days but also a higher extreme number of wet days. The median differences of number of wet days between ERA5L and CERRAL is roughly 5 to 8 days per year, but can reach up to 35 days on some locations. This results in almost 1000 more rainy days during the training period for such locations.

Let us consider the difference in the probability of rain in Figure 3, where a positive value means a higher probability of rain in ERA5L compared to CERRAL for a given location.
<figure><img src="image_4.png" /><figcaption>**Fig. 2:** Yearly mean of number of wet days (¿1mm) for ERA5L and CERRAL per season for the training (a) and validation (b) periods.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Fig. 3:** Boxplot of differences in probability of rain (days with ¿1mm of rain) of ERA5L against CERRAL by season for the training (a) and validation (b) periods.</figcaption></figure>

In Figure 3, we can see a 5% higher probability of rain across all seasons for ERA5L, with some locations with 20 to 40% more rainy days in ERA5L than CERRAL. On the other hand, CERRAL locations have at most around 15% more rain than ERA5L, and this only happens rarely. Results are very similar for the training and validation period which seems to indicate that the evolution in rain probabilities evolved similarly in both datasets between the two periods.

If left as it is, not only the resulting bias-corrected time series will keep the wet bias, but the amount of rain might also be greatly superior to the target time series. The wet or dry bias must then be corrected to reduce the number of bias sources and produce bias-corrected data as trustworthy as possible. Multiple methods are available to correct differences in dry days probability with the threshold adaptation method (Schmidli et al. 2006; Lavaysse et al. 2012), positive correction method (Mao et al. 2015), or a more direct approach (Piani et al. 2010; Vigaud et al. 2013). A short intercomparison was done by Vrac et al. (2016) where a new method in the direct approach family, able to correct both wet and dry bias was also introduced as the Singularity Stochastic Removal method (SSR). The positive threshold method is the most widely used with various thresholds considered by the scientific community (Semenov et al. 1998; Ambrosino et al. 2014; Vaittinada Ayar et al. 2016; Bouvier et al. 2003; Ear et al. 2025). However, it is not able to take into account dry bias or the specificity of each location. Dry bias is much more difficult to correct than wet bias since the correction is two-fold: the intensity needs to be chosen as well as the temporal location of the correction. Randomly adding wet days to the time series
may break temporal correlation structure, precipitation being notably autocorrelated (Khan et al. 2007; Abbott et al. 2016).

The SSR method (Vrac et al. 2016) corrects both wet and dry bias, location-wise. However, for the dry bias correction, randomly chosen dry days are turned into wet days (thus the stochastic part) without any regard to the autocorrelation of daily precipitation data.

While the SSR method does not take autocorrelation into account, it allows for a simple method to correct dry days proportion, with no distinction between adding or removing wet days which could be the case for other methods.

The SSR method applies the following steps:

1. Select a threshold $th$ such that any value above $th$ is considered a wet day and any value below is considered a dry day. This can either be a common threshold or the minimum positive value of all datasets;

2. Set all days below $th$ (null days) to a random uniformly taken between 0 and $th$;

3. Perform the bias correction technique;

4. Set the bias-corrected data lower than $th$ to 0.

In the case where a distribution has been fitted to the wet days' data, the method will be applied as such, considering the quantile mapping for the bias correction method (see Equation (1)).

Let $F_{CERRAL,ref}$, $F_{ERA5L,ref}$ be the fitted models to respectively the reference observation data and reference model data. Let $x_{ERA5L,fut}$ be the to-be-corrected future model data. Let $\alpha_{CERRAL,ref}$ and $\alpha_{ERA5L,ref}$ be the probability of dry days respectively for the reference observed and reference modeled.

Quantile mapping methods are usually only applied to correct wet days distributions (Themeßl et al. 2012; Ajaaj et al. 2016; Lafon et al. 2013). To include the dry day proportion correction, let us define:

$$ \tilde{F}_{ERA5L,ref}(x) = \begin{cases} F_{ERA5L,ref}(x) (1 - \alpha_{ERA5L,ref}) + \alpha_{ERA5L,ref} & x \geq th, \\ \frac{\alpha_{ERA5L,ref}}{th} x & x < th, \end{cases} \quad (2) $$

and

$$ \tilde{F}_{CERRAL,ref}^{-1}(x) = \begin{cases} F_{CERRAL,ref}^{-1}\left(\frac{p-\alpha_{CERRAL,ref}}{1-\alpha_{CERRAL,ref}}\right) & p \geq \alpha_{CERRAL,ref}, \\ \frac{th}{\alpha_{CERRAL,ref}}p & p < \alpha_{CERRAL,ref}. \end{cases} \quad (3) $$

In Equations (2) and (3), any probability lower than $\alpha_{CERRAL,ref}$ will result in a value lower than $th$. The SSR method will both correct the dry days proportion with regards to the reference period, but also keep the evolution of dry days proportion between the reference and future period.
## 2.3 Parametric, semi-parametric and non-parametric models

Multiple distributional models were used in the bias correction procedure to compare their performances. We mostly used the same distributions as in (Ear et al. 2025) but we also added the empirical distribution for comparison. Distributions were fitted using the Maximum Likelihood Estimator (MLE) for all parametric models. Fitting was done location-wise, seasonally and using only wet days (¿1mm). When available, a 1mm location parameter was used, or a 1mm shift was applied to the fitting procedure for distributions without one. A left-censor has been used for the EGPD as in (Naveau et al. 2016), and the threshold has been determined through trials not shown here.

In the present study, the distributions used are the Gamma, Exponentiated Weibull, Extended Generalized Pareto, Stitch-BJ distribution and empirical model.

**Gamma distribution** The Gamma distribution is one of the most used parametric distribution for daily rainfall (Martinez-Villalobos and Neelin 2019; Husak et al. 2007). The distribution is defined by its cumulative distribution function (CDF) as follows:

$$F(x) = \frac{1}{\Gamma(x)} \gamma \left( k, \frac{x}{\theta} \right),$$

with $x > 0, k > 0$ and $\theta > 0$. $\Gamma$ is the gamma function and $\gamma$ is the lower incomplete gamma function.

**Exponentiated Weibull distribution.** The Exponentiated Weibull (ExpW) generalizes the Weibull distribution and while it has historically been used for failure rates and survival data modelisation, (Khan 2018; Mudholkar et al. 1996), the ExpW has also been used for precipitation modeling (Sharma et al. 2022). This distribution is defined by its CDF:

$$F(x) = \left[ 1 - e^{-(x/\lambda)^k} \right]^{\alpha},$$

with $x > 0, k > 0, \alpha > 0$ and $\lambda > 0$. With $\alpha = 1$, we retrieve the original Weibull distribution.

**Extended Generalized Pareto distribution.** Here we recall a family of distributions which has theoretical properties in line with Extreme Value Theory for both the lower and upper heavy-tail behavior (see Naveau et al. (2016)). The Extended Generalized Pareto distribution (EGP) has recently been used in multiple studies for rainfall modeling due to its flexibility and theoretical guarantees (Tencaliec et al. 2019; Rivoire et al. 2021; Haruna et al. 2023). We will only use the type 1 EGPD with the same 3mm left-censoring as used in Ear et al. (2025). The CDF of the type 1 EGPD is:

$$F(x) = \begin{cases} \left(1 - \left(1 + \frac{\xi x}{\sigma}\right)^{-\frac{1}{\xi}}\right)^{\kappa}, & \xi > 0, \\ \left(1 - e^{-x/\sigma}\right)^{\kappa}, & \xi = 0, \end{cases}$$

with $x > 0, \sigma > 0, \xi \ge 0$ and $\kappa > 0$. Notice that for $\xi = 0$, the model matches the previous ExpW model with $k = 1$.

**Empirical model.** The empirical distribution is often used in climate studies for its robustness (Enayati et al. 2021). The CDF can be written as follows:
$$F_n(x) = \frac{1}{n} \sum_{i=1}^{n} \mathbf{1}_{\{x_{(i)} \le x\}}$$

with $x_{(1)}, ..., x_{(n)}$ the order statistics of the considered sample, and n the sample size.

**Stitch-BJ model.** The Stitch-BJ model is a semi-parametric distributional model based on a penalized version of the Berk-Jones statistical test (PBJ) recently introduced in Ear et al. (2025).

The model is a flexible distribution which allows for an automatic stitch between the EGP, ExpW and empirical distributions using the PBJ. For more information, we refer to the original article where details can be found, as well as an analytical form of the resulting final CDF.

<figure><img src="image_6.png" /><figcaption>**Fig. 4:** Example of replacement indexes for the fitting of the Stitch-BJ model for a given location on ERA5L (a) and CERRAL (b) in the reference period (01/01/1985 to 31/12/2009).</figcaption></figure>

To easily understand the working mechanism of the Stitch-BJ, a simple example is given in Figure 4 where the quantile-quantile plot (QQplot) of the fitted EGP and ExpW distributions against the actual quantiles of respectively ERA5L (panel a) and CERRAL (panel b) are shown for a selected location. In both panels, the vertical line indicates an upper and lower rejection index (respectively in blue and green), named $i_u$ and $i_l$. These indexes are derived from the Penalized Berk-Jones statistical test and indicate the point from which the considered parametric distribution is deemed unfit.

On panel (a) of Figure 4, the final distribution can be expressed as :

$$F(x) = \begin{cases} F_{EGP}(x) & \text{if } x < F_{EGP}^{-1}(\frac{i_u}{n}) \\ F_n(x) & \text{else,} \end{cases} \quad (4)$$

with n the number of wet days for the considered time series. The upper tail of the EGP distribution has been replaced with the empirical distribution. Note that while the ExpW seems to be well fitted in panel (a) of Figure 4, errors on the upper tail can be up to 48mm compared to the reference data, which is why it was also rejected on the upper tail.
On panel (b), the lower tail of the EGP has been replaced with the ExpW distribution which resulted in a better fitting distribution. Using the PBJ test, the final stitch distribution can then be expressed as :

$$F(x) = \begin{cases} F_{ExpW}(x) & \text{if } x < F_{EGP}^{-1}(\frac{iu}{n}) \\ F_{EGP}(x) & \text{else.} \end{cases} \quad (5)$$

Note that Equations (4) and (5) are simple illustrations related to specific situations in panel (a) and (b) of Figure 4. Complete general analytical forms with monotonicity guarantees are available in Ear et al. (2025).

Using this methodology, we fitted the Stitch-BJ model to the daily precipitation data of ERA5L and CERRAL. From these fitted models, both a spatial analysis of the stitching type and a global analysis of the stitching proportion can be obtained in Figures 5 and 6.

<figure><img src="image_7.png" /><figcaption>Fig. 5: Replacement map for the Stitch-BJ fitted on the reference period of ERA5L (a, b) and CERRAL (c, d) for DJF and JJA.</figcaption></figure>

For both DJF and JJA season on ERA5L, the vast majority (resp 86.35% and 62.3%) of locations used a pure EGP distribution as shown in Figure 5 panels (a, b). Considering all fully parametric models, stitch between parametric models included, the proportion goes up to around 94% for both seasons. This shows that most locations where able to get a satisfying fit out of purely parametric distributions. The remaining locations used mainly an empirical lower tail replacement (3.6% and 4.75% of the total locations respectively), while almost no locations used a purely empirical distribution. Spatially, most of the locations using the empirical distribution for either the lower or
upper tail replacement are located in historically complex climate regions such as the Cévènnes, the Alps or the Corsica (Vautard et al. 2015; Emmanuel et al. 2017; Cortés-Hernández et al. 2024; Estermann et al. 2025). For the JJA season, these locations are a bit more spread out but still have a large concentration around the Alps and the Cévènnes, with some replacement along the Pyrénées.

Replacement maps for CERRAL in Figure 5 panel (c, d) are similar to those for ERA5L in panels (a, b) with most of the locations for both DJF and JJA seasons using pure EGP distributions (88.9% and 84% respectively). The use of an empirical lower or upper tail is even rarer with respectively 95.39% and 96.5% of locations using a purely parametric distribution. Use of the empirical distribution for the upper tail is mainly located around the Cévènnes for the DJF season, while for JJA locations seems to be spread out on all the study areas. No locations used a fully empirical distribution for both DJF and JJA season while fitted on the CERRAL dataset.

For locations where Stitch-BJ was used, one might want to know how much of the original distributions' the stitching replaced. Overall, the median replacement affected 5% of the distribution for the lower tail, and less than 1% for the upper tail. More details can be found in Figure 6 with boxplots of the percentage of replacement for the lower and upper tail per season, for respectively ERA5L (panels a and c) and CERRAL (panels b and d). The lower replacement index $i_l$ concerned a larger portion of the distribution compared to $i_u$ but its impact is much smaller due to the small intensity of rain affected. Overall, $i_l$ for CERRAL were smaller, both for the median and outliers, than for ERA5L, throughout all seasons. For the upper tail index $i_u$, similar observations can be made with CERRAL upper replacement portion being roughly equal to or lower than for ERA5L.

The Stitch-BJ model is heavily flexible, allowing to automatically replace a portion of the EGP distribution with other models, either parametric or empirical, using p-values results from the PBJ goodness-of-fit statistical tests. This allows for a better fit to precipitation distribution which can be extremely variable depending on the season and the location chosen.

In the next section, bias correction performance of all methods will be assessed using the CERRAL validation period.

## 3 Bias-correction results on the period 2010-2020

In this section, we compare the performance of the newly introduced Stitch-BJ against classical parametric distributions (EGP, ExpW, Gamma) and the empirical distribution (referred to as *emp*).

As explained in Section 2.1, both datasets ERA5L and CERRAL were first mapped on the same regular longitude/latitude grid (approximately 0.05° × 0.05°) and the following periods have been used for training and validation:

* Training period: 01/01/1985 to 31/12/2009;

* Validation period: 01/01/2010 to 31/12/2020;

along with a separation using the meteorological seasons. Only the results of the correction of the validation period are shown in this section as they are similar to the
<figure><img src="image_8.png" /><figcaption>Fig. 6: Boxplots of percentages of the lower (first row) and upper tail (second row) being replaced by the stitching algorithm for ERA5L (a, c) and CERRAL (b, d).</figcaption></figure>

ones from the previous study. In this section, the following metrics were computed on a pixel-by-pixel basis and we noted $q_i^n$ the quantile $\frac{i}{n}$ of the empirical distribution of the target data i.e. $F_n^{-1}(\frac{i}{n})$, and $q_i^F$ the same quantile for a given distribution $F$, i.e. $F^{-1}(\frac{i}{n})$.

**Mean Absolute Error**

$$MAE_F = \frac{1}{n} \sum_{i=1}^{n} |q_i^F - q_i^n|.$$

**Mean Absolute Error over the 95th percentile**

$$MAE95sup_F = \frac{1}{n - \lceil 95\% \rceil \times n} \sum_{i=\lceil 95\% \rceil \times n}^{n} |q_i^F - q_i^n|,$$

where $\lceil \cdot \rceil$ is the ceiling function.

**Root Mean Squared Error**

$$RMSE_F = \sqrt{\frac{1}{n} \sum_{i=1}^{n} (q_i^F - q_i^n)^2}.$$

The number of quantiles $n$ was chosen to be 50 to ensure all locations had enough wet days without giving an artificially high weight to extremes.
Differences of metrics were also computed to increase visibility on the spatial improvements of one distribution over another. For two distributional models F and G respectively, we compute the differences metrics as follows:

$$MAE_{diff_{F,G}} = MAE_{F} - MAE_{G},$$

$$MAE_{95supdiff_{F,G}} = MAE_{95sup_{F}} - MAE_{95sup_{G}},$$

$$RMSE_{diff_{F,G}} = RMSE_{F} - RMSE_{G}.$$

Note that in the following sections, Gamma's error maps will not be shown. We will display only Gamma's results in boxplots due to the very high errors. We do include for some figures the interpolated but non-corrected future model data, referred to as ERA5L in the legends. Moreover, the colour ranges are fixed between maps of the same metrics to represent the 99.9th quantile error of the Stitch-BJ method. Finally, for differences of metrics, the scale is not fixed and may vary from distribution to distribution.

## 3.1 Mean Absolute Error

Differences of MAE for the bias correction performance shown in Figure 7 are quite difficult to assess. Both the spatial pattern and intensity of error are very similar for all chosen distributions. However, some extreme MAE for the EGP and ExpW models can be seen for the DJF season (panels b and c), respectively near the Cévènnes and the Alps region. Those mountainous areas are known to be difficult to model due to their strong orographic features, paired with intense rain events during the SON and DJF seasons. Performances for the JJA season (panels e, f, g and h) are mostly identical which may be caused by the lower amount of rain and fewer extreme events.

<figure><img src="image_9.png" /><figcaption>Fig. 7: Maps of MAE for DJF (first row) and JJA (second row) for the following models: Stitch-BJ panels (panels a,e), EGP (panels b,f), emp (panels c,g) and ExpW (panels d,h).</figcaption></figure>
Differences maps in Figure 8 allow us to see how much the extreme errors are corrected with the Stitch-BJ method, compared to the EGP and ExpW. We refer to *model-season* the error maps corresponding to the differences between Stitch-BJ’s error and *model’s error*, for the given *season*.

For the EGP-DJF (panel a), we see an important increase in performance around the Cévènnes region, with a difference of more than 10mm in terms of MAE, while the EGP is at most 1.46mm of MAE better than the Stitch-BJ. For EGP-JJA (panel d), performances are relatively similar between the Stitch-BJ and the EGP models.

For the ExpW-DJF (panel c), improvements cover a greater area but are less important with at most an improvement of the MAE of 5mm. As for the emp-DJF and emp-JJA (panels b and e), improvements are much less important and performance is relatively equal with no models improving significantly over the other. The only exception may be with the ExpW-JJA (panel f) where the Stitch-BJ was able to improve the MAE of more than 5mm in some areas.

<figure><img src="image_10.png" /><figcaption>**Fig. 8:** Maps of $MAEdiff_{Stitch-BJ,G}$ for DJF (first row) and JJA (second row) with G being the following models: EGP (a,d), emp (b,e) and ExpW (c,f). Note that value ranges are specific to each map.</figcaption></figure>

Boxplots of Figure 9 show how extreme errors are reduced for both the DJF and JJA seasons by using the Stitch-BJ model to replace misfitted portions of the EGP and ExpW distributions. For the DJF season (panel a), the Stitch-BJ was able to improve on the extremes error of both the EGP and ExpW distributions, but fail to reduce the maximum error against the empirical distribution. However, the bulk of the outliers are essentially identical for the Stitch-BJ and empirical distribution. Moreover, in terms of median MAE, the Stitch-BJ and EGP are improving against both the ExpW and empirical models (with an MAE of 0.61, 0.61, 0.64 and 0.67 for the Stitch-BJ, EGP, emp and ExpW respectively). For the JJA season (panel b), conclusions are similar
to the DJF season but at a lesser extent with much less extreme outliers produced by the EGP and ExpW distributions. All models are relatively competitive for the JJA season and the Stitch-BJ and EGP still produce very slightly lower median MAE than other models.

<figure><img src="image_11.png" /><figcaption>**Fig. 9:** Boxplots of MAE for DJF (a) and JJA (b) for all models. Note that the *y*-axis is in log10-scale.</figcaption></figure>

## 3.2 Mean Absolute Error over the 95th percentile

MAE95sup maps in Figure 10 are very similar to the MAE maps of Figure 7. However, for the DJF season, the Stitch-BJ was able to improve against all other shown models significantly in some regions: Cévènnes and Corsica regions for both EGP and ExpW (panels a and c) and Alps and Corsica regions for the empirical model (panel b). However, for the JJA season (panels d, e and f), no significant improvements are easily noticeable.

<figure><img src="image_12.png" /><figcaption>**Fig. 10:** Maps of MAE95sup for DJF (first row) and JJA (second row) for the following models: Stitch-BJ panels (panels a,e), EGP (panels b,f), *emp* (panels c,g) and ExpW (panels d,h).</figcaption></figure>
Differences of MAE95sup in Figure 11 clearly show the advantages of the Stitch-BJ against all other methods for the DJF season, mainly focused on the Cévènes and Corsica region, with errors differences ranging from 30mm to 162mm (panels a, b and c). However, some locations show poor performance from the Stitch-BJ compared to other models, especially in the Eastern Pyrenees where other models produce MAE95sup 17mm to 21mm lower. For the JJA season, performance is equivalent for most models (panels d and e), with only against the ExpW (panel f) where the Stitch-BJ can show higher improvements than the opposite with a maximum improvement of 32mm against 20mm for the ExpW.

<figure><img src="image_13.png" /><figcaption>**Fig. 11:** Maps of MAE95supdiff$_{Stitch-BJ,G}$ for DJF (first row) and JJA (second row) with G being the following models: EGP (a,d), emp (b,e) and ExpW (c,f). Note that value ranges are specific to each map.</figcaption></figure>

In Figure 12, boxplots for the DJF season (panel a) show the clear improvement of the Stitch-BJ against the EGP and ExpW models. At first glance, the ExpW model produces lower outliers than the Stitch-BJ but this is only true for the three most extreme outliers produced by the Stitch-BJ. The rest of the outliers are much more condensed and generally lower than both the ExpW and empirical models’ outliers. Median errors are also slightly lower for the EGP and Stitch-BJ models compared to all other models for the DJF season. Models over the JJA season (panel b) produce very similar extreme outliers, and median errors are also very close and higher than for DJF.
<figure><img src="image_14.png" /><figcaption>**Fig. 12:** Boxplots of MAE95sup for DJF (a) and JJA (b) for all models. Note that the y-axis is in log10-scale.</figcaption></figure>

Additional figures for the RMSE results are available in Appendix A.

## 3.3 Impact of seasonality on performance

In the following figures, all seasons are shown on the same boxplot to assess the performance of the different methods season-wise. For the dry days' probability in Figure 13, we show the difference between the target dry days probability and the modeled one: a positive value means a dry bias, while a negative one implies a wet bias.

<figure><img src="image_15.png" /><figcaption>**Fig. 13:** Multi-season dry days probability difference between CERRAL and Stitch-BJ, EGP, emp, ExpW and uncorrected model ERA5L data.</figcaption></figure>

In Figure 13, the impact of the SSR method is very noticeable for all seasons, with a large reduction in differences in dry days probability for all methods using the SSR correction. The method performed best for the DJF and MAM season, while a small wet bias remains in JJA season and a dry bias in the SON season. Since the underlying model has almost no impact on the SSR correction, all models corrected with the SSR produce similar dry days probability.

For all metrics boxplots in Figures 14 and 15, the *MAEdiff*<sub>Stitch-BJ,mod</sub>, *MAE95supdiff*<sub>Stitch-BJ,mod</sub> and *RMSEdiff*<sub>Stitch-BJ,mod</sub> (see Section 3) are shown
with the EGP, ExpW and empirical distributions as competing models. A positive value means the competing model performed better while a negative one means the Stitch-BJ performed better.

<figure><img src="image_16.png" /><figcaption>**Fig. 14:** Multi-season MAE (a), MAE95sup (b) and RMSE (c) difference of Stitch-BJ, against *emp*, ExpW and uncorrected ERA5L. Note that outliers are not shown in these boxplots for readability.</figcaption></figure>

In Figure 14, we have seasonal boxplots of respectively the MAE, MAE95sup and RMSE differences in panels a, b and c. For the MAE differences in panel (a), the median differences are always negative and the first quartile is consistently lower than the third quartile for all seasons for the empirical and ExpW models. This means that the Stitch-BJ resulted in a lower median MAE for all seasons compared to all the other competing distributions. However, improvements on MAE are marginal with the highest median improvement being for the SON season over the empirical distribution with an improvement of approximately 0.05mm. As for the uncorrected model ERA5L, improvements are more noticeable with a median MAE differences being much lower than the others. We have similar observations in panel (b), with an improvement over the MAE95sup for all seasons over the empirical distribution. The ExpW model produces a median MAE95sup almost identical to the Stitch-BJ model for all seasons.
This may be due to the ExpW being a flexible model as well, being able to model both light and heavy tails distributions. Results are constant throughout the seasons, with JJA and SON producing the most variations when comparing with the empirical distribution. The empirical and uncorrected model were both outperformed in terms of median MAE95sup by the proposed model.

Similarly to Figure 14 (panels a and b), an improvement can be seen over the empirical method for the median differences in panel c. A small improvement on the MAM season is noticeable against the ExpW, while performances are similar for all other seasons. The first and third quartiles are once again asymmetric, favouring the Stitch-BJ method slightly over the two others.

<figure><img src="image_17.png" /><figcaption>Fig. 15: Multi-season RMSE difference between Stitch-BJ and EGP.</figcaption></figure>

In Figure 15, we show the RMSE difference between the Stitch-BJ and EGP models, including outliers in the boxplot. We separated the EGP from the other distributions from Figure 19 to show the extent of the errors of the EGP without impacting the y-axis scale. For the MAM and JJA season, differences are quite symmetrical. This is expected due to the low replacement rate of the EGP for these seasons (see Figure 5). However, for DJF and SON seasons, improvements provided by the Stitch-BJ model over the EGP model are much more noticeable. Many locations produced RMSE 100 to 600mm higher than the Stitch-BJ respective location. This shows that even though most of the locations used the EGP to model the upper tail, the model cannot be used on all the study area at the risk of producing extremely high errors. This emphasizes the importance of the considered stitch model.

## 3.4 Local analysis on a selected location

While Figure 4 gave an idea of how the stitching procedure worked, let's take a corrected location for a finer analysis. In Figure 16, we have quantile-quantile plots of the fitting procedure on both ERA5L and CERRAL data in panels (a, b), and the resulting bias-corrected data in panel (c). In the PBJ test, p-values $k_i$ are produced and are used in the Stitch-BJ model to detect misfitted quantiles. When they cross a given
threshold (red dashed line), we consider the corresponding quantile misfitted and we reject it. In the first two panels, dark green and light green stars correspond to the $k_i$ and rejected $k_i$ respectively. The green and blue vertical lines correspond to the lower and upper rejection index $i_l$ and $i_u$ introduced in Section 2.3, and any quantiles lower (resp. higher) than the index are replaced with ones from another distribution.

<figure><img src="image_18.png" /><figcaption>**Fig. 16:** Local analysis for location (Lon 9.3 , Lat 41.96). Panel (a): Fitting performance of the EGP and ExpW on ERA5L. Panel (b): Fitting performance of the EGP and ExpW on CERRAL. Panel (c): Bias correction performance for Stitch-BJ, EGP, ExpW and empirical model.</figcaption></figure>
In Figure 16 (panel a), we have a local analysis where a full ExpW model has been used on ERA5L data. While the ExpW performed globally worse on high quantiles, the extreme errors of the EGP distribution (50mm on the last quantile) resulted in the selection of the ExpW distribution as the final distributional model. Note that while rejected $k_i$ can be seen in the figure, since we completely replaced the distribution with an ExpW, no lower and upper stitching indexes are shown with vertical lines, meaning that the adopted model is a pure ExpW for this ERA5L’s location. In Figure 16 (panel b), a stitch distribution has been used for this location for CERRAL data. This specific stitch model is composed by an empirical distribution in the lower tail, an EGP distribution in the centre of the law and an ExpW distribution in the upper tail. The notation associated is EMP-EGP-ExpW (see Figure 5).

This resulted in the bias-correction qqplot in panel (c), where the Stitch-BJ is the best fitting model over all contenders. Moreover, we can see the impact of the stitch: if left as is, the EGP produced a large overestimation of almost 1000mm for the highest quantile.

To put into perspective, for this location, the MAE is at respectively 1.1, 3.3, 4.1 and 25.6mm for the Stitch-BJ, empirical, ExpW and EGP models. When considering the MAE95sup which is more focused on the extremes, the errors are 6.3, 29.5, 17 and 440 mm respectively. For differences metrics as shown in Figures 8, 11 and 18, the MAE differences are then respectively -2.2, -3 and -26.5mm against the empirical, ExpW and EGP models. MAE95sup differences exhibit similar results with -23.2, -10.6 and -433.5 mm respectively.

The analysis of Figure 16 shows how the stitching procedure can prevent extremely large errors from a misfitted EGP, and even improve upon the bias-corrected data using and ExpW or the empirical distribution.

# 4 Conclusion and discussion

In this study, our objective was to seasonally bias correct daily precipitation data from ERA5-Land with CERRA-Land data using out-of-sample validation. We applied the Stitch-BJ (Ear et al. 2025) model along with other distributions (see Section 2.3) combined with the Singularity Stochastic Removal method from Vrac et al. (2016) for dry days probability correction. This allowed to show the potential of parametric and semi-parametric bias correction methods against the classical empirical model.

We first showed the necessity of introducing a seasonal separation of our datasets in Section 2.1 and correcting the probability of dry days in Section 2.2. If left as is, residual bias would be left in the final result, which may worsen the performance of all models. To correct the dry-day probabilities of our models, we adapted the SSR methodto apply it to semi and fully-parametric models. This allowed us to obtain well-calibrated dry-day probabilities compared to uncorrected models, with medians of differences with the target dry-day probabilities for the validation period of less than 0.01 across all seasons.

In Section 3, we show the bias correction performance of all models for the DJF and JJA season across multiple metrics. We showed that the Stitch-BJ model produced performance similar to the empirical model in terms of MAE, with a slight reduction
in the median error across the study area. For more extreme focused metrics such as the MAE95sup or the RMSE, the Stitch-BJ model was able to show a reduction in extreme errors against all models for the DJF season, while performance were mostly equivalent for the JJA season. Again, median errors were lower for the Stitch-BJ compared to the empirical model, and error outliers were globally lower, or similar to the empirical model's.

Comparing all seasons together, we do not see an impact on the relative perfor-
mance of models against the Stitch-BJ, with the latter outperforming the ExpW and
empirical models for all seasons.

While this study shows promising results for the Stitch-BJ and demonstrates its flexibility compared to other models, performances could be improved by considering other distributions. All metrics from this study were also solely based on quantiles of wet days intensity distribution, and on the proportion of dry days. Interesting results may arise from studying the impact of such models on climate indicators, which often consider the temporality of such variables.

Moreover, the present analysis shows a high sensitivity to extreme events. Indeed, depending on the user, one could want to produce the worst-case scenario and would rather produce extremes too high than potentially miss them, and in this case, the presented models can be used as such. However, some users might want to be conservative and avoid sending alarming signals and results as much as possible. In this case, multiple methods exist to reduce the sensitivity of the fitting procedure to the few highest data:

* The few highest data points may be removed, either completely from the dataset (right-truncation) or only the information on their intensity (right-censoring), since outliers are known to cause misfitting in some cases as in Berg et al. (2024). Some variation can be found where extreme outliers are not corrected and left as is in (Gutjahr and Heinemann 2013);

* Adding a more robust model for the upper tail modelisation such as the EQM-LIN from (Holthuijzen et al. 2022) which uses a linear correction for extremes above a selected threshold. This allows the smoothening of the upper tail and ignores very high outliers produced by the very few last data quantiles.

Finally, this study's validation and training periods were all historical periods that could be compared with station data. However, climate projections such as data from CMIP6 also need to be bias-corrected. Such data are influenced by green gas emission scenarios which may greatly impact the probability of more extreme events. Such settings may be better suited for parametric distributions for their ability to extrapolate to new extremes and would be an interesting future study to pursue, as in Andrade-Velázquez and Montero-Martínez (2023) or Enyew et al. (2024).

A RMSE results for the bias correction of ERA5L on
the 2010-2020 period using CERRAL data

Maps from Figure ?? are very similar to the one for the MAE95sup in Figure 10, with
mostly similar results for the DJF season except for the Cévènnes, Corsica and Alps
region where the Stitch-BJ may perform better than the other models. For the JJA
season, performance are very similar among all models.

<figure><img src="image_19.png" /><figcaption>**Fig. 17:** Maps of RMSE diff$_{Stitch-BJ,G}$ for DJF (first row) and JJA (second row) with G being the following models: EGP (a,d), emp (b,e) and ExpW (c,f). Note that value ranges are specific to each map.</figcaption></figure>

In Figure 18, spatial patterns and intensity of errors match the ones from Figure 11 for both seasons. The Stitch-BJ method offers an advantage over the DJF season for all methods, while the JJA season shows very close performance for all models. With similar observations to previous figures, conclusion from RMSE boxplots of Figure 19 are almost identical to Figure 12 with a slight improvement over median error for the Stitch-BJ and EGP models, and an improvement on extremes for the Stitch-BJ compared to the EGP and ExpW models for both seasons.

<figure><img src="image_20.png" /><figcaption>**Fig. 19:** Boxplots of RMSE for DJF (a) and JJA (b) for all models. Note that the y-axis is in log10-scale.</figcaption></figure>

**Author contributions:** Conceptualization, P.E., E.D., T.L., A.L and M.T.; methodology, P.E., E.D. and T.L.; software, P.E.; validation, P.E., E.D. and T.L.;
<figure><img src="image_21.png" /><figcaption>**Fig. 18:** Maps of RMSEdiff<sub>Stitch-BJ,G</sub> for DJF (first row) and JJA (second row) with G being the following models: EGP (a,d), emp (b,e) and ExpW (c,f). Note that value ranges are specific to each map.</figcaption></figure>

formal analysis, P.E.; investigation, P.E.; resources, P.E.; data curation, P.E.; writing–original draft preparation, P.E.; writing–review and editing, P.E., E.D. and T.L.; visualization, P.E.; supervision, E.D., T.L. and A.L.; project administration, E.D., T.L., A.L. and M.T.; funding acquisition, E.D., T.L., A.L. and M.T. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work has been partially supported by the French government through the CIFRE funding (CIFRE no 2022/0519) managed by the National Association of Research and Technology and the 3IA Côte d'Azur Investments in the Future project managed by the National Research Agency (ANR-19-P3IA-0002).

**Data availability:** All data, material, and programming codes used in this study are available upon request. ERA5-Land and CERRA-Land datasets analyzed in the current study are available on the Copernicus Climate Change Service (C3S) Climate Data Store.

**Conflict of interests:** The authors declare no conflicts of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

## References

Ambrosino, C., Chandler, R.E., Todd, M.C.: Rainfall-derived growing season characteristics for agricultural impact assessments in South Africa. Theoretical and Applied Climatology **115**(3), 411–426 (2014) https://doi.org/10.1007/s00704-013-0896-y. Accessed 2025-01-21

Argüeso, D., Evans, J.P., Fita, L.: Precipitation bias correction of very high resolution
regional climate models. Hydrology and Earth System Sciences **17**(11), 4379-4388
(2013) https://doi.org/10.5194/hess-17-4379-2013. Accessed 2024-10-29

Ajaaj, A.A., Mishra, A.K., Khan, A.A.: Comparison of BIAS correction techniques for GPCC rainfall data in semi-arid climate. Stochastic Environmental Research and Risk Assessment **30**(6), 1659–1675 (2016) https://doi.org/10.1007/s00477-015-1155-9. Accessed 2024-04-25

Abbott, T.H., Stechmann, S.N., Neelin, J.D.: Long temporal autocorrelations in tropical precipitation data and spike train prototypes. Geophysical Research Letters **43**(21), 11472–11480 (2016) https://doi.org/10.1002/2016GL071282. Accessed 2024-10-29

Alfieri, L., Thielen, J.: A European precipitation index for extreme rain-storm and flash flood early warning. Meteorological Applications **22**(1), 3–13 (2015) https://doi.org/10.1002/met.1328. Accessed 2023-05-17

Andrade-Velázquez, M., Montero-Martínez, M.J.: Statistical Downscaling of Precipitation in the South and Southeast of Mexico. Climate **11**(9), 186 (2023) https://doi.org/10.3390/cli11090186. Accessed 2025-01-13

Berg, P., Bosshard, T., Bozhinova, D., Bärring, L., Löw, J., Nilsson, C., Strandberg, G., Södling, J., Thuresson, J., Wilcke, R., Yang, W.: Robust handling of extremes in quantile mapping – “Murder your darlings”. Geoscientific Model Development **17**(22), 8173–8179 (2024) https://doi.org/10.5194/gmd-17-8173-2024. Accessed 2025-01-21

Bador, M., Boé, J., Terray, L., Alexander, L.V., Baker, A., Bellucci, A., Haarsma, R., Koenigk, T., Moine, M.-P., Lohmann, K., Putrasahan, D.A., Roberts, C., Roberts, M., Scoccimarro, E., Schiemann, R., Seddon, J., Senan, R., Valcke, S., Vanniere, B.: Impact of Higher Spatial Atmospheric Resolution on Precipitation Extremes Over Land in Global Climate Models. Journal of Geophysical Research: Atmospheres **125**(13), 2019–032184 (2020) https://doi.org/10.1029/2019JD032184. Accessed 2025-01-21

Bouvier, C., Cisneros, L., Dominguez, R., Laborde, J.-P., Lebel, T.: Generating rainfall fields using principal components (PC) decomposition of the covariance matrix: a case study in Mexico City. Journal of Hydrology **278**(1), 107–120 (2003) https://doi.org/10.1016/S0022-1694(03)00122-7. Accessed 2025-01-21

Byun, K., Hamlet, A.F.: An improved empirical quantile mapping approach for bias correction of extreme values in climate model simulations. Environmental Research Letters **20**(1), 014041 (2024) https://doi.org/10.1088/1748-9326/ad9b3d. Accessed 2025-01-21

Braunstein, S.L.: How large a sample is needed for the maximum likelihood estimator to be approximately Gaussian? Journal of Physics A: Mathematical and General
25(13), 3813 (1992) https://doi.org/10.1088/0305-4470/25/13/027. Accessed 2025-01-21

Boé, J., Terray, L., Habets, F., Martin, E.: Statistical and dynamical downscaling of the Seine basin climate for hydro-meteorological studies. International Journal of Climatology **27**(12), 1643–1655 (2007) https://doi.org/10.1002/joc.1602. Accessed 2023-09-18

Chen, D., Dai, A., Hall, A.: The Convective-To-Total Precipitation Ratio and the “Drizzling” Bias in Climate Models. Journal of Geophysical Research: Atmospheres **126**(16), 2020–034198 (2021) https://doi.org/10.1029/2020JD034198. Accessed 2024-07-08

Cortés-Hernández, V.E., Caillaud, C., Bellon, G., Brisson, E., Alias, A., Lucas-Picher, P.: Evaluation of the convection permitting regional climate model CNRM-AROME on the orographically complex island of Corsica. Climate Dynamics **62**(6), 4673–4696 (2024) https://doi.org/10.1007/s00382-024-07232-z. Accessed 2025-01-21

Casanueva, A., Kotlarski, S., Herrera, S., Fernández, J., Gutiérrez, J.M., Boberg, F., Colette, A., Christensen, O.B., Goergen, K., Jacob, D., Keuler, K., Nikulin, G., Teichmann, C., Vautard, R.: Daily precipitation statistics in a EURO-CORDEX RCM ensemble: added value of raw and bias-corrected high-resolution simulations. Climate Dynamics **47**(3), 719–737 (2016) https://doi.org/10.1007/s00382-015-2865-x. Accessed 2024-12-03

Chaouche, K., Neppel, L., Dieulin, C., Pujol, N., Ladouche, B., Martin, E., Salas, D., Caballero, Y.: Analyses of precipitation, temperature and evapotranspiration in a French Mediterranean region in the context of climate change. Comptes Rendus Geoscience **342**(3), 234–243 (2010) https://doi.org/10.1016/j.crte.2010.02.001. Accessed 2025-01-21

Derdour, S., Ghenim, A.N., Megnounif, A., Tangang, F., Chung, J.X., Ayoub, A.B.: Bias Correction and Evaluation of Precipitation Data from the CORDEX Regional Climate Model for Monitoring Climate Change in the Wadi Chemora Basin (Northeastern Algeria). Atmosphere **13**(11), 1876 (2022) https://doi.org/10.3390/atmos13111876. Accessed 2025-01-13

Déqué, M.: Frequency of precipitation and temperature extremes over France in an anthropogenic scenario: Model results and statistical correction according to observed values. Global and Planetary Change **57**(1), 16–26 (2007) https://doi.org/10.1016/j.gloplacha.2006.11.030. Accessed 2023-02-06

Enayati, M., Bozorg-Haddad, O., Bazrafshan, J., Hejabi, S., Chu, X.: Bias correction capabilities of quantile mapping methods for rainfall and temperature variables. Journal of Water and Climate Change **12**(2), 401–419 (2021) https://doi.org/10.2166/wcc.2020.261. Accessed 2024-12-03
Ear, P., Di Bernardino, E., Laloë, T., Troin, M., Lambert, A.: A semi-parametric distribution stitch based on the Berk-Jones test for French daily precipitation bias correction. Stochastic Environmental Research and Risk Assessment (2025). to appear. Accessed 2024-10-29

Emmanuel, I., Payrastre, O., Andrieu, H., Zuber, F.: A method for assessing the influence of rainfall spatial variability on hydrograph modeling. First case study in the Cevennes Region, southern France. Journal of Hydrology **555**, 314–322 (2017) https://doi.org/10.1016/j.jhydrol.2017.10.011. Accessed 2025-01-21

Estermann, R., Rajczak, J., Velasquez, P., Lorenz, R., Schär, C.: Projections of Heavy Precipitation Characteristics Over the Greater Alpine Region Using a Kilometer-Scale Climate Model Ensemble. Journal of Geophysical Research: Atmospheres **130**(2), 2024–040901 (2025) https://doi.org/10.1029/2024JD040901. Accessed 2025-01-21

Enyew, F.B., Sahlu, D., Tarekegn, G.B., Hama, S., Debele, S.E.: Performance Evaluation of CMIP6 Climate Model Projections for Precipitation and Temperature in the Upper Blue Nile Basin, Ethiopia. Climate **12**(11), 169 (2024) https://doi.org/10.3390/cli12110169. Accessed 2025-01-13

Fosser, G., Gaetani, M., Kendon, E.J., Adinolfi, M., Ban, N., Belušić, D., Caillaud, C., Careto, J.A.M., Coppola, E., Demory, M.-E., Vries, H., Dobler, A., Feldmann, H., Goergen, K., Lenderink, G., Pichelli, E., Schär, C., Soares, P.M.M., Somot, S., Tölle, M.H.: Convection-permitting climate models offer more certain extreme rainfall projections. npj Climate and Atmospheric Science **7**(1), 1–10 (2024) https://doi.org/10.1038/s41612-024-00600-w. Accessed 2025-01-21

Gutowski, W.J., Decker, S.G., Donavon, R.A., Pan, Z., Arritt, R.W., Takle, E.S.: Temporal-Spatial Scales of Observed and Simulated Precipitation in Central U.S. Climate. Journal of Climate (2003). Accessed 2025-02-28

Gutjahr, O., Heinemann, G.: Comparing precipitation bias correction methods for high-resolution regional climate simulations using COSMO-CLM. Theoretical and Applied Climatology **114**(3), 511–529 (2013) https://doi.org/10.1007/s00704-013-0834-z. Accessed 2023-07-18

Gutiérrez, J.M., Maraun, D., Widmann, M., Huth, R., Hertig, E., Benestad, R., Roessler, O., Wibig, J., Wilcke, R., Kotlarski, S., San Martín, D., Herrera, S., Bedia, J., Casanueva, A., Manzanas, R., Iturbide, M., Vrac, M., Dubrovsky, M., Ribalaygua, J., Pórtoles, J., Räty, O., Räisänen, J., Hingray, B., Raynaud, D., Casado, M.J., Ramos, P., Zerenner, T., Turco, M., Bosshard, T., Štěpánek, P., Bartholy, J., Pongracz, R., Keller, D.E., Fischer, A.M., Cardoso, R.M., Soares, P.M.M., Czernecki, B., Pagé, C.: An intercomparison of a large ensemble of statistical downscaling methods over Europe: Results from the VALUE perfect predictor cross-validation experiment. International Journal of Climatology **39**(9), 3750–3785 (2019) https://doi.org/10.1002/joc.5462. Accessed 2025-02-25
Guo, C., Ning, N., Guo, H., Tian, Y., Bao, A., De Maeyer, P.: Does ERA5-Land Effectively Capture Extreme Precipitation in the Yellow River Basin? Atmosphere 15(10), 1254 (2024) https://doi.org/10.3390/atmos15101254. Accessed 2025-02-25

Holthuijzen, M., Beckage, B., Clemins, P.J., Higdon, D., Winter, J.M.: Robust bias-correction of precipitation extremes using a novel hybrid empirical quantile-mapping method: Advantages of a linear correction for extremes. Theoretical and Applied Climatology 149(1-2), 863–882 (2022) https://doi.org/10.1007/s00704-022-04035-2. Accessed 2023-01-23

Haruna, A., Blanchet, J., Favre, A.-C.: Modeling Intensity-Duration-Frequency Curves for the Whole Range of Non-Zero Precipitation: A Comparison of Models. Water Resources Research 59(6), 2022-033362 (2023) https://doi.org/10.1029/2022WR033362. Accessed 2025-01-21

Henckes, P., Knaut, A., Obermüller, F., Frank, C.: The benefit of long-term high resolution wind data for electricity system analysis. Energy 143, 934–942 (2018) https://doi.org/10.1016/j.energy.2017.10.049. Accessed 2023-05-03

Husak, G.J., Michaelsen, J.C., Funk, C.C.: Use of the Gamma distribution to represent monthly rainfall in Africa for drought monitoring applications. International Journal of Climatology 27(7), 935–944 (2007) https://doi.org/10.1002/joc.1441. Accessed 2024-06-12

Joly, D., Brossard, T., Cardot, H., Cavailhes, J., Hilal, M., Wavresky, P.: Les types de climats en France, une construction spatiale. Cybergeo: European Journal of Geography (2010) https://doi.org/10.4000/cybergeo.23155. Accessed 2023-08-30

Katiraie-Boroujerdy, P.-S., Rahnamay Naeini, M., Akbari Asanjan, A., Chavoshian, A., Hsu, K.-l., Sorooshian, S.: Bias Correction of Satellite-Based Precipitation Estimations Using Quantile Mapping Approach in Different Climate Regions of Iran. Remote Sensing 12(13), 2102 (2020) https://doi.org/10.3390/rs12132102. Accessed 2025-02-24

Khan, S.A.: Exponentiated Weibull regression for time-to-event data. Lifetime Data Analysis 24(2), 328–354 (2018) https://doi.org/10.1007/s10985-017-9394-3. Accessed 2024-01-29

Khan, S., Kuhn, G., Ganguly, A.R., Erickson III, D.J., Ostrouchov, G.: Spatio-temporal variability of daily and weekly precipitation extremes in South America. Water Resources Research 43(11) (2007) https://doi.org/10.1029/2006WR005384. Accessed 2024-10-29

Lafon, T., Dadson, S., Buys, G., Prudhomme, C.: Bias correction of daily precipitation simulated by a regional climate model: a comparison of methods. International Journal of Climatology 33(6), 1367–1381 (2013) https://doi.org/10.1002/joc.3518. Accessed 2023-05-03
Li, B., Huang, Y., Du, L., Wang, D.: Bias Correction for Precipitation Simulated by RegCM4 over the Upper Reaches of the Yangtze River Based on the Mixed Distribution Quantile Mapping Method. Atmosphere 12(12), 1566 (2021) https://doi.org/10.3390/atmos12121566. Accessed 2023-03-28

Langousis, A., Mamalakis, A., Deidda, R., Marrocu, M.: Assessing the relative effec-
tiveness of statistical downscaling and distribution mapping in reproducing rainfall
statistics based on climate model results. Water Resources Research 52(1), 471–494
(2016) https://doi.org/10.1002/2015WR017556. Accessed 2024-10-17

Lavaysse, C., Vrac, M., Drobinski, P., Lengaigne, M., Vischel, T.: Statistical downscaling of the French Mediterranean climate: assessment for present and projection in an anthropogenic scenario. Natural Hazards and Earth System Sciences 12(3), 651–670 (2012) https://doi.org/10.5194/nhess-12-651-2012. Accessed 2025-01-21

Maraun, D.: Bias Correction, Quantile Mapping, and Downscaling: Revisiting the Inflation Issue (2013) https://doi.org/10.1175/JCLI-D-12-00821.1. Accessed 2025-02-24

Mamalakis, A., Langousis, A., Deidda, R., Marrocu, M.: A parametric approach for simultaneous bias correction and high-resolution downscaling of climate model rainfall. Water Resources Research 53(3), 2149–2170 (2017) https://doi.org/10.1002/2016WR019578. Accessed 2024-10-17

Muñoz-Sabater, J., Dutra, E., Agustí-Panareda, A., Albergel, C., Arduini, G., Bal-
samo, G., Boussetta, S., Choulga, M., Harrigan, S., Hersbach, H., Martens, B.,
Miralles, D.G., Piles, M., Rodríguez-Fernández, N.J., Zsoter, E., Buontempo, C.,
Thépaut, J.-N.: ERA5-Land: a state-of-the-art global reanalysis dataset for land
applications. Earth System Science Data 13(9), 4349–4383 (2021) https://doi.org/
10.5194/essd-13-4349-2021. Accessed 2023-05-02

Mudholkar, G.S., Srivastava, D.K., Kollia, G.D.: A Generalization of the Weibull Dis-
tribution with Application to the Analysis of Survival Data. Journal of the American
Statistical Association 91(436), 1575–1583 (1996) https://doi.org/10.2307/2291583
. Accessed 2024-01-29

Michelangeli, P.-A., Vrac, M., Loukos, H.: Probabilistic downscaling approaches: Application to wind cumulative distribution functions. Geophysical Research Letters 36(11) (2009) https://doi.org/10.1029/2009GL038401. Accessed 2022-10-21

Mao, G., Vogl, S., Laux, P., Wagner, S., Kunstmann, H.: Stochastic bias correction of dynamically downscaled precipitation fields for Germany through Copula-based integration of gridded observation data. Hydrology and Earth System Sciences 19(4), 1787–1806 (2015) https://doi.org/10.5194/hess-19-1787-2015. Accessed 2025-01-21
Martinez-Villalobos, C., Neelin, J.D.: Why Do Precipitation Intensities Tend to Follow Gamma Distributions? Journal of the Atmospheric Sciences 76(11), 3611–3631 (2019) https://doi.org/10.1175/JAS-D-18-0343.1. Accessed 2023-03-28

Naveau, P., Huser, R., Ribereau, P., Hannart, A.: Modeling jointly low, moderate, and heavy rainfall intensities without a threshold selection. Water Resources Research 52(4), 2753–2769 (2016) https://doi.org/10.1002/2015WR018552. Accessed 2023-09-20

Pelosi, A.: Performance of the Copernicus European Regional Reanalysis (CERRA) dataset as proxy of ground-based agrometeorological data. Agricultural Water Management 289, 108556 (2023) https://doi.org/10.1016/j.agwat.2023.108556. Accessed 2025-01-21

Prein, A.F., Gobiet, A., Truhetz, H., Keuler, K., Goergen, K., Teichmann, C., Fox Maule, C., Meijgaard, E., Déqué, M., Nikulin, G., Vautard, R., Colette, A., Kjellström, E., Jacob, D.: Precipitation in the EURO-CORDEX 0.11° and 0.44° simulations: high resolution, high benefits? Climate Dynamics 46(1), 383–412 (2016) https://doi.org/10.1007/s00382-015-2589-y. Accessed 2023-05-03

Piani, C., Haerter, J.O., Coppola, E.: Statistical bias correction for daily precipitation in regional climate models over Europe. Theoretical and Applied Climatology 99(1), 187–192 (2010) https://doi.org/10.1007/s00704-009-0134-9. Accessed 2023-01-03

Reiter, P., Gutjahr, O., Schefczyk, L., Heinemann, G., Casper, M.: Does applying quantile mapping to subsamples improve the bias correction of daily precipitation? International Journal of Climatology 38(4), 1623–1633 (2018) https://doi.org/10.1002/joc.5283. Accessed 2025-02-24

Rivoire, P., Martius, O., Naveau, P.: A Comparison of Moderate and Extreme ERA-5 Daily Precipitation With Two Observational Data Sets. Earth and Space Science 8(4), 2020–001633 (2021) https://doi.org/10.1029/2020EA001633. Accessed 2023-10-13

Sangati, M., Borga, M.: Influence of rainfall spatial resolution on flash flood modelling. Natural Hazards and Earth System Sciences 9(2), 575–584 (2009) https://doi.org/10.5194/nhess-9-575-2009. Accessed 2023-05-17

Semenov, M.A., Brooks, R.J., Barrow, E.M., Richardson, C.W.: Comparison of the WGEN and LARS-WG stochastic weather generators for diverse climates. Climate Research 10(2), 95–107 (1998) https://doi.org/10.3354/cr010095. Accessed 2025-01-21

Strohmenger, L., Collet, L., Andréassian, V., Corre, L., Rousset, F., Thirel, G.: Köppen–Geiger climate classification across France based on an ensemble of high-resolution climate projections. Comptes Rendus. Géoscience 356(G1), 67–82 (2024) https://doi.org/10.5802/crgeos.263. Accessed 2025-01-21
Schmidli, J., Frei, C., Vidale, P.L.: Downscaling from GCM precipitation: a benchmark for dynamical and statistical downscaling methods. International Journal of Climatology **26**(5), 679–689 (2006) https://doi.org/10.1002/joc.1287. Accessed 2023-01-03

Sharma, V.K., Singh, S.V., Shekhawat, K.: Exponentiated Teissier distribution with increasing, decreasing and bathtub hazard functions. Journal of Applied Statistics **49**(2), 371–393 (2022) https://doi.org/10.1080/02664763.2020.1813694. Accessed 2024-01-29

Shayeghi, A., Ziveh, A.R., Bakhtar, A., Teymoori, J., Hanel, M., Vargas Godoy, M.R., Markonis, Y., AghaKouchak, A.: Assessing drought impacts on groundwater and agriculture in Iran using high-resolution precipitation and evapotranspiration products. Journal of Hydrology **631**, 130828 (2024) https://doi.org/10.1016/j.jhydrol.2024.130828. Accessed 2024-12-03

Trentini, L., Dal Gesso, S., Venturini, M., Guerrini, F., Calmanti, S., Petitta, M.: A Novel Bias Correction Method for Extreme Events. Climate **11**(1), 3 (2023) https://doi.org/10.3390/cli11010003. Accessed 2023-05-03

Tencaliec, P., Favre, A.-C., Naveau, P., Prieur, C., Nicolet, G.: Flexible semiparametric Generalized Pareto modeling of the entire range of rainfall amount. Environmetrics **31**(2), 2582–1 (2019) https://doi.org/10.1002/env.2582. Accessed 2023-10-12

Themeßl, M.J., Gobiet, A., Heinrich, G.: Empirical-statistical downscaling and error correction of regional climate models and its impact on the climate change signal. Climatic Change **112**(2), 449–468 (2012) https://doi.org/10.1007/s10584-011-0224-4. Accessed 2023-02-08

Vaittinada Ayar, P., Vrac, M., Bastin, S., Carreau, J., Déqué, M., Gallardo, C.: Inter-comparison of statistical and dynamical downscaling models under the EURO- and MED-CORDEX initiative framework: present climate evaluations. Climate Dynamics **46**(3), 1301–1329 (2016) https://doi.org/10.1007/s00382-015-2647-5. Accessed 2025-01-21

Verrelle, A., Glinton, M., Bazile, E., Moigne, P.L.: CERRA-Land: A new land surface reanalysis at 5.5 km resolution over Europe. Technical Report EMS2021-492, Copernicus Meetings (June 2021). https://meetingorganizer.copernicus.org/EMS2021/EMS2021-492.html. Accessed 2024-03-21

Velasquez, P., Messmer, M., Raible, C.C.: A new bias-correction method for precipitation over complex terrain suitable for different climate states: a case study using WRF (version 3.8.1). Geoscientific Model Development **13**(10), 5007–5027 (2020) https://doi.org/10.5194/gmd-13-5007-2020. Accessed 2023-02-23

Vrac, M., Noël, T., Vautard, R.: Bias correction of precipitation through Singularity Stochastic Removal: Because occurrences matter. Journal of Geophysical Research:
Atmospheres **121**(10), 5237-5258 (2016) https://doi.org/10.1002/2015JD024511. Accessed 2022-11-04

Vigaud, N., Vrac, M., Caballero, Y.: Probabilistic downscaling of GCM scenarios over southern India. International Journal of Climatology **33**(5), 1248–1263 (2013) https://doi.org/10.1002/joc.3509. Accessed 2025-01-21

Vautard, R., Yiou, P., Oldenborgh, G.-J.v., Lenderink, G., Thao, S., Ribes, A., Plan-
ton, S., Dubuisson, B., Soubeyroux, J.-M.: Extreme Fall 2014 Precipitation in
the Cévennes Mountains. Bulletin of the American Meteorological Society (2015)
https://doi.org/10.1175/BAMS-D-15-00088.1. Accessed 2025-01-21

Wang, C., Zhang, L., Lee, S.-K., Wu, L., Mechoso, C.R.: A global perspective on CMIP5 climate model biases. Nature Climate Change **4**(3), 201–205 (2014) https://doi.org/10.1038/nclimate2118. Accessed 2024-11-14

Xu, H., Xu, C.-Y., Sælthun, N.R., Zhou, B., Xu, Y.: Evaluation of reanalysis and satellite-based precipitation datasets in driving hydrological models in a humid region of Southern China. Stochastic Environmental Research and Risk Assessment **29**(8), 2003–2020 (2015) https://doi.org/10.1007/s00477-014-1007-z. Accessed 2024-07-09

Şan, M., Nacar, S., Kankal, M., Bayram, A.: Daily precipitation performances of regression-based statistical downscaling models in a basin with mountain and semi-arid climates. Stochastic Environmental Research and Risk Assessment **37**(4), 1431–1455 (2023) https://doi.org/10.1007/s00477-022-02345-5. Accessed 2024-07-08