Repository: sktime/sktime
Pull request: [BUG] Fix `predict_residuals` internal data type conversion
Fixes https://github.com/sktime/sktime/issues/4766

The `predict_residuals` internal calculation of residuals could break due to different data types for the observed and the predicted, as in some cases the observed would be obtained from `self._y` which can be data type coerced.

The "optimal" solution would be to not do any coercion or conversion, but for now this should be fixed by converting to the predicted which has the data type consistent with the `fit` input.

Problem statement:
[BUG] `predict_residuals(y=None)` does not work if forecaster is fitted with a Series.
**Describe the bug**
I fitted a simple pipeline using a `pd.Series` as input for fit. Calling `predict_residuals()` lead to an exception.

**To Reproduce**
``` python
from sktime.datasets import load_longley
from sktime.forecasting.model_selection import temporal_train_test_split
from sktime.forecasting.sarimax import SARIMAX
from sktime.transformations.series.difference import Differencer

y, X = load_longley()
y_train, y_test, X_train, X_test = temporal_train_test_split(y, X)
pipe = Differencer() * SARIMAX()
pipe.fit(y=y_train, X=X_train, fh=[1, 2, 3, 4])
result = pipe.predict_residuals()
```
This raises the following exception

```
TypeError: y must have same type, dims, index as expected predict return. expected type <class 'pandas.core.series.Series'>, but found <class 'pandas.co
re.frame.DataFrame'>
```

Note the following command that should be identical to the above example works:

```
y, X = load_longley()
y_train, y_test, X_train, X_test = temporal_train_test_split(y, X)
pipe = Differencer() * SARIMAX()
pipe.fit(y=y_train, X=X_train, fh=[1, 2, 3, 4])
result = pipe.predict_residuals(y=y_train) # The difference is that here that y_train is passed to the parameter y.
```

**Expected behavior**
The default predict_residuals should also work if the estimator is fitted with a pd.Series. 

**Additional context**

The reason of this bug is that `self.y` is used in the default version, which is a pd.DataFrame. However, the return type of the `predict` is a `pd.Series`. 
I think the underlying is some inconsistencies in the fit method of `_base` line 334 to 338, where y_inner is passed to `_update_y_X`, which is a DataFrame and not a Series anymore. 

```python
X_inner, y_inner = self._check_X_y(X=X, y=y)

# set internal X/y to the new X/y
# this also updates cutoff from y
self._update_y_X(y_inner, X_inner)

```


I see two ways of fixing this bug:
* ensuring that self.y is of the same type, which is passed to fit and not a transformed version, would require either changing the return type of `check_X_y` or to pass `y` instead of `y_inner` to `_update_y_X`
* retransforming self.y to the type of predict or vise versa.

Any considerations or preferences on that?

**Versions**
sktime 19.2
python 3.9.16 (windows)

Hints from issue discussion:
Confirmed on current main, python 3.11, windows

Commit message: bugfix

Patch to review:
diff --git a/sktime/forecasting/base/_base.py b/sktime/forecasting/base/_base.py
index e09e0ff28b7..c25c19d2827 100644
--- a/sktime/forecasting/base/_base.py
+++ b/sktime/forecasting/base/_base.py
@@ -1129,10 +1129,7 @@ def predict_residuals(self, y=None, X=None):
         y_pred = self.predict(fh=fh, X=X)
 
         if not type(y_pred) == type(y):
-            raise TypeError(
-                "y must have same type, dims, index as expected predict return. "
-                f"expected type {type(y_pred)}, but found {type(y)}"
-            )
+            y = convert_to(y, self._y_mtype_last_seen)
 
         y_res = y - y_pred
 
diff --git a/sktime/forecasting/base/tests/test_base_bugs.py b/sktime/forecasting/base/tests/test_base_bugs.py
index 6a92d303792..3d556c52f1f 100644
--- a/sktime/forecasting/base/tests/test_base_bugs.py
+++ b/sktime/forecasting/base/tests/test_base_bugs.py
@@ -10,8 +10,10 @@
     ForecastingGridSearchCV,
 )
 from sktime.forecasting.reconcile import ReconcilerForecaster
+from sktime.forecasting.sarimax import SARIMAX
 from sktime.forecasting.trend import PolynomialTrendForecaster
 from sktime.transformations.hierarchical.aggregate import Aggregator
+from sktime.transformations.series.difference import Differencer
 from sktime.utils._testing.hierarchical import _make_hierarchical
 from sktime.utils.validation._dependencies import _check_estimator_deps
 
@@ -58,3 +60,21 @@ def test_heterogeneous_get_fitted_params():
 
     reconciler.fit(y_agg)
     reconciler.get_fitted_params()  # triggers an error pre-fix
+
+
+@pytest.mark.skipif(
+    not _check_estimator_deps(SARIMAX, severity="none"),
+    reason="skip test if required soft dependency not available",
+)
+def test_predict_residuals_conversion():
+    """Regression test for bugfix #4766, related to predict_residuals internal type."""
+    from sktime.datasets import load_longley
+    from sktime.forecasting.model_selection import temporal_train_test_split
+
+    y, X = load_longley()
+    y_train, y_test, X_train, X_test = temporal_train_test_split(y, X)
+    pipe = Differencer() * SARIMAX()
+    pipe.fit(y=y_train, X=X_train, fh=[1, 2, 3, 4])
+    result = pipe.predict_residuals()
+
+    assert type(result) == type(y_train)


Code review comments:
diff_hunk: @@ -58,3 +60,21 @@ def test_heterogeneous_get_fitted_params():
 
     reconciler.fit(y_agg)
     reconciler.get_fitted_params()  # triggers an error pre-fix
+
+
+@pytest.mark.skipif(
+    not _check_estimator_deps(SARIMAX, severity="none"),
+    reason="skip test if required soft dependency not available",
+)
+def test_predict_residuals_conversion():
+    """Regression test for bugfix #4766, related to predict_residuals internal type."""
+    from sktime.datasets import load_longley
+    from sktime.forecasting.model_selection import temporal_train_test_split
+
+    y, X = load_longley()
+    y_train, y_test, X_train, X_test = temporal_train_test_split(y, X)
+    pipe = Differencer() * SARIMAX()
path: sktime/forecasting/base/tests/test_base_bugs.py
text: @user1:
Not Blocking, but should we replace SARIMAX with NaiveForecaster? With that change, this test would also be executed if the soft dependency is not available. 


@author:
Hm, perhaps - although `NaiveForecaster` ignores `X`. The issue could also be specific to the inner types (I didn't check), so I simply went with a known failure.

Did you ccheck whether `NaiveForecaster` also causes the issue (pre-PR)?

@user1:
I tested it a minute ago and 

```python
   from sktime.datasets import load_longley
   from sktime.forecasting.model_selection import temporal_train_test_split
   
   y, X = load_longley()
   y_train, y_test, X_train, X_test = temporal_train_test_split(y, X)
   pipe = Differencer() * NaiveForecaster()
   pipe.fit(y=y_train, X=X_train, fh=[1, 2, 3, 4])
   result = pipe.predict_residuals()
   
```

results in the same error.

@author:
ok, replaced with `NaiveForecaster` now.