Repository: scikit-learn/scikit-learn
Pull request: Fix ValueError on HistGradientBoostingRegressor with least_absolute_deviation loss function fit call with sample_weight parameter
#### Reference Issues/PRs
Fixes #19400.


Problem statement:
HistGradientBoostingRegressor with 'least_absolute_deviation' loss function and sample_weight raises ValueError: `indices` and `arr` must have the same number of dimensions
#### Describe the bug
`HistGradientBoostingRegressor` with `least_absolute_deviation` loss function raises ValueError on `fit` call with `sample_weight` parameter.

#### Steps/Code to Reproduce
```python
import numpy as np
from sklearn.experimental import enable_hist_gradient_boosting
from sklearn.ensemble import HistGradientBoostingRegressor

n = 500000
x = np.random.uniform(-1, 1, [n, 3])
y = np.random.uniform(-1, 1, n)
sample_weight = np.random.uniform(0, 1, n)
gb = HistGradientBoostingRegressor(loss='least_absolute_deviation')
gb.fit(x, y, sample_weight=sample_weight)
```

#### Expected Results
No error is thrown.

#### Actual Results
```
Traceback (most recent call last):
  File "/snap/pycharm-professional/230/plugins/python/helpers/pydev/pydevd.py", line 1477, in _exec
    pydev_imports.execfile(file, globals, locals)  # execute the script
  File "/snap/pycharm-professional/230/plugins/python/helpers/pydev/_pydev_imps/_pydev_execfile.py", line 18, in execfile
    exec(compile(contents+"\n", file, 'exec'), glob, loc)
  File "/home/vadim/Загрузки/tt.py", line 11, in <module>
    gb.fit(x, y, sample_weight=sample_weight)
  File "/home/vadim/projects/monza/venv/lib/python3.8/site-packages/sklearn/ensemble/_hist_gradient_boosting/gradient_boosting.py", line 466, in fit
    self._loss.update_leaves_values(grower, y_train,
  File "/home/vadim/projects/monza/venv/lib/python3.8/site-packages/sklearn/ensemble/_hist_gradient_boosting/loss.py", line 264, in update_leaves_values
    median_res = _weighted_percentile(y_true[indices]
  File "/home/vadim/projects/monza/venv/lib/python3.8/site-packages/sklearn/utils/stats.py", line 43, in _weighted_percentile
    sorted_weights = _take_along_axis(sample_weight, sorted_idx, axis=0)
  File "/home/vadim/projects/monza/venv/lib/python3.8/site-packages/sklearn/utils/fixes.py", line 172, in _take_along_axis
    return np.take_along_axis(arr=arr, indices=indices, axis=axis)
  File "<__array_function__ internals>", line 5, in take_along_axis
  File "/home/vadim/projects/monza/venv/lib/python3.8/site-packages/numpy/lib/shape_base.py", line 170, in take_along_axis
    return arr[_make_along_axis_idx(arr_shape, indices, axis)]
  File "/home/vadim/projects/monza/venv/lib/python3.8/site-packages/numpy/lib/shape_base.py", line 34, in _make_along_axis_idx
    raise ValueError(
ValueError: `indices` and `arr` must have the same number of dimensions
```

#### Versions
Output from `sklearn.show_versions()`:
```System:
    python: 3.8.5 (default, Jul 28 2020, 12:59:40)  [GCC 9.3.0]
executable: /home/vadim/projects/monza/venv/bin/python
   machine: Linux-5.8.0-41-generic-x86_64-with-glibc2.29

Python dependencies:
          pip: 20.3.3
   setuptools: 51.3.3
      sklearn: 0.24.0
        numpy: 1.19.5
        scipy: 1.6.0
       Cython: None
       pandas: 1.2.0
   matplotlib: 3.3.3
       joblib: 1.0.0
threadpoolctl: 2.1.0

Built with OpenMP: True
```

#### Possible bug location
It seems like bug in `sklearn.ensemble._hist_gradient_boosting.loss.py` file `LeastAbsoluteDeviation.update_leaves_values` method, line 264:
```
median_res = _weighted_percentile(y_true[indices]
                                                  - raw_predictions[indices],
                                                  sample_weight=sample_weight,  # -> sample_weight[indices]
                                                  percentile=50)
```

