Repository: scikit-learn/scikit-learn
Pull request: FIX ravel prediction of `PLSRegression` when fitted on 1d `y`
#### Reference Issues/PRs

Fixes: #26549.

#### What does this implement/fix? Explain your changes.

This implements https://github.com/scikit-learn/scikit-learn/issues/26549#issuecomment-1594874351 for `PLSRegression`. Please let me know if the alternative (https://github.com/scikit-learn/scikit-learn/issues/26549#issuecomment-1594874351) is more desired, or if I should implement this for other multi-output regressors as well in this PR.


Problem statement:
PLSRegression not working with VotingRegressor
Hi, 

Unfortunately, PLSRegression does not work in VotingRegressor. Likely related to PLSRegression returning predictions in shape of (n_samples, 1) instead of (n_samples,) like other single-target regressors. 

thanks!

```python
import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.cross_decomposition import PLSRegression
from sklearn.ensemble import VotingRegressor

r1 = LinearRegression()
r2 = PLSRegression()

X = np.array([[1, 1], [2, 4], [3, 9], [4, 16], [5, 25], [6, 36]])
y = np.array([2, 6, 12, 20, 30, 42])

r1.fit(X, y).predict(X).shape
```
```
(6,)
```

```python
r2.fit(X, y).predict(X).shape
```
```
(6, 1)
```

```python
er = VotingRegressor([('lr', r1), ('plsr', r2)])
er.fit(X, y).predict(X)
```
```python-traceback
Traceback (most recent call last):
  Cell In[25], line 1
    er.fit(X, y).predict(X)
  File ~/code/scikit-learn/sklearn/ensemble/_voting.py:625 in predict
    return np.average(self._predict(X), axis=1, weights=self._weights_not_none)
  File ~/code/scikit-learn/sklearn/ensemble/_voting.py:69 in _predict
    return np.asarray([est.predict(X) for est in self.estimators_]).T
ValueError: setting an array element with a sequence. The requested array has an inhomogeneous shape after 2 dimensions. The detected shape was (2, 6) + inhomogeneous part.
```

Hints from issue discussion:
`PLSRegresson `supports multi-target regression, hence the predicition is 2d-array.

As a proposed solution, I suggest creating a custom subclass of `PLSRegression`, which we can call `PLS1D`. This subclass overrides the `predict `method of `PLSRegression `to ensure it returns a 1-dimensional array, making it compatible with the `VotingRegressor.` This will effectively handle the case where `n_components=1`.

Here is the relevant code for the proposed PLS1D subclass:

```
import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.cross_decomposition import PLSRegression
from sklearn.ensemble import VotingRegressor



class PLS1D(PLSRegression):
    def __init__(self, n_components=1, scale=True, max_iter=500, tol=1e-06, copy=True):
        if n_components > 1:
            raise ValueError("PLS1D can only handle n_components=1")
        super().__init__(n_components=n_components, scale=scale, max_iter=max_iter, tol=tol, copy=copy)
    
    def predict(self, X):
        y_pred = super().predict(X)
        return y_pred.ravel()


r1 = LinearRegression()
r2 = PLS1D()

X = np.array([[1, 1], [2, 4], [3, 9], [4, 16], [5, 25], [6, 36]])
y = np.array([2, 6, 12, 20, 30, 42])
er = VotingRegressor([('lr', r1), ('plsr', r2)])
print(er.fit(X, y).predict(X))

r1.fit(X, y).predict(X).shape
r2.fit(X, y).predict(X).shape
```

Output
`[ 1.21007429  6.09332284 12.55111269 20.58344384 30.1903163  41.37173005]`
I think that we might want `PLSRegression` to behave consistently with `LinearRegression` and `Ridge`, that is, when `y` is provided as shape `(n_samples, 1)` at `fit`, then `predict` should output that same shape, but when `y` is passed as shape `(n_samples,)` then `predict` would automatically ravel its output.

