issue_comments

The alignment of data and weights is not done automatically. So I agree this would be ideal, but I'll need some guidance to make it happen without the align call.

This is the error I get python ValueError: operand to apply_ufunc has required core dimensions ['dim_0', 'dim_1'], but some of these dimensions are absent on an input variable: ['dim_1']

@mathause I'm not 100% confident the methods other than linear work as expected, so I suggest we do not expose method until a robust test suite can confirm results are valid.

I have added a test to verify that using equal weights with the different interpolation methods that this PR supports would work (i.e. would yield the same results as np.quantile, with the corresponding methods). It is skipped however because the method argument is not currently exposed in the API (it would be in future work, ideally).

We have local changes that expose method in the public API, and compare results against numpy for the case with uniform weights to confirm results are identical.

Do you want us to push those changes here ?

@mathause This PR goes beyond what is currently implemented in numpy. For now, all weighted quantiles PR on numpy are more or less based on "linear" method (method 7) and none have been merged. I plan to work on integrating weights with the other interpolation methods but don't have the time right now. I'll probably pick some ideas from here.

As for the numerics here, everything looks good to me. The only limitations I can see are: - This only handles sampling weights, which is fine I guess. - Some interpolation methods are missing, they can be added later. - ~A nan_weighted_quantile could also be interesting to add~

Correct, it's not exposed yet because I don't have the bandwidth to create tests for all the different methods.

The equal weight case could be tested against numpy fairly easily though.

@mathause Should this PR be amended to account for #6108 ?

FYI, weighted quantiles topic will be discussed in numpy's triage meeting of today (17:00 UTC). I'm not a maintainer but I'm sure you are welcomed to join in if you are interested. Meeting information: https://hackmd.io/68i_JvOYQfy9ERiHgXMPvg

I had the same interrogation. My guess is that the DataArray.quantile and DataArrayWeighted.quantile should have similar signatures.

The way I see it is that another PR would bring the DataArray.quantile signature up to date with numpy 1.22, which includes the 9 h types, and harmonize DataArrayWeighted.quantile with it.

@huard's latest commit modifies the algorithm so that it uses Kish's effective sample size, as described in the blog where the algorithm comes from: https://aakinshin.net/posts/kish-ess-weighted-quantiles/, which seems to solve the problem mentioned by @mathause.

Also he adds support for the interpolation types 4 to 9 (those that share a common way of computing Qp, as described here: https://en.wikipedia.org/wiki/Quantile#Estimating_quantiles_from_a_sample

@mathause About this:

I did some tries and got an unexpected result:

```python data = xr.DataArray([0, 1, 2, 3]) weights = xr.DataArray([1, 0, 1, 0]) data.weighted(weights).quantile([0.75])

np.quantile([0, 2], 0.75) ```

Can you double-check? Or do I misunderstand something?

My latest commit should fix (and test) this.

@mathause Thanks for the many excellent suggestions! After having removed the for loop the way you suggested, I tried to address this:

The algorithm is quite clever but it multiplies all elements (except 2) with 0 - this could maybe be sped up by only using the relevant elements.

At first I thought that something like this could work:

python w = np.diff(v) nz = np.nonzero(w) d2 = np.tile(data, (n, 1)) r = w[nz] * d2[nz] r = r[::2] + r[1::2]

The problem however is that it turns out that w's rows sometimes have one element only, instead of two (which is when an h coincides exactly with a weight value, instead of lying between two). Given that difficulty, my impression is that it's not really solvable, or at least not in a way that would result in a more efficient version.