Hints from issue discussion:
Thanks for the report and for the reproducible example, I can reproduce.

The fix is simply to update

```py
                median_res = _weighted_percentile(y_true[indices]
                                                  - raw_predictions[indices],
                                                  sample_weight=sample_weight,
                                                  percentile=50
```

to use `sample_weight[indices]`.

Would you like to submit a PR? We'd need a small non-regression test as well, your snippet above would be perfect (with a lower number of samples).

Commit message: fix code style error

Patch to review:
diff --git a/sklearn/ensemble/_hist_gradient_boosting/loss.py b/sklearn/ensemble/_hist_gradient_boosting/loss.py
index f256408bf01fb..4bbf59dc01088 100644
--- a/sklearn/ensemble/_hist_gradient_boosting/loss.py
+++ b/sklearn/ensemble/_hist_gradient_boosting/loss.py
@@ -261,10 +261,11 @@ def update_leaves_values(self, grower, y_true, raw_predictions,
                 median_res = np.median(y_true[indices]
                                        - raw_predictions[indices])
             else:
-                median_res = _weighted_percentile(y_true[indices]
-                                                  - raw_predictions[indices],
-                                                  sample_weight=sample_weight,
-                                                  percentile=50)
+                median_res = _weighted_percentile(
+                    y_true[indices] - raw_predictions[indices],
+                    sample_weight=sample_weight[indices],
+                    percentile=50
+                )
             leaf.value = grower.shrinkage * median_res
             # Note that the regularization is ignored here
 
diff --git a/sklearn/ensemble/_hist_gradient_boosting/tests/test_gradient_boosting.py b/sklearn/ensemble/_hist_gradient_boosting/tests/test_gradient_boosting.py
index c501125059c8f..dc394c096bfd9 100644
--- a/sklearn/ensemble/_hist_gradient_boosting/tests/test_gradient_boosting.py
+++ b/sklearn/ensemble/_hist_gradient_boosting/tests/test_gradient_boosting.py
@@ -203,6 +203,18 @@ def test_least_absolute_deviation():
     assert gbdt.score(X, y) > .9
 
 
+def test_least_absolute_deviation_sample_weight():
+    # non regression test for issue #19400
+    # make sure no error is thrown during fit of
+    # HistGradientBoostingRegressor with least_absolute_deviation loss function
+    n = 100
+    X = np.random.uniform(-1, 1, [n, 2])
+    y = np.random.uniform(-1, 1, n)
+    sample_weight = np.random.uniform(0, 1, n)
+    gbdt = HistGradientBoostingRegressor(loss='least_absolute_deviation')
+    gbdt.fit(X, y, sample_weight=sample_weight)
+
+
 @pytest.mark.parametrize('y', [([1., -2., 0.]), ([0., 0., 0.])])
 def test_poisson_y_positive(y):
     # Test that ValueError is raised if either one y_i < 0 or sum(y_i) <= 0.


Code review comments:
diff_hunk: @@ -203,6 +203,18 @@ def test_least_absolute_deviation():
     assert gbdt.score(X, y) > .9
 
 
+def test_least_absolute_deviation_sample_weight():
+    # non regression test for issue #19400
+    # make sure no error is thrown during fit of
+    # HistGradientBoostingRegressor with least_absolute_deviation loss function
+    n = 100
+    X = np.random.uniform(-1, 1, [n, 2])
+    y = np.random.uniform(-1, 1, n)
+    sample_weight = np.random.uniform(0, 1, n)
path: sklearn/ensemble/_hist_gradient_boosting/tests/test_gradient_boosting.py
text: @user1:
Small nit: we like to set the RNG (even if it doesn't really matter here) and use more descriptive names


```suggestion
	rng = np.random.RandomState(0)
    n_samples = 100
    X = rng.uniform(-1, 1, size=(n_samples, 2))
    y = rng.uniform(-1, 1, size=n_samples)
    sample_weight = rng.uniform(0, 1, size=n_samples)
```

@author:
OK, I'll take note