What's the difference in using more nquantiles?

[Joaogmr472]

Setup Information

• Xclim version:
• Python version:
• Operating System:

Context

I'm using the xclim module applied for bias correction (xclim.sdba.adjustment) for daily data from the historical period (7305 times) and future from the CMIP6 models.
In the training section (sdba.QuantileDeltaMapping.train) the option of nquantiles is offered. In the description of the library it is shown that the default value is 20, but it does not explain why. I did tests with higher values, up to the time length of the historical data.
My question is: What is the difference between setting the value of nquantiles to 20 and leaving nquantiles=len(hist.time)?

Code of Conduct

• I agree to follow this project's Code of Conduct

[aulemahal]

nquantiles is the number of points xclim uses to approximate the CDF of the timeseries. During the train step, for each gridpoint and group (if you're using that option, the timeseries's distribution is summarized into a array of N quantiles. In the case of QuantileDeltaMapping, this means that there will be N adjustment factors per grid point.

As this number is essentially the resolution to which the distribution is resolved, it's usually a matter of performance : the more points you have, the heavier the computation is. I am not sure of this, but I would fear creating situations with overfitting if you used "all" quantiles, meaning you're merely sorting the timeseries, not computing any quantile per se.

But I have to dig a bit more in the literature to provide a better answer.

However, as I hinted above, the use of group makes it so that passing nquantiles=len(hist.time) doesn't make much sense because you'll have more quantiles than then length of each sub-timeseries. Would you be interested in an nquantiles='all' option to explicitly use all quantiles? This could be implemented, and because we would sort the timeseries rather than find a distribution, it should be faster than nquantiles=len(hist.time).

EDIT: I see that Cannon's own R package uses all quantiles per default... I guess that nulls my overfitting fear and makes the nquantiles='all' option even more interesting. See https://search.r-project.org/CRAN/refmans/MBC/html/QDM.html

[aulemahal]

Hi @Joaogmr472,

It took some time before I had time to come back to this issue. I discussed it with others and did some quick tests and I have a few conclusions :

One issue with taking as many quantiles as values is indeed overfitting, but particularly on the extremes. An example I actually encountered with xclim is related to issue #1015. If we were to use the sorted timeseries as I suggested above, it would assign the quantiles 0 and 1 to the minimum and maximum of the series. In the adjust step, this has the effect of clipping the scenario to the range of the reference, at least within the reference period.

The default way xclim chooses its quantiles (after this issue was solved) is to divide the [0, 1] range in N + 1 parts (where N = nquantiles) and to use only the inner boundaries. But as N increases, we get close to this "clipping" problem of the extremes. And this is one reason why the default nquantiles is kept relatively small. (FYI, this is done here)

In a recent project I used 50, which mean the min and max quantiles are 0.01 and 0.99. It seems to have performed well without overfitting the extremes (but this is hard to exhaustively analyze...). Also I was using DetrendedQuantileMapping.

However, this is not heavily documented, but you can pass an array to nquantiles. So you could construct an array of quantiles with a high "resolution" in the middle of the distribution and a coarser one on in the tails, avoiding quantiles too close to 0 or 1. Also, on the contrary, you could also use this to imitate using the CDF of the timeseries directlly (instead of approximating it), this was my suggestion in the previous comment when I was speculating about nquantiles='all'. It would go like: nquantiles=np.linspace(0, 1, num=len(hist.time), endpoint=True).

The last technique would still be slower than my proposition, and would be inexact when grouping is used. I do not have time to implement something to support this nquantiles='all', but if you do, we're open to PRs and we will happily help you for your contribution!

I hope that answers your question better. I know it's been a long time, but hopefully, this will help other people as well...