I think other multi-output regressors such as `RandomForestRegressor` and other tree-based regressors should also do the same.
We might also accept a PR to make `VotingRegressor` regressor do such raveling automatically when all inputs have a mix of `(n_samples,)` and `(n_samples, 1)` shapes.
I would like to make a PR for the first suggestion (i.e., automatically ravel the prediction when `Y.ndim == 1` when fitting). I often needed to ravel manually when using scikit-learn so I think this is desirable, but I do think this is a breaking change.

Commit message: changelof entry updated

Patch to review:
diff --git a/doc/whats_new/v1.3.rst b/doc/whats_new/v1.3.rst
index 7bdbe7841f0d4..cdcde6561ff30 100644
--- a/doc/whats_new/v1.3.rst
+++ b/doc/whats_new/v1.3.rst
@@ -276,6 +276,13 @@ Changelog
   :func:`covariance.graphical_lasso_path`, and :class:`covariance.GraphicalLassoCV`.
   :pr:`26033` by :user:`Genesis Valencia <genvalen>`.
 
+:mod:`sklearn.cross_decomposition`
+..................................
+
+- |Enhancement| :class:`cross_decomposition.PLSRegression` now automatically ravels
+  the output of `predict` if fitted with one dimensional `y`.
+  :pr:`26602` by :user:`Yao Xiao <Charlie-XIAO>`.
+
 :mod:`sklearn.datasets`
 .......................
 
diff --git a/sklearn/cross_decomposition/_pls.py b/sklearn/cross_decomposition/_pls.py
index da395d8f060fb..815e1804f10de 100644
--- a/sklearn/cross_decomposition/_pls.py
+++ b/sklearn/cross_decomposition/_pls.py
@@ -236,7 +236,10 @@ def fit(self, X, Y):
             Y, input_name="Y", dtype=np.float64, copy=self.copy, ensure_2d=False
         )
         if Y.ndim == 1:
+            self._predict_1d = True
             Y = Y.reshape(-1, 1)
+        else:
+            self._predict_1d = False
 
         n = X.shape[0]
         p = X.shape[1]
@@ -467,8 +470,8 @@ def predict(self, X, copy=True):
         # Normalize
         X -= self._x_mean
         X /= self._x_std
-        Ypred = X @ self.coef_.T
-        return Ypred + self.intercept_
+        Ypred = X @ self.coef_.T + self.intercept_
+        return Ypred.ravel() if self._predict_1d else Ypred
 
     def fit_transform(self, X, y=None):
         """Learn and apply the dimension reduction on the train data.
diff --git a/sklearn/cross_decomposition/tests/test_pls.py b/sklearn/cross_decomposition/tests/test_pls.py
index 8f4840c9b9f21..add67cc950a3a 100644
--- a/sklearn/cross_decomposition/tests/test_pls.py
+++ b/sklearn/cross_decomposition/tests/test_pls.py
@@ -622,3 +622,17 @@ def test_pls_set_output(Klass):
     assert isinstance(y_trans, np.ndarray)
     assert isinstance(X_trans, pd.DataFrame)
     assert_array_equal(X_trans.columns, est.get_feature_names_out())
+
+
+def test_pls_regression_fit_1d_y():
+    """Check that when fitting with 1d y, prediction should also be 1d.
+
+    Non-regression test for Issue #26549.
+    """
+    X = np.array([[1, 1], [2, 4], [3, 9], [4, 16], [5, 25], [6, 36]])
+    y = np.array([2, 6, 12, 20, 30, 42])
+    expected = y.copy()
+
+    plsr = PLSRegression().fit(X, y)
+    y_pred = plsr.predict(X)
+    assert_allclose(y_pred, expected)


Code review comments:
diff_hunk: @@ -276,6 +276,13 @@ Changelog
   :func:`covariance.graphical_lasso_path`, and :class:`covariance.GraphicalLassoCV`.
   :pr:`26033` by :user:`Genesis Valencia <genvalen>`.
 
+:mod:`sklearn.cross_decomposition`
+..................................
+
+- |Enhancement| :class:`cross_decomposition.PLSRegression` now automatically ravels
path: doc/whats_new/v1.3.rst
text: @user1:
```suggestion
- |Fix| :class:`cross_decomposition.PLSRegression` now automatically ravels
```
This seems more like a fix than an enhancement.