Repository: sktime/sktime
Pull request #5657 (2023-12-24 14:31:13)

[ENH] Replace `isinstance` by `object_type` tag based checks
<!--
Welcome to sktime, and thanks for contributing!
Please have a look at our contribution guide:
https://www.sktime.net/en/latest/get_involved/contributing.html
-->

#### Reference Issues/PRs
Without an issue.

Discussed in #5656


#### What does this implement/fix? Explain your changes.
Replace `isinstance` or `issubclass` based checks by `object_type` tag.

#### Does your contribution introduce a new dependency? If yes, which one?

No

#### What should a reviewer concentrate their feedback on?

<!-- This section is particularly useful if you have a pull request that is still in development. You can guide the reviews to focus on the parts that are ready for their comments. We suggest using bullets (indicated by * or -) and filled checkboxes [x] here -->

#### Did you add any tests for the change?
No
#### Any other comments?
<!--
We value all user contributions, no matter how small or complex they are. If you have any questions, feel free to post
in the dev-chat channel on the sktime discord https://discord.com/invite/54ACzaFsn7. If we are slow to review (>3 working days), likewise feel free to ping us on discord. Thank you for your understanding during the review process.
-->

#### PR checklist
<!--
Please go through the checklist below. Please feel free to remove points if they are not applicable.
-->

##### For all contributions
- [ ] I've added myself to the [list of contributors](https://github.com/sktime/sktime/blob/main/CONTRIBUTORS.md) with any new badges I've earned :-)
  How to: add yourself to the [all-contributors file](https://github.com/sktime/sktime/blob/main/.all-contributorsrc) in the `sktime` root directory (not the `CONTRIBUTORS.md`). Common badges: `code` - fixing a bug, or adding code logic. `doc` - writing or improving documentation or docstrings. `bug` - reporting or diagnosing a bug (get this plus `code` if you also fixed the bug in the PR).`maintenance` - CI, test framework, release.
  See here for [full badge reference](https://allcontributors.org/docs/en/emoji-key)
- [x] The PR title starts with either [ENH], [MNT], [DOC], or [BUG]. [BUG] - bugfix, [MNT] - CI, test framework, [ENH] - adding or improving code, [DOC] - writing or improving documentation or docstrings.

<!--
Thanks for contributing!
-->



Patch:
diff --git a/README.md b/README.md
index 538099695..a435ea1b1 100644
--- a/README.md
+++ b/README.md
@@ -4,7 +4,7 @@
 
 > A unified interface for machine learning with time series
 
-:rocket: **Version 0.24.1 out now!** [Check out the release notes here](https://www.sktime.net/en/latest/changelog.html).
+:rocket: **Version 0.25.0 out now!** [Check out the release notes here](https://www.sktime.net/en/latest/changelog.html).
 
 sktime is a library for time series analysis in Python. It provides a unified interface for multiple time series learning tasks. Currently, this includes time series classification, regression, clustering, annotation, and forecasting. It comes with [time series algorithms](https://www.sktime.net/en/stable/estimator_overview.html) and [scikit-learn] compatible tools to build, tune and validate time series models.
 
@@ -169,8 +169,8 @@ flexible choice of soft dependencies is unavailable via `conda`)
 ``` python
 from sktime.datasets import load_airline
 from sktime.forecasting.base import ForecastingHorizon
-from sktime.forecasting.model_selection import temporal_train_test_split
 from sktime.forecasting.theta import ThetaForecaster
+from sktime.split import temporal_train_test_split
 from sktime.performance_metrics.forecasting import mean_absolute_percentage_error
 
 y = load_airline()
diff --git a/build_tools/run_blogposts.sh b/build_tools/run_blogposts.sh
new file mode 100755
index 000000000..95e00b97b
--- /dev/null
+++ b/build_tools/run_blogposts.sh
@@ -0,0 +1,11 @@
+#!/bin/bash
+
+# Script to run all blog post notebooks.
+set -euxo pipefail
+
+CMD="jupyter nbconvert --to notebook --inplace --execute --ExecutePreprocessor.timeout=600"
+
+for notebook in examples/blog_posts/*.ipynb; do
+  echo "Running: $notebook"
+  $CMD "$notebook"
+done
diff --git a/docs/source/api_reference/transformations.rst b/docs/source/api_reference/transformations.rst
index 25513ffd2..0bbc48474 100644
--- a/docs/source/api_reference/transformations.rst
+++ b/docs/source/api_reference/transformations.rst
@@ -665,6 +665,19 @@ Bootstrap transformations
     SplitterBootstrapTransformer
     STLBootstrapTransformer
 
+Panel-to-Series transformers
+----------------------------
+
+These transformers create a single series from a panel.
+
+.. currentmodule:: sktime.transformations.merger
+
+.. autosummary::
+    :toctree: auto_generated/
+    :template: class.rst
+
+    Merger
+
 Outlier detection, changepoint detection
 ----------------------------------------
 
diff --git a/docs/source/changelog.rst b/docs/source/changelog.rst
index 4006f3e90..2fa48c3d4 100644
--- a/docs/source/changelog.rst
+++ b/docs/source/changelog.rst
@@ -19,6 +19,340 @@ For upcoming changes and next releases, see our `milestones <https://github.com/
 For our long-term plan, see our :ref:`roadmap`.
 
 
+Version 0.25.0 - 2023-12-26
+---------------------------
+
+Release with base class updates and scheduled deprecations:
+
+* framework support for multioutput classifiers, regressors
+  (:pr:`5408`, :pr:`5651`, :pr:`5662`) :user:`Vasudeva-bit`, :user:`fkiraly`
+* framework support for panel-to-series transformers (:pr:`5351`) :user:`benHeid`
+* scheduled deprecations
+
+For last larger feature update, see 0.24.2.
+
+Core interface changes
+~~~~~~~~~~~~~~~~~~~~~~
+
+Time series classification and regression
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* the base class framework now supports multioutput classifiers or regressors.
+  All classifiers and regressors are now able to make multioutput predictions,
+  including all third party classifiers and regressors.
+  A multioutput ``y`` can now be passed, in the form of a 2D ``np.ndarray`` or
+  ``pd.DataFrame``, with one column per output.
+  The ``predict`` method will then return a predicted output of the same type.
+  To retain downwards compatibility, ``predict`` will always return a 1D ``np.ndarray``
+  for univariate outputs, this is currently not subject to deprecation.
+
+* Genuinely multioutput classifiers and regressors are labelled with the new
+  tag ``capability:multioutput`` being ``True``.
+  All other classifiers and regressors broadcast by column of ``y``,
+  and a parallelization backend can be selected via ``set_config``,
+  by setting the ``backend:parallel`` and ``backend:parallel:params`` configuration
+  flags, see the ``set_config`` docstring for details.
+  Broadcasting extends automatically to all existing third party classifiers
+  and regressors via base class inheritance once ``sktime`` is updated,
+  the estimator classes themselves do not need to be updated.
+
+* classifiers and regressors now have a tag ``y_inner_mtype``, this allows extenders
+  to specify an internal ``mtype``, of ``Table`` scitype.
+  The mtype specified i the tag is the guaranteed
+  mtype of ``y`` seen in the private ``_fit`` method.
+  The default is the same as previously
+  implicit, the ``numpy1D`` mtype.
+  Therefore, third party classifiers and regressors do not need to be updated,
+  and should be fully upwards compatible.
+
+Transformations
+^^^^^^^^^^^^^^^
+
+* the base class framework now supports transformations that aggregate ``Panel`` data
+  to ``Series`` data, i.e., panel-to-series transformers, e.g., averaging.
+  Such transformers are identified by the tags
+  ``scitype:transform-input`` being ``"Panel"``,
+  and ``scitype:transform-output`` being ``"Series"``.
+  An example is ``Merger``.
+
+Deprecations and removals
+~~~~~~~~~~~~~~~~~~~~~~~~~
+
+Benchmarking, Metrics, Splitters
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* time series splitters, i.e., descendants of ``BaseSplitter``, have moved from
+  ``sktime.forecasting.model_selection`` to ``sktime.split``.
+  They are no longer available in the old location
+  ``sktime.forecasting.model_selection``, since 0.25.0.
+  Forecasting tuners are still present in ``sktime.forecasting.model_selection``,
+  and their locationn is not subject to deprecation.
+
+* in forecasting ``evaluate``, the order of columns in the return data frame
+  has changed. Users should consult the docstring of ``evaluate`` for details.
+
+* in forecasting ``evaluate``, the ``compute`` argument was removed,
+  after deprecation in 0.24.0.
+  Its purpose was to distinguish lazy or eager evaluation in
+  the ``dask`` parallelization backend.
+  To switch between lazy and eager evaluation, users should instead
+  select ``dask`` or ``dask_lazy`` via the ``backend`` parameter.
+
+* in forecasting ``evaluate``, ``kwargs`` are deprecated, removal has been
+  moved to 0.26.0. Users should pass backend parameters via the ``backend_params``
+  parameter instead.
+
+
+Contents
+~~~~~~~~
+
+* [ENH] Multioutput capability for all time series classifiers and regressors, broadcasting and tag (:pr:`5408`) :user:`Vasudeva-bit`
+* [ENH] Support for panel-to-series transformers, merger transformation (:pr:`5351`) :user:`benHeid`
+* [ENH] allow object ``dtype``-s in ``pandas`` based ``Table`` mtype-s (:pr:`5651`) :user:`fkiraly`
+* [ENH] intermediate base class for panel tasks - classification, regression (:pr:`5662`) :user:`fkiraly`
+* [MNT] CI element to test blogpost notebooks (:pr:`5663`) :user:`fkiraly`, :user:`yarnabrina`
+* [MNT] 0.25.0 deprecations and change actions (:pr:`5613`) :user:`fkiraly`
+
+Contributors
+~~~~~~~~~~~~
+
+:user:`benHeid`,
+:user:`fkiraly`,
+:user:`Vasudeva-bit`,
+:user:`yarnabrina`
+
+Version 0.24.2 - 2023-12-24
+---------------------------
+
+Highlights
+~~~~~~~~~~
+
+* ``FunctionParamFitter`` for custom parameter switching, e.g., applying forecaster or transformer
+  conditional on instance properties (:pr:`5630`) :user:`tpvasconcelos`
+* ``calibration_plot`` for probabilistic forecasts (:pr:`5632`) :user:`benHeid`
+* ``prophet`` based piecewise linear trend forecaster (:pr:`5592`) :user:`sbuse`
+* new transformer: dilation mapping (:pr:`5557`) :user:`fspinna`
+* custom ``joblib`` backends are now supported in parallelization via ``set_config`` (:pr:`5537`) :user:`fkiraly`
+
+Dependency changes
+~~~~~~~~~~~~~~~~~~
+
+* ``dask`` (data container and parallelization back-end) bounds have been updated to ``<2023.12.2``.
+* ``holidays`` (transformations soft dependency) bounds have been updated to ``>=0.29,<0.40``.
+
+Core interface changes
+~~~~~~~~~~~~~~~~~~~~~~
+
+Forecasting
+^^^^^^^^^^^
+
+* ``fit_predict`` now allows specification of ``X_pred`` argument for ``predict``.
+  If passed, ``X_pred`` is used as ``X`` in ``predict``, instead of ``X``.
+  This is useful for forecasters that expect ``X`` to be subset to the
+  forecasting horizon.
+* custom ``joblib`` backends for hierarchical and multivariate forecast broadcasting
+  are now supported. To use a custom ``joblib`` backend, use ``set_config`` to
+  set the ``backend:parallel`` configuration flag to ``"joblib"``,
+  and set the ``backend`` parameter in the ``dict`` set via ``backend:parallel:params``
+  to the name of the custom ``joblib`` backend. Further bakcend parameters
+  can be passed in the same ``dict``. See docstring of ``set_config`` for details.
+
+Time series classification
+^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* In ``SimpleRNNClassifier``, the ``num_epochs`` parameter is deprecated and has been
+  renamed to ``n_epochs``. ``num_epochs`` can be used until ``sktime`` 0.25.last,
+  but will be removed in ``sktime`` 0.26.0. A deprecation warning is raised if
+  ``num_epochs`` is used.
+
+Time series regression
+^^^^^^^^^^^^^^^^^^^^^^
+
+* In ``SimpleRNNRegressor``, the ``num_epochs`` parameter is deprecated and has been
+  renamed to ``n_epochs``. ``num_epochs`` can be used until ``sktime`` 0.25.last,
+  but will be removed in ``sktime`` 0.26.0. A deprecation warning is raised if
+  ``num_epochs`` is used.
+
+Transformations
+^^^^^^^^^^^^^^^
+
+* custom ``joblib`` backends for hierarchical and multivariate transformer broadcasting
+  are now supported. To use a custom ``joblib`` backend, use ``set_config`` to
+  set the ``backend:parallel`` configuration flag to ``"joblib"``,
+  and set the ``backend`` parameter in the ``dict`` set via ``backend:parallel:params``
+  to the name of the custom ``joblib`` backend. Further bakcend parameters
+  can be passed in the same ``dict``. See docstring of ``set_config`` for details.
+
+Enhancements
+~~~~~~~~~~~~
+
+BaseObject and base framework
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* [ENH] improved error messages for input checks in base classes (:pr:`5510`) :user:`fkiraly`
+* [ENH] support for custom ``joblib`` backends in parallelization (:pr:`5537`) :user:`fkiraly`
+* [ENH] consistent use of ``np.ndarray`` for mtype tags (:pr:`5648`) :user:`fkiraly`
+* [ENH] set output format parameter in ``sktime`` internal ``check_is_mtype`` calls to silence deprecation warnings (:pr:`5563`) :user:`benHeid`
+
+Benchmarking, Metrics, Splitters
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* [ENH] cutoff and forecasting horizon ``loc`` based splitter (:pr:`5575`) :user:`fkiraly`
+* [ENH] enable tag related registry tests for ``splitter`` estimator type (:pr:`5576`) :user:`fkiraly`
+
+Data types, checks, conversions
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* [ENH] ``sklearn`` facing coercion utility for ``pd.DataFrame``, to ``str`` columns (:pr:`5550`) :user:`fkiraly`
+* [ENH] ``deep_equals`` - clearer return on diffs from ``dtypes`` and ``index``, relaxation of ``MultiIndex`` equality check (:pr:`5560`) :user:`fkiraly`
+* [ENH] Uniformization of ``pandas`` index types in mtypes (:pr:`5561`) :user:`fkiraly`
+* [ENH] ``n_features`` and ``feature_names`` metadata field for time series mtypes (:pr:`5596`) :user:`fkiraly`
+
+Forecasting
+^^^^^^^^^^^
+
+* [ENH] expected forecast prediction index utility in ``ForecastingHorizon`` (:pr:`5501`) :user:`fkiraly`
+* [ENH] refactor index generation in reducers to use ``ForecastingHorizon`` method (:pr:`5539`) :user:`fkiraly`
+* [ENH] fix index name check for reduction forecasters (:pr:`5543`) :user:`fkiraly`
+* [ENH] forecaster ``fit_predict`` with ``X_pred`` argument for ``predict`` (:pr:`5562`) :user:`fkiraly`
+* [ENH] refactor ``DirectReductionForecaster``to use ``sklearn`` input coercion utility (:pr:`5581`) :user:`fkiraly`
+* [ENH] export and test ``DirectReductionForecaster`` (:pr:`5582`) :user:`fkiraly`
+* [ENH] ``prophet`` based piecewise linear trend forecaster (:pr:`5592`) :user:`sbuse`
+* [ENH] Add ``fit_kwargs`` to ``Prophet`` (:pr:`5597`) :user:`tpvasconcelos`
+* [ENH] ``Croston`` test parameters - integer smoothing parameter (:pr:`5608`) :user:`NguyenChienFelix33`
+* [ENH] ``prophet`` adapter - safer handling of ``fit_kwargs`` (:pr:`5622`) :user:`fkiraly`
+
+Parameter estimation and hypothesis testing
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* [ENH] Add new ``FunctionParamFitter`` parameter estimator (:pr:`5630`) :user:`tpvasconcelos`
+
+Time series annotation
+^^^^^^^^^^^^^^^^^^^^^^
+* [ENH] Change ``GGS`` to inherit from ``BaseSeriesAnnotator`` (:pr:`5315`) :user:`Alex-JG3`
+
+Time series classification
+^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* [ENH] enable testing ``MrSQM`` for persistence in ``nsfa>0`` case after upstream bugfix (:pr:`5171`) :user:`fkiraly`
+* [ENH] ``num_epochs`` renamed to ``n_epochs`` in ``SimpleRNNClassifier`` and ``SimpleRNNRegressor`` (:pr:`5607`) :user:`aeyazadil`
+
+Time series clustering
+^^^^^^^^^^^^^^^^^^^^^^
+
+* [ENH] enable tag related registry tests for ``clusterer`` estimator type (:pr:`5576`) :user:`fkiraly`
+
+Transformations
+^^^^^^^^^^^^^^^
+
+* [ENH] dilation mapping transformer (:pr:`5557`) :user:`fspinna`
+* [ENH] second test parameter set for ``TSFreshRelevantFeatureExtractor`` (:pr:`5623`) :user:`fkiraly`
+
+Visualization
+^^^^^^^^^^^^^
+
+* [ENH] Add ``calibration_plot`` for probabilistic forecasts (:pr:`5632`) :user:`benHeid`
+
+Test framework
+^^^^^^^^^^^^^^
+
+* [ENH] reactivate and fix ``test_multiprocessing_idempotent`` (:pr:`5573`) :user:`fkiraly`
+* [ENH] test class register, refactor ``check_estimator`` test gathering to central location (:pr:`5574`) :user:`fkiraly`
+* [ENH] conditional testing of objects - test if covering test class has changed (:pr:`5579`) :user:`fkiraly`
+
+
+Fixes
+~~~~~
+
+BaseObject and base framework
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* [BUG] fix ``scitype`` ``coerce_to_list`` parameter, add test coverage (:pr:`5578`) :user:`fkiraly`
+
+Data types, checks, conversions
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+* [BUG] Fix typos in mtype tags ``np.ndarray``, from erroneous ``nd.array`` (:pr:`5645`) :user:`yarnabrina`
+
+Forecasting
+^^^^^^^^^^^
+
+* [BUG] in ``ARCH``, fix ``str`` coercion of ``pd.Series`` name (:pr:`5407`) :user:`Vasudeva-bit`
+* [BUG] in reduced regressor, copy or truncate ``X`` if it does not fit the forecasting horizon (:pr:`5542`) :user:`benHeid`
+* [BUG] pass correct level argument from ``StatsForecastBackAdapter`` to ``statsforecast`` (:pr:`5587`) :user:`sd2k`
+* [BUG] fix ``HierarchyEnsembleForecaster`` returned unexpected predictions if data had only one hierarchy level and forecasters specified by node (:pr:`5615`) :user:`VyomkeshVyas`
+* [BUG] fix loss of time zone attribute in ``ForecastingHorizon.to_absolute`` (:pr:`5628`) :user:`fkiraly`
+* [BUG] change index match to integer in ``_StatsModelsAdapter`` predict (:pr:`5642`) :user:`ciaran-g`
+
+Transformations
+^^^^^^^^^^^^^^^
+
+* [BUG] ``TsFreshFeatureExtractor`` - correct wrong forwarded parameter name ``profiling`` (:pr:`5600`) :user:`sssilvar`
+* [BUG] Correct inference of ``TransformerPipeline`` output type tag (:pr:`5625`) :user:`fkiraly`
+
+Visualization
+^^^^^^^^^^^^^
+
+* [BUG] Fix multiple figures created by ``plot_windows`` (:pr:`5636`) :user:`benHeid`
+
+
+Maintenance
+~~~~~~~~~~~
+
+* [MNT] CI Modifications (:pr:`5498`) :user:`yarnabrina`
+* [MNT] rename variables in base (:pr:`5502`) :user:`yarnabrina`
+* [MNT] addressing various ``pandas`` related deprecations (:pr:`5583`) :user:`fkiraly`
+* [MNT] Update pre commit hooks (:pr:`5646`) :user:`yarnabrina`
+* [MNT] [Dependabot](deps-dev): Update ``pytest-xdist`` requirement from ``<3.4,>=3.3`` to ``>=3.3,<3.5`` (:pr:`5551`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``dask`` requirement from ``<2023.7.1`` to ``<2023.11.1`` (:pr:`5552`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``dask`` requirement from ``<2023.11.1`` to ``<2023.12.2`` (:pr:`5629`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``holidays`` requirement from ``<0.36,>=0.29`` to ``>=0.29,<0.37`` (:pr:`5538`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``holidays`` requirement from ``<0.37,>=0.29`` to ``>=0.29,<0.38`` (:pr:`5565`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``holidays`` requirement from ``<0.38,>=0.29`` to ``>=0.29,<0.40`` (:pr:`5637`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``sphinx-gallery`` requirement from ``<0.15.0`` to ``<0.16.0`` (:pr:`5566`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``pytest-xdist`` requirement from ``<3.5,>=3.3`` to ``>=3.3,<3.6`` (:pr:`5567`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps-dev): Update ``pycatch22`` requirement from ``<0.4.4`` to ``<0.4.5`` (:pr:`5542`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps): Bump actions/download-artifact from 3 to 4 (:pr:`5627`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps): Bump actions/setup-python from 4 to 5 (:pr:`5605`) :user:`dependabot[bot]`
+* [MNT] [Dependabot](deps): Bump actions/upload-artifact from 3 to 4 (:pr:`5626`) :user:`dependabot[bot]`
+
+Documentation
+~~~~~~~~~~~~~
+
+* [DOC] splitter full API reference page (:pr:`5577`) :user:`fkiraly`
+* [DOC] Correct ReST syntax in "RocketClassifier" (:pr:`5564`) :user:`rahulporuri`
+* [DOC] Added notebook accompanying Joanna Lenczuk's blog post for testing (:pr:`5604`) :user:`onyekaugochukwu`, :user:`joanlenczuk`
+* [DOC] Remove extra parameter in docstring with incorrect definition (:pr:`5617`) :user:`wayneadams`
+* [DOC] fix and complete ``YfromX`` docstring (:pr:`5593`) :user:`fkiraly`
+* [DOC] fix typo in ``AA_datatypes_and_datasets.ipynb`` panel data loading example (:pr:`5594`) :user:`fkiraly`
+* [DOC] forecasting ``evaluate`` utility - improved algorithm description in docstring #5603  (:pr:`5603`) :user:`adamkells`
+* [DOC] add explanation about fit/transform instance linking behaviour of rocket transformers (:pr:`5621`) :user:`fkiraly`
+* [DOC] Adjust ``FunctionTransformer``'s docstring (:pr:`5634`) :user:`tpvasconcelos`
+* [DOC] fixed typo in ``pytest.mark.skipif`` (:pr:`5640`) :user:`yarnabrina`
+
+Contributors
+~~~~~~~~~~~~
+
+:user:`adamkells`,
+:user:`aeyazadil`,
+:user:`Alex-JG3`,
+:user:`benHeid`,
+:user:`ciaran-g`,
+:user:`fkiraly`,
+:user:`fspinna`,
+:user:`joanlenczuk`,
+:user:`NguyenChienFelix33`,
+:user:`onyekaugochukwu`,
+:user:`rahulporuri`,
+:user:`sbuse`,
+:user:`sd2k`,
+:user:`sssilvar`,
+:user:`tpvasconcelos`,
+:user:`Vasudeva-bit`,
+:user:`VyomkeshVyas`,
+:user:`wayneadams`,
+:user:`yarnabrina`
+
 Version 0.24.1 - 2023-11-05
 ---------------------------
 
@@ -343,7 +677,7 @@ Benchmarking, Metrics, Splitters
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
 * time series splitters, i.e., descendants of ``BaseSplitter``, have moved from
-  ``sktime.forecasting.model_selection`` to ``sktime.`split``.
+  ``sktime.forecasting.model_selection`` to ``sktime.split``.
   The old location ``model_selection`` is deprecated and will be removed in 0.25.0.
   Until 0.25.0, it is still available but will raise an informative warning message.
 
@@ -412,7 +746,7 @@ Transformations
 * [ENH] add proper ``inverse_transform`` to ``STLTransformer`` (:pr:`5300`) :user:`fkiraly`
 * [ENH] ``joblib`` and ``dask`` backends in broadcasting of estimators in multivariate or hierarchical case - part 1, ``VectorizedDF.vectorize_est`` (:pr:`5267`) :user:`fkiraly`
 * [ENH] ``joblib`` and ``dask`` backends in broadcasting of estimators in multivariate or hierarchical case - part 2, base class config (:pr:`5301`) :user:`fkiraly`
-* [ENH] Refactor of `DateTimeFeatures` tests to `pytest` fixtures (:pr:`5397`) :user:`adamkells`
+* [ENH] Refactor of ``DateTimeFeatures`` tests to ``pytest`` fixtures (:pr:`5397`) :user:`adamkells`
 
 Testing framework
 ^^^^^^^^^^^^^^^^^
diff --git a/docs/source/get_started.rst b/docs/source/get_started.rst
index cce05ec13..5295cf594 100644
--- a/docs/source/get_started.rst
+++ b/docs/source/get_started.rst
@@ -75,9 +75,9 @@ Forecasting
 
     >>> from sktime.datasets import load_airline
     >>> from sktime.forecasting.base import ForecastingHorizon
-    >>> from sktime.forecasting.model_selection import temporal_train_test_split
     >>> from sktime.forecasting.theta import ThetaForecaster
     >>> from sktime.performance_metrics.forecasting import mean_absolute_percentage_error
+    >>> from sktime.split import temporal_train_test_split
 
     >>> y = load_airline()
     >>> y_train, y_test = temporal_train_test_split(y)
diff --git a/docs/source/user_guide/introduction.rst b/docs/source/user_guide/introduction.rst
index 672dca8a6..70f6c7fd8 100644
--- a/docs/source/user_guide/introduction.rst
+++ b/docs/source/user_guide/introduction.rst
@@ -31,8 +31,8 @@ can simply write:
     from sktime.datasets import load_airline
     from sktime.forecasting.compose import make_reduction
     from sklearn.ensemble import RandomForestRegressor
-    from sktime.forecasting.model_selection import temporal_train_test_split
     from sktime.performance_metrics.forecasting import MeanAbsolutePercentageError
+    from sktime.split import temporal_train_test_split
 
     y = load_airline()
     y_train, y_test = temporal_train_test_split(y)
diff --git a/extension_templates/classification.py b/extension_templates/classification.py
index 532743cd5..884d7f77f 100644
--- a/extension_templates/classification.py
+++ b/extension_templates/classification.py
@@ -71,9 +71,11 @@ class MyTimeSeriesClassifier(BaseClassifier):
     # these are the default values, only add if different to these.
     _tags = {
         "X_inner_mtype": "numpy3D",  # which type do _fit/_predict accept, usually
-        # this is either "numpy3D" or "nested_univ" (nested pd.DataFrame). Other
+        "y_inner_mtype": "numpy1D",  # which type do _fit/_predict return, usually
+        # this is either "numpy3D", "pd-multiindex" or "nested_univ" (nested df). Other
         # types are allowable, see datatypes/panel/_registry.py for options.
-        "capability:multivariate": False,
+        "capability:multivariate": False,  # ability to handle multivariate X
+        "capability:multioutput": False,  # ability to predict multiple columns in y
         "capability:unequal_length": False,
         "capability:missing_values": False,
         "capability:train_estimate": False,
@@ -129,13 +131,19 @@ class MyTimeSeriesClassifier(BaseClassifier):
         ----------
         X : guaranteed to be of a type in self.get_tag("X_inner_mtype")
             if self.get_tag("X_inner_mtype") = "numpy3D":
-                3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
-            if self.get_tag("X_inner_mtype") = "nested_univ":
-                pd.DataFrame with each column a dimension, each cell a pd.Series
+            3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            if self.get_tag("X_inner_mtype") = "pd-multiindex:":
+            pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
             for list of other mtypes, see datatypes.SCITYPE_REGISTER
             for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of int, of shape [n_instances] - class labels for fitting
-            indices correspond to instance indices in X
+        y : guaranteed to be of a type in self.get_tag("y_inner_mtype")
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            if self.get_tag("capaility:multioutput") = False, guaranteed to be 1D
+            if self.get_tag("capaility:multioutput") = True, guaranteed to be 2D
 
         Returns
         -------
@@ -152,31 +160,32 @@ class MyTimeSeriesClassifier(BaseClassifier):
         #   3. read from self in _fit,  4. pass to interfaced_model.fit in _fit
 
     # todo: implement this, mandatory
-    def _predict(self, X) -> np.ndarray:
+    def _predict(self, X):
         """Predict labels for sequences in X.
 
         private _predict containing the core logic, called from predict
 
-        State required:
-            Requires state to be "fitted".
-
-        Accesses in self:
-            Fitted model attributes ending in "_"
-
         Parameters
         ----------
         X : guaranteed to be of a type in self.get_tag("X_inner_mtype")
             if self.get_tag("X_inner_mtype") = "numpy3D":
-                3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
-            if self.get_tag("X_inner_mtype") = "nested_univ":
-                pd.DataFrame with each column a dimension, each cell a pd.Series
+            3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            if self.get_tag("X_inner_mtype") = "pd-multiindex:":
+            pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
             for list of other mtypes, see datatypes.SCITYPE_REGISTER
             for specifications, see examples/AA_datatypes_and_datasets.ipynb
 
         Returns
         -------
-        y : 1D np.array of int, of shape [n_instances] - predicted class labels
+        y : should be of mtype in self.get_tag("y_inner_mtype")
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            predicted class labels
             indices correspond to instance indices in X
+            if self.get_tag("capaility:multioutput") = False, should be 1D
+            if self.get_tag("capaility:multioutput") = True, should be 2D
         """
 
         # implement here
@@ -200,9 +209,11 @@ class MyTimeSeriesClassifier(BaseClassifier):
         ----------
         X : guaranteed to be of a type in self.get_tag("X_inner_mtype")
             if self.get_tag("X_inner_mtype") = "numpy3D":
-                3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
-            if self.get_tag("X_inner_mtype") = "nested_univ":
-                pd.DataFrame with each column a dimension, each cell a pd.Series
+            3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            if self.get_tag("X_inner_mtype") = "pd-multiindex:":
+            pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
             for list of other mtypes, see datatypes.SCITYPE_REGISTER
             for specifications, see examples/AA_datatypes_and_datasets.ipynb
 
diff --git a/pyproject.toml b/pyproject.toml
index 0d86f8947..f2c3b7235 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -1,6 +1,6 @@
 [project]
 name = "sktime"
-version = "0.24.1"
+version = "0.25.0"
 description = "A unified framework for machine learning with time series"
 readme = "README.md"
 keywords = [
@@ -203,7 +203,7 @@ transformations = [
   "holidays>=0.29,<0.40",
   "mne<1.6,>=1.5",
   'numba<0.59,>=0.53; python_version < "3.12"',
-  "pycatch22<0.4.4,>=0.4",
+  "pycatch22>=0.4,<0.4.5",
   "pykalman-bardo<0.10,>=0.9.7",
   "statsmodels<0.15,>=0.12.1",
   'stumpy<1.13,>=1.5.1; python_version < "3.12"',
diff --git a/sktime/__init__.py b/sktime/__init__.py
index 41607cda3..5ca449788 100644
--- a/sktime/__init__.py
+++ b/sktime/__init__.py
@@ -1,6 +1,6 @@
 """sktime."""
 
-__version__ = "0.24.1"
+__version__ = "0.25.0"
 
 __all__ = ["show_versions"]
 
diff --git a/sktime/base/__init__.py b/sktime/base/__init__.py
index 4a30db2d0..6cf431c83 100644
--- a/sktime/base/__init__.py
+++ b/sktime/base/__init__.py
@@ -6,10 +6,12 @@ __author__ = ["mloning", "RNKuhns", "fkiraly"]
 __all__ = [
     "BaseObject",
     "BaseEstimator",
+    "BasePanelMixin",
     "_HeterogenousMetaEstimator",
     "load",
 ]
 
 from sktime.base._base import BaseEstimator, BaseObject
+from sktime.base._base_panel import BasePanelMixin
 from sktime.base._meta import _HeterogenousMetaEstimator
 from sktime.base._serialize import load
diff --git a/sktime/base/_base.py b/sktime/base/_base.py
index 6ff54472b..677b35257 100644
--- a/sktime/base/_base.py
+++ b/sktime/base/_base.py
@@ -78,7 +78,15 @@ class BaseObject(_BaseObject):
     Extends skbase BaseObject with additional features.
     """
 
-    _config = {"warnings": "on"}
+    _config = {
+        "warnings": "on",
+        "backend:parallel": None,  # parallelization backend for broadcasting
+        #  {None, "dask", "loky", "multiprocessing", "threading"}
+        #  None: no parallelization
+        #  "loky", "multiprocessing" and "threading": uses `joblib` Parallel loops
+        #  "dask": uses `dask`, requires `dask` package in environment
+        "backend:parallel:params": None,  # params for parallelization backend,
+    }
 
     _config_doc = {
         "display": """
@@ -101,6 +109,36 @@ class BaseObject(_BaseObject):
             * "on" = will raise warnings from sktime
             * "off" = will not raise warnings from sktime
         """,
+        "backend:parallel": """
+        backend:parallel : str, optional, default="None"
+            backend to use for parallelization when broadcasting/vectorizing, one of
+
+            - "None": executes loop sequentally, simple list comprehension
+            - "loky", "multiprocessing" and "threading": uses ``joblib.Parallel``
+            - "joblib": custom and 3rd party ``joblib`` backends, e.g., ``spark``
+            - "dask": uses ``dask``, requires ``dask`` package in environment
+        """,
+        "backend:parallel:params": """
+        backend:parallel:params : dict, optional, default={} (no parameters passed)
+            additional parameters passed to the parallelization backend as config.
+            Valid keys depend on the value of ``backend:parallel``:
+
+            - "None": no additional parameters, ``backend_params`` is ignored
+            - "loky", "multiprocessing" and "threading": default ``joblib`` backends
+              any valid keys for ``joblib.Parallel`` can be passed here, e.g.,
+              ``n_jobs``, with the exception of ``backend`` which is directly
+              controlled by ``backend``.
+              If ``n_jobs`` is not passed, it will default to ``-1``, other parameters
+              will default to ``joblib`` defaults.
+            - "joblib": custom and 3rd party ``joblib`` backends,
+              e.g., ``spark``. Any valid keys for ``joblib.Parallel``
+              can be passed here, e.g., ``n_jobs``,
+              ``backend`` must be passed as a key of ``backend_params`` in this case.
+              If ``n_jobs`` is not passed, it will default to ``-1``, other parameters
+              will default to ``joblib`` defaults.
+            - "dask": any valid keys for ``dask.compute`` can be passed,
+              e.g., ``scheduler``
+        """,
     }
 
     def __init__(self):
diff --git a/sktime/base/_base_panel.py b/sktime/base/_base_panel.py
new file mode 100644
index 000000000..9cffa2a41
--- /dev/null
+++ b/sktime/base/_base_panel.py
@@ -0,0 +1,461 @@
+# copyright: sktime developers, BSD-3-Clause License (see LICENSE file)
+"""Intermediate base class mixin with common functionality for panel tasks.
+
+Inherits from BaseEstimator, descendents are:
+
+BaseClassifier
+BaseRegressor
+"""
+
+__author__ = ["fkiraly"]
+__all__ = ["BasePanelMixin"]
+
+
+import numpy as np
+import pandas as pd
+
+from sktime.base import BaseEstimator
+from sktime.datatypes import (
+    MTYPE_LIST_PANEL,
+    MTYPE_LIST_TABLE,
+    VectorizedDF,
+    check_is_error_msg,
+    check_is_scitype,
+    convert,
+)
+from sktime.utils.warnings import warn
+
+
+class BasePanelMixin(BaseEstimator):
+    """Abstract base class for time series panel tasks, e.g., classifiers, regressors.
+
+    The class contains boilerplate for checks and data conversions,
+    which are common to multiple descendants, most importantly,
+    BaseClassifier and BaseRegressor.
+    """
+
+    # convenience constant to control which metadata of input data
+    # are regularly retrieved in input checks
+    METADATA_REQ_IN_CHECKS = [
+        "n_instances",
+        "has_nans",
+        "is_univariate",
+        "is_equal_length",
+    ]
+
+    # attribute name where vectorized estimators are stored
+    VECTORIZATION_ATTR = "estimators_"  # e.g., classifiers_, regressors_
+
+    # used in error messages
+    TASK = "panel data tasks"  # e.g., classification, regression
+    EST_TYPE = "estimator"  # e.g., classifier, regressor
+    EST_TYPE_PLURAL = "estimators"  # e.g., classifiers, regressors
+
+    def _vectorize(self, methodname, **kwargs):
+        """Vectorized/iterated loop over methods of BaseClassifier, BaseRegressor.
+
+        Stores one estimator per loop index.
+        """
+        # retrieve data arguments
+        X = kwargs.pop("X", None)
+        y = kwargs.pop("y", None)
+
+        # add some common arguments to kwargs
+        kwargs["rowname_default"] = self.EST_TYPE_PLURAL
+        kwargs["colname_default"] = self.EST_TYPE_PLURAL
+        kwargs["backend"] = self.get_config()["backend:parallel"]
+        kwargs["backend_params"] = self.get_config()["backend:parallel:params"]
+
+        if methodname == "fit":
+            self._yvec = y
+
+            ests_ = self._yvec.vectorize_est(self, method="clone", **kwargs)
+            ests_fit = self._yvec.vectorize_est(
+                ests_,
+                method=methodname,
+                args={"y": y},
+                X=X,
+                **kwargs,
+            )
+            setattr(self, self.VECTORIZATION_ATTR, ests_fit)
+            return self
+        else:  # methodname == "predict" or methodname == "predict_proba":
+            ests_ = getattr(self, self.VECTORIZATION_ATTR)
+            y_pred = self._yvec.vectorize_est(
+                ests_,
+                method=methodname,
+                X=X,
+                **kwargs,
+            )
+            y_pred = pd.DataFrame(
+                {str(i): y_pred[col].values[0] for i, col in enumerate(y_pred.columns)}
+            )
+            return y_pred
+
+    def _fit_predict_boilerplate(self, X, y, cv, change_state, method):
+        """Boilerplate logic for fit_predict and fit_predict_proba."""
+        from sklearn.model_selection import KFold
+
+        if isinstance(cv, int):
+            random_state = getattr(self, "random_state", None)
+            cv = KFold(cv, random_state=random_state, shuffle=True)
+
+        if change_state:
+            self.reset()
+            est = self
+        else:
+            est = self.clone()
+
+        if cv is None:
+            return getattr(est.fit(X, y), method)(X)
+        elif change_state:
+            self.fit(X, y)
+
+        # we now know that cv is an sklearn splitter
+        X, y = self._internal_convert(X, y)
+        X_metadata = self._check_input(
+            X, y, return_metadata=self.METADATA_REQ_IN_CHECKS
+        )
+        X_mtype = X_metadata["mtype"]
+        # Check this estimator can handle characteristics
+        self._check_capabilities(X_metadata)
+
+        # handle single class case
+        if len(self._class_dictionary) == 1:
+            return self._single_class_y_pred(X)
+
+        # Convert data to format easily usable for applying cv
+        if isinstance(X, np.ndarray):
+            X = convert(
+                X,
+                from_type=X_mtype,
+                to_type="numpy3D",
+                as_scitype="Panel",
+                store_behaviour="freeze",
+            )
+        else:
+            X = convert(
+                X,
+                from_type=X_mtype,
+                to_type="nested_univ",
+                as_scitype="Panel",
+                store_behaviour="freeze",
+            )
+
+        if method == "predict_proba":
+            y_pred = np.empty([len(y), len(np.unique(y))])
+        else:
+            y_pred = np.empty_like(y)
+        y_pred[:] = -1
+        if isinstance(X, np.ndarray):
+            for tr_idx, tt_idx in cv.split(X):
+                X_train = X[tr_idx]
+                X_test = X[tt_idx]
+                y_train = y[tr_idx]
+                fitted_est = self.clone().fit(X_train, y_train)
+                y_pred[tt_idx] = getattr(fitted_est, method)(X_test)
+        else:
+            for tr_idx, tt_idx in cv.split(X):
+                X_train = X.iloc[tr_idx]
+                X_test = X.iloc[tt_idx]
+                y_train = y[tr_idx]
+                fitted_est = self.clone().fit(X_train, y_train)
+                y_pred[tt_idx] = getattr(fitted_est, method)(X_test)
+
+        return y_pred
+
+    def _check_convert_X_for_predict(self, X):
+        """Input checks, capability checks, repeated in all predict/score methods.
+
+        Parameters
+        ----------
+        X : any object (to check/convert)
+            should be of a supported Panel mtype or 2D numpy.ndarray
+
+        Returns
+        -------
+        X: an object of a supported Panel mtype, numpy3D if X was a 2D numpy.ndarray
+
+        Raises
+        ------
+        ValueError if X is of invalid input data type, or there is not enough data
+        ValueError if the capabilities in self._tags do not handle the data.
+        """
+        X = self._internal_convert(X)
+        X_metadata = self._check_input(X, return_metadata=self.METADATA_REQ_IN_CHECKS)
+        X_mtype = X_metadata["mtype"]
+        # Check that estimator can handle characteristics
+        self._check_capabilities(X_metadata)
+        # Convert data as dictated by the estimator tags
+        X = self._convert_X(X, X_mtype=X_mtype)
+
+        return X
+
+    def _check_capabilities(self, X_metadata):
+        """Check whether this estimator can handle the data characteristics.
+
+        Parameters
+        ----------
+        missing : boolean, does the data passed to fit contain missing values?
+        multivariate : boolean, does the data passed to fit contain missing values?
+        unequal : boolea, do the time series passed to fit have variable lengths?
+
+        Raises
+        ------
+        ValueError if the capabilities in self._tags do not handle the data.
+        """
+        missing = X_metadata["has_nans"]
+        multivariate = not X_metadata["is_univariate"]
+        unequal = not X_metadata["is_equal_length"]
+
+        allow_multivariate = self.get_tag("capability:multivariate")
+        allow_missing = self.get_tag("capability:missing_values")
+        allow_unequal = self.get_tag("capability:unequal_length")
+
+        self_name = type(self).__name__
+
+        # identify problems, mismatch of capability and inputs
+        problems = []
+        if missing and not allow_missing:
+            problems += ["missing values"]
+        if multivariate and not allow_multivariate:
+            problems += ["multivariate series"]
+        if unequal and not allow_unequal:
+            problems += ["unequal length series"]
+
+        # construct error message
+        problems_and = " and ".join(problems)
+        problems_or = " or ".join(problems)
+        msg = (
+            f"Data seen by {self_name} instance has {problems_and}, "
+            f"but this {self_name} instance cannot handle {problems_or}. "
+            f"Calls with {problems_or} may result in error or unreliable results."
+        )
+
+        # raise exception or warning with message
+        # if self is composite, raise a warning, since passing could be fine
+        #   see discussion in PR 2366 why
+        if len(problems) > 0:
+            if self.is_composite():
+                warn(msg, obj=self)
+            else:
+                raise ValueError(msg)
+
+    def _convert_X(self, X, X_mtype):
+        """Convert equal length series from DataFrame to numpy array or vice versa.
+
+        Parameters
+        ----------
+        X : input data for the estimator, any Panel mtype
+        X_mtype : str, a Panel mtype string, e.g., "pd_multiindex", "numpy3D"
+
+        Returns
+        -------
+        X : input X converted to type in "X_inner_mtype" tag
+            usually a pd.DataFrame (nested) or 3D np.ndarray
+            Checked and possibly converted input data
+        """
+        inner_type = self.get_tag("X_inner_mtype")
+        # convert pd.DataFrame
+        X = convert(
+            X,
+            from_type=X_mtype,
+            to_type=inner_type,
+            as_scitype="Panel",
+        )
+        return X
+
+    def _check_y(self, y=None):
+        """Check and coerce X/y for fit/transform functions.
+
+        Parameters
+        ----------
+        y : pd.DataFrame, pd.Series or np.ndarray
+
+        Returns
+        -------
+        y : object of sktime compatible time series type
+            can be Series, Panel, Hierarchical
+        y_metadata : dict
+            metadata of y, retured by check_is_scitype
+        """
+        if y is None:
+            return None
+
+        capa_multioutput = self.get_tag("capability:multioutput")
+        y_inner_mtype = self.get_tag("y_inner_mtype")
+
+        y_valid, y_msg, y_metadata = check_is_scitype(
+            y, "Table", return_metadata=["is_univariate"]
+        )
+
+        if not y_valid:
+            allowed_msg = (
+                f"In {self.TASK}, y must be of a supported type, "
+                f"for instance 1D or 2D numpy arrays, pd.DataFrame, or pd.Series. "
+                f"Allowed compatible mtype format specifications are:"
+                f" {MTYPE_LIST_TABLE} ."
+            )
+            check_is_error_msg(
+                y_msg, var_name="y", allowed_msg=allowed_msg, raise_exception=True
+            )
+
+        y_uni = y_metadata["is_univariate"]
+        y_mtype = y_metadata["mtype"]
+
+        requires_vectorization = not capa_multioutput and not y_uni
+
+        if requires_vectorization:
+            y_df = convert(
+                y,
+                from_type=y_mtype,
+                to_type="pd_DataFrame_Table",
+                as_scitype="Table",
+                store=self._converter_store_y,
+            )
+            y_vec = VectorizedDF([y_df], iterate_cols=True)
+            return y_vec, y_metadata
+
+        y_inner = convert(
+            y,
+            from_type=y_mtype,
+            to_type=y_inner_mtype,
+            as_scitype="Table",
+            store=self._converter_store_y,
+        )
+
+        return y_inner, y_metadata
+
+    def _convert_output_y(self, y):
+        """Convert output y to original format.
+
+        Parameters
+        ----------
+        y : np.ndarray or pd.DataFrame
+
+        Returns
+        -------
+        y : np.ndarray or pd.DataFrame
+        """
+        # for consistency with legacy behaviour:
+        # output is coerced to numpy1D in case of univariate output
+        if not self._y_metadata["is_univariate"]:
+            output_mtype = self._y_metadata["mtype"]
+            converter_store = self._converter_store_y
+        else:
+            output_mtype = "numpy1D"
+            converter_store = None
+
+        y = convert(
+            y,
+            from_type=self.get_tag("y_inner_mtype"),
+            to_type=output_mtype,
+            as_scitype="Table",
+            store=converter_store,
+            store_behaviour="freeze",
+        )
+        return y
+
+    def _check_input(self, X, y=None, enforce_min_instances=1, return_metadata=True):
+        """Check whether input X and y are valid formats with minimum data.
+
+        Raises a ValueError if the input is not valid.
+
+        Parameters
+        ----------
+        X : check whether conformant with any sktime Panel mtype specification
+        y : check whether a pd.Series or np.array
+        enforce_min_instances : int, optional (default=1)
+            check there are a minimum number of instances.
+        return_metadata : bool, str, or list of str
+            metadata fields to return with X_metadata, input to check_is_scitype
+
+        Returns
+        -------
+        metadata : dict with metadata for X returned by datatypes.check_is_scitype
+
+        Raises
+        ------
+        ValueError
+            If y or X is invalid input data type, or there is not enough data
+        """
+        # Check X is valid input type and recover the data characteristics
+        X_valid, msg, X_metadata = check_is_scitype(
+            X, scitype="Panel", return_metadata=return_metadata
+        )
+        # raise informative error message if X is in wrong format
+        allowed_msg = (
+            f"Allowed scitypes for {self.EST_TYPE_PLURAL} are Panel mtypes, "
+            f"for instance a pandas.DataFrame with MultiIndex and last(-1) "
+            f"level an sktime compatible time index. "
+            f"Allowed compatible mtype format specifications are: {MTYPE_LIST_PANEL} ."
+        )
+        if not X_valid:
+            check_is_error_msg(
+                msg, var_name="X", allowed_msg=allowed_msg, raise_exception=True
+            )
+
+        n_cases = X_metadata["n_instances"]
+        if n_cases < enforce_min_instances:
+            raise ValueError(
+                f"Minimum number of cases required is {enforce_min_instances} but X "
+                f"has : {n_cases}"
+            )
+
+        # Check y if passed
+        if y is not None:
+            # Check y valid input
+            if not isinstance(y, (pd.Series, pd.DataFrame, np.ndarray)):
+                raise ValueError(
+                    "y must be a np.array or a pd.Series or pd.DataFrame, but found ",
+                    f"type: {type(y)}",
+                )
+            # Check matching number of labels
+            n_labels = y.shape[0]
+            if n_cases != n_labels:
+                raise ValueError(
+                    f"Mismatch in number of cases. Number in X = {n_cases} nos in y = "
+                    f"{n_labels}"
+                )
+            if isinstance(y, np.ndarray):
+                if y.ndim > 2:
+                    raise ValueError(
+                        f"np.ndarray y must be 1-dimensional or 2-dimensional, "
+                        f"but found {y.ndim} dimensions"
+                    )
+            # warn if only a single class label is seen
+            # this should not raise exception since this can occur by train subsampling
+            if len(np.unique(y)) == 1:
+                warn(
+                    "only single label seen in y passed to "
+                    f"fit of {self.EST_TYPE} {type(self).__name__}",
+                    obj=self,
+                )
+
+        return X_metadata
+
+    def _internal_convert(self, X, y=None):
+        """Convert X and y if necessary as a user convenience.
+
+        Convert X to a 3D numpy array if already a 2D and convert y into an 1D numpy
+        array if passed as a Series.
+
+        Parameters
+        ----------
+        X : an object of a supported Panel mtype, or 2D numpy.ndarray
+        y : np.ndarray or pd.Series
+
+        Returns
+        -------
+        X: an object of a supported Panel mtype, numpy3D if X was a 2D numpy.ndarray
+        y: np.ndarray
+        """
+        if isinstance(X, np.ndarray):
+            # Temporary fix to insist on 3D numpy. For univariate problems,
+            # most panel estimators simply convert back to 2D. This squeezing should be
+            # done here, but touches a lot of files, so will get this to work first.
+            if X.ndim == 2:
+                X = X.reshape(X.shape[0], 1, X.shape[1])
+        if y is None:
+            return X
+        return X, y
diff --git a/sktime/classification/base.py b/sktime/classification/base.py
index 7d8110f49..8f218969c 100644
--- a/sktime/classification/base.py
+++ b/sktime/classification/base.py
@@ -23,25 +23,17 @@ __all__ = [
 __author__ = ["mloning", "fkiraly", "TonyBagnall", "MatthewMiddlehurst"]
 
 import time
-from abc import ABC, abstractmethod
 
 import numpy as np
-import pandas as pd
-
-from sktime.base import BaseEstimator
-from sktime.datatypes import (
-    MTYPE_LIST_PANEL,
-    check_is_error_msg,
-    check_is_scitype,
-    convert_to,
-)
+
+from sktime.base import BasePanelMixin
+from sktime.datatypes import VectorizedDF, check_is_scitype, convert
 from sktime.utils.sklearn import is_sklearn_transformer
 from sktime.utils.validation import check_n_jobs
 from sktime.utils.validation._dependencies import _check_estimator_deps
-from sktime.utils.warnings import warn
 
 
-class BaseClassifier(BaseEstimator, ABC):
+class BaseClassifier(BasePanelMixin):
     """Abstract base class for time series classifiers.
 
     The base classifier specifies the methods and method signatures that all
@@ -60,7 +52,9 @@ class BaseClassifier(BaseEstimator, ABC):
     _tags = {
         "object_type": "classifier",  # type of object
         "X_inner_mtype": "numpy3D",  # which type do _fit/_predict, support for X?
+        "y_inner_mtype": "numpy1D",  # which type do _fit/_predict, support for y?
         #    it should be either "numpy3D" or "nested_univ" (nested pd.DataFrame)
+        "capability:multioutput": False,  # whether classifier supports multioutput
         "capability:multivariate": False,
         "capability:unequal_length": False,
         "capability:missing_values": False,
@@ -81,6 +75,14 @@ class BaseClassifier(BaseEstimator, ABC):
         "is_equal_length",
     ]
 
+    # attribute name where vectorized estimators are stored
+    VECTORIZATION_ATTR = "classifiers_"  # e.g., classifiers_, regressors_
+
+    # used in error messages
+    TASK = "classification"  # e.g., classification, regression
+    EST_TYPE = "classifier"  # e.g., classifier, regressor
+    EST_TYPE_PLURAL = "classifiers"  # e.g., classifiers, regressors
+
     def __init__(self):
         # reserved attributes written to in fit
         self.classes_ = []  # classes seen in y, unique labels
@@ -93,6 +95,8 @@ class BaseClassifier(BaseEstimator, ABC):
         # required for compatibility with some sklearn interfaces
         # i.e. CalibratedClassifierCV
         self._estimator_type = "classifier"
+        self._is_vectorized = False
+        self._converter_store_y = {}
 
         super().__init__()
         _check_estimator_deps(self)
@@ -136,48 +140,71 @@ class BaseClassifier(BaseEstimator, ABC):
             return NotImplemented
 
     def fit(self, X, y):
-        """Fit time series classifier to training data.
+        """
+        Fit time series classifier to training data.
+
+        State change:
+            Changes state to "fitted".
+
+        Writes to self:
+            Sets self.is_fitted to True.
+            Sets fitted model attributes ending in "_".
 
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of int, of shape [n_instances] - class labels for fitting
-            indices correspond to instance indices in X
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
+        y : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            supported sktime types: np.ndarray (1D, 2D), pd.Series, pd.DataFrame
 
         Returns
         -------
         self : Reference to self.
-
-        Notes
-        -----
-        Changes state by creating a fitted model that updates attributes
-        ending in "_" and sets is_fitted flag to True.
         """
         # reset estimator at the start of fit
         self.reset()
 
+        # fit timer start
         start = int(round(time.time() * 1000))
+
+        # check and convert y for multioutput vectorization
+        y, y_metadata = self._check_y(y)
+        self._y_metadata = y_metadata
+        self._is_vectorized = isinstance(y, VectorizedDF)
+
+        if self._is_vectorized:
+            self._vectorize("fit", X=X, y=y)
+            # fit timer end
+            self.fit_time_ = int(round(time.time() * 1000)) - start
+            # this should happen last: fitted state is set to True
+            self._is_fitted = True
+            return self
+
+        # no vectorization needed, proceed with normal fit
+
         # convenience conversions to allow user flexibility:
         # if X is 2D array, convert to 3D, if y is Series, convert to numpy
         X, y = self._internal_convert(X, y)
-        X_metadata = self._check_classifier_input(
+        X_metadata = self._check_input(
             X, y, return_metadata=self.METADATA_REQ_IN_CHECKS
         )
-        missing = X_metadata["has_nans"]
-        multivariate = not X_metadata["is_univariate"]
-        unequal = not X_metadata["is_equal_length"]
+        X_mtype = X_metadata["mtype"]
         self._X_metadata = X_metadata
 
         # Check this classifier can handle characteristics
-        self._check_capabilities(missing, multivariate, unequal)
+        self._check_capabilities(X_metadata)
 
         # remember class labels
         self.classes_ = np.unique(y)
@@ -194,7 +221,7 @@ class BaseClassifier(BaseEstimator, ABC):
             return self
 
         # Convert data as dictated by the classifier tags
-        X = self._convert_X(X)
+        X = self._convert_X(X, X_mtype)
         multithread = self.get_tag("capability:multithreading")
         if multithread:
             try:
@@ -208,32 +235,42 @@ class BaseClassifier(BaseEstimator, ABC):
         self._fit(X, y)
         self.fit_time_ = int(round(time.time() * 1000)) - start
 
-        # this should happen last
+        # this should happen last: fitted state is set to True
         self._is_fitted = True
         return self
 
-    def predict(self, X) -> np.ndarray:
+    def predict(self, X):
         """Predicts labels for sequences in X.
 
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
 
         Returns
         -------
-        y : 1D np.array of int, of shape [n_instances] - predicted class labels
-            indices correspond to instance indices in X
+        y_pred : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            predicted class labels
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            1D np.npdarray, if y univariate (one dimension)
+            otherwise, same type as y passed in fit
         """
         self.check_is_fitted()
 
+        # enter vectorized mode if needed
+        if self._is_vectorized:
+            return self._vectorize("predict", X=X)
+
         # boilerplate input checks for predict-like methods
         X = self._check_convert_X_for_predict(X)
 
@@ -241,33 +278,42 @@ class BaseClassifier(BaseEstimator, ABC):
         if len(self._class_dictionary) == 1:
             return self._single_class_y_pred(X, method="predict")
 
-        # call internal _predict_proba
-        return self._predict(X)
+        # call internal _predict, convert output
+        y_pred_inner = self._predict(X)
+        y_pred = self._convert_output_y(y_pred_inner)
+        return y_pred
 
-    def predict_proba(self, X) -> np.ndarray:
+    def predict_proba(self, X):
         """Predicts labels probabilities for sequences in X.
 
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
 
         Returns
         -------
-        y : 2D array of shape [n_instances, n_classes] - predicted class probabilities
-            1st dimension indices correspond to instance indices in X
-            2nd dimension indices correspond to possible labels (integers)
-            (i, j)-th entry is predictive probability that i-th instance is of class j
+        y_pred : 2D np.array of int, of shape [n_instances, n_classes]
+            predicted class label probabilities
+            0-th indices correspond to instance indices in X
+            1-st indices correspond to class index, in same order as in self.classes_
+            entries are predictive class probabilities, summing to 1
         """
         self.check_is_fitted()
 
+        # enter vectorized mode if needed
+        if self._is_vectorized:
+            return self._vectorize("predict_proba", X=X)
+
+        self.check_is_fitted()
+
         # boilerplate input checks for predict-like methods
         X = self._check_convert_X_for_predict(X)
 
@@ -278,7 +324,7 @@ class BaseClassifier(BaseEstimator, ABC):
         # call internal _predict_proba
         return self._predict_proba(X)
 
-    def fit_predict(self, X, y, cv=None, change_state=True) -> np.ndarray:
+    def fit_predict(self, X, y, cv=None, change_state=True):
         """Fit and predict labels for sequences in X.
 
         Convenience method to produce in-sample predictions and
@@ -292,17 +338,22 @@ class BaseClassifier(BaseEstimator, ABC):
 
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of int, of shape [n_instances] - class labels for fitting
-            indices correspond to instance indices in X
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
+        y : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            supported sktime types: np.ndarray (1D, 2D), pd.Series, pd.DataFrame
         cv : None, int, or sklearn cross-validation object, optional, default=None
             None : predictions are in-sample, equivalent to fit(X, y).predict(X)
             cv : predictions are equivalent to fit(X_train, y_train).predict(X_test)
@@ -320,9 +371,14 @@ class BaseClassifier(BaseEstimator, ABC):
 
         Returns
         -------
-        y : 1D np.array of int, of shape [n_instances] - predicted class labels
-            indices correspond to instance indices in X
-            if cv is passed, -1 indicates entries not seen in union of test sets
+        y_pred : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            predicted class labels
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            1D np.npdarray, if y univariate (one dimension)
+            otherwise, same type as y passed in fit
         """
         return self._fit_predict_boilerplate(
             X=X, y=y, cv=cv, change_state=change_state, method="predict"
@@ -349,14 +405,12 @@ class BaseClassifier(BaseEstimator, ABC):
 
         # we now know that cv is an sklearn splitter
         X, y = self._internal_convert(X, y)
-        X_metadata = self._check_classifier_input(
+        X_metadata = self._check_input(
             X, y, return_metadata=self.METADATA_REQ_IN_CHECKS
         )
-        missing = X_metadata["has_nans"]
-        multivariate = not X_metadata["is_univariate"]
-        unequal = not X_metadata["is_equal_length"]
+        X_mtype = X_metadata["mtype"]
         # Check this classifier can handle characteristics
-        self._check_capabilities(missing, multivariate, unequal)
+        self._check_capabilities(X_metadata)
 
         # handle single class case
         if len(self._class_dictionary) == 1:
@@ -364,15 +418,17 @@ class BaseClassifier(BaseEstimator, ABC):
 
         # Convert data to format easily usable for applying cv
         if isinstance(X, np.ndarray):
-            X = convert_to(
+            X = convert(
                 X,
+                from_type=X_mtype,
                 to_type="numpy3D",
                 as_scitype="Panel",
                 store_behaviour="freeze",
             )
         else:
-            X = convert_to(
+            X = convert(
                 X,
+                from_type=X_mtype,
                 to_type="nested_univ",
                 as_scitype="Panel",
                 store_behaviour="freeze",
@@ -400,25 +456,36 @@ class BaseClassifier(BaseEstimator, ABC):
 
         return y_pred
 
-    def fit_predict_proba(self, X, y, cv=None, change_state=True) -> np.ndarray:
+    def fit_predict_proba(self, X, y, cv=None, change_state=True):
         """Fit and predict labels probabilities for sequences in X.
 
         Convenience method to produce in-sample predictions and
         cross-validated out-of-sample predictions.
 
+        Writes to self, if change_state=True:
+            Sets self.is_fitted to True.
+            Sets fitted model attributes ending in "_".
+
+        Does not update state if change_state=False.
+
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of int, of shape [n_instances] - class labels for fitting
-            indices correspond to instance indices in X
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
+        y : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            supported sktime types: np.ndarray (1D, 2D), pd.Series, pd.DataFrame
         cv : None, int, or sklearn cross-validation object, optional, default=None
             None : predictions are in-sample, equivalent to fit(X, y).predict(X)
             cv : predictions are equivalent to fit(X_train, y_train).predict(X_test)
@@ -434,10 +501,11 @@ class BaseClassifier(BaseEstimator, ABC):
 
         Returns
         -------
-        y : 2D array of shape [n_instances, n_classes] - predicted class probabilities
-            1st dimension indices correspond to instance indices in X
-            2nd dimension indices correspond to possible labels (integers)
-            (i, j)-th entry is predictive probability that i-th instance is of class j
+        y_pred : 2D np.array of int, of shape [n_instances, n_classes]
+            predicted class label probabilities
+            0-th indices correspond to instance indices in X
+            1-st indices correspond to class index, in same order as in self.classes_
+            entries are predictive class probabilities, summing to 1
         """
         return self._fit_predict_boilerplate(
             X=X, y=y, cv=cv, change_state=change_state, method="predict_proba"
@@ -460,17 +528,22 @@ class BaseClassifier(BaseEstimator, ABC):
 
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.ndarray of int, of shape [n_instances] - class labels (ground truth)
-            indices correspond to instance indices in X
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
+        y : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            supported sktime types: np.ndarray (1D, 2D), pd.Series, pd.DataFrame
 
         Returns
         -------
@@ -506,35 +579,37 @@ class BaseClassifier(BaseEstimator, ABC):
         """
         return super().get_test_params(parameter_set=parameter_set)
 
-    @abstractmethod
     def _fit(self, X, y):
         """Fit time series classifier to training data.
 
         Abstract method, must be implemented.
 
+        Writes to self:
+            Sets fitted model attributes ending in "_".
+
         Parameters
         ----------
         X : guaranteed to be of a type in self.get_tag("X_inner_mtype")
             if self.get_tag("X_inner_mtype") = "numpy3D":
-                3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
-            if self.get_tag("X_inner_mtype") = "nested_univ":
-                pd.DataFrame with each column a dimension, each cell a pd.Series
+            3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            if self.get_tag("X_inner_mtype") = "pd-multiindex:":
+            pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
             for list of other mtypes, see datatypes.SCITYPE_REGISTER
             for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of int, of shape [n_instances] - class labels for fitting
-            indices correspond to instance indices in X
+        y : guaranteed to be of a type in self.get_tag("y_inner_mtype")
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            if self.get_tag("capaility:multioutput") = False, guaranteed to be 1D
+            if self.get_tag("capaility:multioutput") = True, guaranteed to be 2D
 
         Returns
         -------
-        self :
-            Reference to self.
-
-        Notes
-        -----
-        Changes state by creating a fitted model that updates attributes
-        ending in "_" and sets is_fitted flag to True.
+        self : Reference to self.
         """
-        ...
+        raise NotImplementedError("abstract method")
 
     def _predict(self, X) -> np.ndarray:
         """Predicts labels for sequences in X.
@@ -545,16 +620,23 @@ class BaseClassifier(BaseEstimator, ABC):
         ----------
         X : guaranteed to be of a type in self.get_tag("X_inner_mtype")
             if self.get_tag("X_inner_mtype") = "numpy3D":
-                3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
-            if self.get_tag("X_inner_mtype") = "nested_univ":
-                pd.DataFrame with each column a dimension, each cell a pd.Series
+            3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            if self.get_tag("X_inner_mtype") = "pd-multiindex:":
+            pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
             for list of other mtypes, see datatypes.SCITYPE_REGISTER
             for specifications, see examples/AA_datatypes_and_datasets.ipynb
 
         Returns
         -------
-        y : 1D np.array of int, of shape [n_instances] - predicted class labels
+        y : should be of mtype in self.get_tag("y_inner_mtype")
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            predicted class labels
             indices correspond to instance indices in X
+            if self.get_tag("capaility:multioutput") = False, should be 1D
+            if self.get_tag("capaility:multioutput") = True, should be 2D
         """
         y_proba = self._predict_proba(X)
         y_pred = y_proba.argmax(axis=1)
@@ -572,9 +654,11 @@ class BaseClassifier(BaseEstimator, ABC):
         ----------
         X : guaranteed to be of a type in self.get_tag("X_inner_mtype")
             if self.get_tag("X_inner_mtype") = "numpy3D":
-                3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
-            if self.get_tag("X_inner_mtype") = "nested_univ":
-                pd.DataFrame with each column a dimension, each cell a pd.Series
+            3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            if self.get_tag("X_inner_mtype") = "pd-multiindex:":
+            pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
             for list of other mtypes, see datatypes.SCITYPE_REGISTER
             for specifications, see examples/AA_datatypes_and_datasets.ipynb
 
@@ -592,211 +676,3 @@ class BaseClassifier(BaseEstimator, ABC):
             dists[i, self._class_dictionary[preds[i]]] = 1
 
         return dists
-
-    def _check_convert_X_for_predict(self, X):
-        """Input checks, capability checks, repeated in all predict/score methods.
-
-        Parameters
-        ----------
-        X : any object (to check/convert)
-            should be of a supported Panel mtype or 2D numpy.ndarray
-
-        Returns
-        -------
-        X: an object of a supported Panel mtype, numpy3D if X was a 2D numpy.ndarray
-
-        Raises
-        ------
-        ValueError if X is of invalid input data type, or there is not enough data
-        ValueError if the capabilities in self._tags do not handle the data.
-        """
-        X = self._internal_convert(X)
-        X_metadata = self._check_classifier_input(
-            X, return_metadata=self.METADATA_REQ_IN_CHECKS
-        )
-        missing = X_metadata["has_nans"]
-        multivariate = not X_metadata["is_univariate"]
-        unequal = not X_metadata["is_equal_length"]
-        # Check this classifier can handle characteristics
-        self._check_capabilities(missing, multivariate, unequal)
-        # Convert data as dictated by the classifier tags
-        X = self._convert_X(X)
-
-        return X
-
-    def _check_capabilities(self, missing, multivariate, unequal):
-        """Check whether this classifier can handle the data characteristics.
-
-        Parameters
-        ----------
-        missing : boolean, does the data passed to fit contain missing values?
-        multivariate : boolean, does the data passed to fit contain missing values?
-        unequal : boolea, do the time series passed to fit have variable lengths?
-
-        Raises
-        ------
-        ValueError if the capabilities in self._tags do not handle the data.
-        """
-        allow_multivariate = self.get_tag("capability:multivariate")
-        allow_missing = self.get_tag("capability:missing_values")
-        allow_unequal = self.get_tag("capability:unequal_length")
-
-        self_name = type(self).__name__
-
-        # identify problems, mismatch of capability and inputs
-        problems = []
-        if missing and not allow_missing:
-            problems += ["missing values"]
-        if multivariate and not allow_multivariate:
-            problems += ["multivariate series"]
-        if unequal and not allow_unequal:
-            problems += ["unequal length series"]
-
-        # construct error message
-        problems_and = " and ".join(problems)
-        problems_or = " or ".join(problems)
-        msg = (
-            f"Data seen by {self_name} instance has {problems_and}, "
-            f"but this {self_name} instance cannot handle {problems_or}. "
-            f"Calls with {problems_or} may result in error or unreliable results."
-        )
-
-        # raise exception or warning with message
-        # if self is composite, raise a warning, since passing could be fine
-        #   see discussion in PR 2366 why
-        if len(problems) > 0:
-            if self.is_composite():
-                warn(msg, obj=self)
-            else:
-                raise ValueError(msg)
-
-    def _convert_X(self, X):
-        """Convert equal length series from DataFrame to numpy array or vice versa.
-
-        Parameters
-        ----------
-        self : this classifier
-        X : pd.DataFrame or np.ndarray. Input attribute data
-
-        Returns
-        -------
-        X : input X converted to type in "X_inner_mtype" tag
-                usually a pd.DataFrame (nested) or 3D np.ndarray
-            Checked and possibly converted input data
-        """
-        inner_type = self.get_tag("X_inner_mtype")
-        # convert pd.DataFrame
-        X = convert_to(
-            X,
-            to_type=inner_type,
-            as_scitype="Panel",
-        )
-        return X
-
-    def _check_classifier_input(
-        self, X, y=None, enforce_min_instances=1, return_metadata=True
-    ):
-        """Check whether input X and y are valid formats with minimum data.
-
-        Raises a ValueError if the input is not valid.
-
-        Parameters
-        ----------
-        X : check whether conformant with any sktime Panel mtype specification
-        y : check whether a pd.Series or np.array
-        enforce_min_instances : int, optional (default=1)
-            check there are a minimum number of instances.
-        return_metadata : bool, str, or list of str
-            metadata fields to return with X_metadata, input to check_is_scitype
-
-        Returns
-        -------
-        metadata : dict with metadata for X returned by datatypes.check_is_scitype
-
-        Raises
-        ------
-        ValueError
-            If y or X is invalid input data type, or there is not enough data
-        """
-        # Check X is valid input type and recover the data characteristics
-        X_valid, msg, X_metadata = check_is_scitype(
-            X, scitype="Panel", return_metadata=return_metadata
-        )
-        # raise informative error message if X is in wrong format
-        allowed_msg = (
-            f"Allowed scitypes for classifiers are Panel mtypes, "
-            f"for instance a pandas.DataFrame with MultiIndex and last(-1) "
-            f"level an sktime compatible time index. "
-            f"Allowed compatible mtype format specifications are: {MTYPE_LIST_PANEL} ."
-        )
-        if not X_valid:
-            check_is_error_msg(
-                msg, var_name="X", allowed_msg=allowed_msg, raise_exception=True
-            )
-
-        n_cases = X_metadata["n_instances"]
-        if n_cases < enforce_min_instances:
-            raise ValueError(
-                f"Minimum number of cases required is {enforce_min_instances} but X "
-                f"has : {n_cases}"
-            )
-
-        # Check y if passed
-        if y is not None:
-            # Check y valid input
-            if not isinstance(y, (pd.Series, np.ndarray)):
-                raise ValueError(
-                    f"y must be a np.array or a pd.Series, but found type: {type(y)}"
-                )
-            # Check matching number of labels
-            n_labels = y.shape[0]
-            if n_cases != n_labels:
-                raise ValueError(
-                    f"Mismatch in number of cases. Number in X = {n_cases} nos in y = "
-                    f"{n_labels}"
-                )
-            if isinstance(y, np.ndarray):
-                if y.ndim > 1:
-                    raise ValueError(
-                        f"np.ndarray y must be 1-dimensional, "
-                        f"but found {y.ndim} dimensions"
-                    )
-            # warn if only a single class label is seen
-            # this should not raise exception since this can occur by train subsampling
-            if len(np.unique(y)) == 1:
-                warn(
-                    "only single class label seen in y passed to "
-                    f"fit of classifier {type(self).__name__}",
-                    obj=self,
-                )
-
-        return X_metadata
-
-    def _internal_convert(self, X, y=None):
-        """Convert X and y if necessary as a user convenience.
-
-        Convert X to a 3D numpy array if already a 2D and convert y into an 1D numpy
-        array if passed as a Series.
-
-        Parameters
-        ----------
-        X : an object of a supported Panel mtype, or 2D numpy.ndarray
-        y : np.ndarray or pd.Series
-
-        Returns
-        -------
-        X: an object of a supported Panel mtype, numpy3D if X was a 2D numpy.ndarray
-        y: np.ndarray
-        """
-        if isinstance(X, np.ndarray):
-            # Temporary fix to insist on 3D numpy. For univariate problems,
-            # most classifiers simply convert back to 2D. This squeezing should be
-            # done here, but touches a lot of files, so will get this to work first.
-            if X.ndim == 2:
-                X = X.reshape(X.shape[0], 1, X.shape[1])
-        if y is not None and isinstance(y, pd.Series):
-            # y should be a numpy array, although we allow Series for user convenience
-            y = pd.Series.to_numpy(y)
-        if y is None:
-            return X
-        return X, y
diff --git a/sktime/classification/deep_learning/rnn.py b/sktime/classification/deep_learning/rnn.py
index cbac0badb..369ab4c0e 100644
--- a/sktime/classification/deep_learning/rnn.py
+++ b/sktime/classification/deep_learning/rnn.py
@@ -11,8 +11,13 @@ from sklearn.utils import check_random_state
 from sktime.classification.deep_learning.base import BaseDeepClassifier
 from sktime.networks.rnn import RNNNetwork
 from sktime.utils.validation._dependencies import _check_dl_dependencies
+from sktime.utils.warnings import warn
 
 
+# todo: 0.26.0 - please remove num_epochs parameter and related logic,
+# because parameter num_epochs usage is deprecated and will be renamed
+# to n_epochs.  Also n_epochs should be moved to become the first
+# argument in the __init__ and super().__init__.
 class SimpleRNNClassifier(BaseDeepClassifier):
     """Simple recurrent neural network.
 
@@ -51,7 +56,7 @@ class SimpleRNNClassifier(BaseDeepClassifier):
 
     def __init__(
         self,
-        num_epochs=100,
+        num_epochs=None,
         batch_size=1,
         units=6,
         callbacks=None,
@@ -63,10 +68,26 @@ class SimpleRNNClassifier(BaseDeepClassifier):
         activation="sigmoid",
         use_bias=True,
         optimizer=None,
+        n_epochs=100,
     ):
         _check_dl_dependencies(severity="error")
         super().__init__()
+        # todo: 0.26.0 - remove this, replace by
+        # self._n_epochs = n_epochs
+        # Deprecated Parameter Handling
+        if num_epochs is not None:
+            warn(
+                "In SimpleRNNClassifier, the parameter 'num_epochs' is deprecated "
+                "and will be removed in v0.26.0. It will be renamed to n_epochs. "
+                "To avoid this warning, update your code to use 'n_epochs'. ",
+                DeprecationWarning,
+                stacklevel=2,
+            )
+            self._n_epochs = num_epochs
+        else:
+            self._n_epochs = n_epochs
         self.num_epochs = num_epochs
+        # end remove
         self.batch_size = batch_size
         self.verbose = verbose
         self.units = units
@@ -80,6 +101,7 @@ class SimpleRNNClassifier(BaseDeepClassifier):
         self.optimizer = optimizer
         self.history = None
         self._network = RNNNetwork(random_state=random_state, units=units)
+        self.n_epochs = n_epochs
 
     def build_model(self, input_shape, n_classes, **kwargs):
         """Construct a compiled, un-trained, keras model that is ready for training.
@@ -189,7 +211,7 @@ class SimpleRNNClassifier(BaseDeepClassifier):
             X,
             y_onehot,
             batch_size=self.batch_size,
-            epochs=self.num_epochs,
+            epochs=self._n_epochs,
             verbose=self.verbose,
             callbacks=self.callbacks_,
         )
@@ -218,7 +240,7 @@ class SimpleRNNClassifier(BaseDeepClassifier):
         """
         params1 = {}
         params2 = {
-            "num_epochs": 50,
+            "n_epochs": 50,
             "batch_size": 2,
             "units": 5,
             "use_bias": False,
diff --git a/sktime/classification/dictionary_based/_boss.py b/sktime/classification/dictionary_based/_boss.py
index cfcfb6ca9..88889a05a 100644
--- a/sktime/classification/dictionary_based/_boss.py
+++ b/sktime/classification/dictionary_based/_boss.py
@@ -666,6 +666,7 @@ class IndividualBOSS(BaseClassifier):
         new_boss.n_classes_ = self.n_classes_
         new_boss.classes_ = self.classes_
         new_boss._class_dictionary = self._class_dictionary
+        new_boss._y_metadata = self._y_metadata
         new_boss._is_fitted = True
 
         return new_boss
diff --git a/sktime/classification/early_classification/_probability_threshold.py b/sktime/classification/early_classification/_probability_threshold.py
index 0464bfd3e..4697bf9de 100644
--- a/sktime/classification/early_classification/_probability_threshold.py
+++ b/sktime/classification/early_classification/_probability_threshold.py
@@ -20,7 +20,7 @@ from sktime.classification.interval_based import CanonicalIntervalForest
 from sktime.utils.validation.panel import check_X
 
 
-# TODO: fix this in 0.25.0
+# TODO: fix this in 0.26.0
 # base class should have been changed to BaseEarlyClassifier
 class ProbabilityThresholdEarlyClassifier(BaseClassifier):
     """Probability Threshold Early Classifier.
diff --git a/sktime/classification/early_classification/base.py b/sktime/classification/early_classification/base.py
index 335cb38b7..f0702c5af 100644
--- a/sktime/classification/early_classification/base.py
+++ b/sktime/classification/early_classification/base.py
@@ -55,7 +55,9 @@ class BaseEarlyClassifier(BaseEstimator, ABC):
     _tags = {
         "object_type": "early_classifier",  # type of object
         "X_inner_mtype": "numpy3D",  # which type do _fit/_predict, support for X?
+        "y_inner_mtype": "numpy1D",  # which type do _fit/_predict, support for y?
         #    it should be either "numpy3D" or "nested_univ" (nested pd.DataFrame)
+        "capability:multioutput": False,  # whether classifier supports multioutput
         "capability:multivariate": False,
         "capability:unequal_length": False,
         "capability:missing_values": False,
@@ -71,6 +73,14 @@ class BaseEarlyClassifier(BaseEstimator, ABC):
         "is_equal_length",
     ]
 
+    # attribute name where vectorized estimators are stored
+    VECTORIZATION_ATTR = "classifiers_"  # e.g., classifiers_, regressors_
+
+    # used in error messages
+    TASK = "early classification"  # e.g., classification, regression
+    EST_TYPE = "early classifier"  # e.g., classifier, regressor
+    EST_TYPE_PLURAL = "early classifiers"  # e.g., classifiers, regressors
+
     def __init__(self):
         self.classes_ = []
         self.n_classes_ = 0
@@ -89,6 +99,8 @@ class BaseEarlyClassifier(BaseEstimator, ABC):
         """
         self.state_info = None
 
+        self._converter_store_y = {}
+
         super().__init__()
 
     def fit(self, X, y):
@@ -587,42 +599,53 @@ class BaseEarlyClassifier(BaseEstimator, ABC):
         _check_convert_X_for_predict = BaseClassifier._check_convert_X_for_predict
         return _check_convert_X_for_predict(self, X)
 
-    def _check_capabilities(self, missing, multivariate, unequal):
+    def _check_capabilities(self, X_metadata):
         """Check whether this classifier can handle the data characteristics.
 
         Parameters
         ----------
-        missing : boolean, does the data passed to fit contain missing values?
-        multivariate : boolean, does the data passed to fit contain missing values?
-        unequal : boolea, do the time series passed to fit have variable lengths?
+        X_metadata : dict with metadata for X returned by datatypes.check_is_scitype
 
         Raises
         ------
         ValueError if the capabilities in self._tags do not handle the data.
         """
         _check_capabilities = BaseClassifier._check_capabilities
-        return _check_capabilities(self, missing, multivariate, unequal)
+        return _check_capabilities(self, X_metadata)
 
-    def _convert_X(self, X):
+    def _convert_X(self, X, X_mtype):
         """Convert equal length series from DataFrame to numpy array or vice versa.
 
         Parameters
         ----------
-        self : this classifier
-        X : pd.DataFrame or np.ndarray. Input attribute data
+        X : input data for the classifier
+        X_mtype : str, a Panel mtype string, e.g., "pd_multiindex", "numpy3D"
 
         Returns
         -------
         X : input X converted to type in "X_inner_mtype" tag
-                usually a pd.DataFrame (nested) or 3D np.ndarray
+            usually a pd.DataFrame (nested) or 3D np.ndarray
             Checked and possibly converted input data
         """
         _convert_X = BaseClassifier._convert_X
-        return _convert_X(self, X)
+        return _convert_X(self, X, X_mtype)
+
+    def _check_y(self, y=None):
+        """Check and coerce X/y for fit/transform functions.
+
+        Parameters
+        ----------
+        y : pd.DataFrame, pd.Series or np.ndarray
+
+        Returns
+        -------
+        y : object of sktime compatible time series type
+            can be Series, Panel, Hierarchical
+        """
+        _check_y = BaseClassifier._check_y
+        return _check_y(self, y)
 
-    def _check_classifier_input(
-        self, X, y=None, enforce_min_instances=1, return_metadata=True
-    ):
+    def _check_input(self, X, y=None, enforce_min_instances=1, return_metadata=True):
         """Check whether input X and y are valid formats with minimum data.
 
         Raises a ValueError if the input is not valid.
@@ -645,10 +668,8 @@ class BaseEarlyClassifier(BaseEstimator, ABC):
         ValueError
             If y or X is invalid input data type, or there is not enough data
         """
-        _check_classifier_input = BaseClassifier._check_classifier_input
-        return _check_classifier_input(
-            self, X, y, enforce_min_instances, return_metadata
-        )
+        _check_input = BaseClassifier._check_input
+        return _check_input(self, X, y, enforce_min_instances, return_metadata)
 
     def _internal_convert(self, X, y=None):
         """Convert X and y if necessary as a user convenience.
diff --git a/sktime/datatypes/_table/_check.py b/sktime/datatypes/_table/_check.py
index 6b25cf333..4a1d96d88 100644
--- a/sktime/datatypes/_table/_check.py
+++ b/sktime/datatypes/_table/_check.py
@@ -65,11 +65,6 @@ def check_pddataframe_table(obj, return_metadata=False, var_name="obj"):
     if _req("has_nans", return_metadata):
         metadata["has_nans"] = obj.isna().values.any()
 
-    # check that no dtype is object
-    if "object" in obj.dtypes.values:
-        msg = f"{var_name} should not have column of 'object' dtype"
-        return _ret(False, msg, None, return_metadata)
-
     return _ret(True, None, metadata, return_metadata)
 
 
@@ -92,11 +87,6 @@ def check_pdseries_table(obj, return_metadata=False, var_name="obj"):
     if _req("n_instances", return_metadata):
         metadata["n_instances"] = len(index)
 
-    # check that dtype is not object
-    if "object" == obj.dtypes:
-        msg = f"{var_name} should not be of 'object' dtype"
-        return _ret(False, msg, None, return_metadata)
-
     # check whether index is equally spaced or if there are any nans
     #   compute only if needed
     if _req("has_nans", return_metadata):
diff --git a/sktime/forecasting/all/__init__.py b/sktime/forecasting/all/__init__.py
index 2d72d78b8..7ac4e1a8e 100644
--- a/sktime/forecasting/all/__init__.py
+++ b/sktime/forecasting/all/__init__.py
@@ -11,13 +11,6 @@ import pandas as pd
 from sktime.datasets import load_airline, load_longley, load_lynx, load_shampoo_sales
 from sktime.forecasting.base import ForecastingHorizon
 from sktime.forecasting.model_evaluation import evaluate
-from sktime.forecasting.model_selection import (
-    CutoffSplitter,
-    ExpandingWindowSplitter,
-    SingleWindowSplitter,
-    SlidingWindowSplitter,
-    temporal_train_test_split,
-)
 from sktime.performance_metrics.forecasting import (
     GeometricMeanRelativeAbsoluteError,
     GeometricMeanRelativeSquaredError,
@@ -58,6 +51,13 @@ from sktime.performance_metrics.forecasting import (
     relative_loss,
 )
 from sktime.registry import all_estimators
+from sktime.split import (
+    CutoffSplitter,
+    ExpandingWindowSplitter,
+    SingleWindowSplitter,
+    SlidingWindowSplitter,
+    temporal_train_test_split,
+)
 from sktime.transformations.series.detrend import Deseasonalizer, Detrender
 from sktime.utils.plotting import plot_series
 
diff --git a/sktime/forecasting/base/_base.py b/sktime/forecasting/base/_base.py
index ed0d26ff9..f22c39e6b 100644
--- a/sktime/forecasting/base/_base.py
+++ b/sktime/forecasting/base/_base.py
@@ -100,48 +100,17 @@ class BaseForecaster(BaseEstimator):
     }
 
     # configs and default config values
+    # see set_config documentation for details
     _config = {
         "backend:parallel": None,  # parallelization backend for broadcasting
         #  {None, "dask", "loky", "multiprocessing", "threading"}
         #  None: no parallelization
         #  "loky", "multiprocessing" and "threading": uses `joblib` Parallel loops
+        #  "joblib": uses custom joblib backend, set via `joblib_backend` tag
         #  "dask": uses `dask`, requires `dask` package in environment
         "backend:parallel:params": None,  # params for parallelization backend
     }
 
-    _config_doc = {
-        "backend:parallel": """
-        backend:parallel : str, optional, default="None"
-            backend to use for parallelization when broadcasting/vectorizing, one of
-
-            - "None": executes loop sequentally, simple list comprehension
-            - "loky", "multiprocessing" and "threading": uses ``joblib.Parallel``
-            - "joblib": custom and 3rd party ``joblib`` backends, e.g., ``spark``
-            - "dask": uses ``dask``, requires ``dask`` package in environment
-        """,
-        "backend:parallel:params": """
-        backend:parallel:params : dict, optional, default={} (no parameters passed)
-            additional parameters passed to the parallelization backend as config.
-            Valid keys depend on the value of ``backend:parallel``:
-
-            - "None": no additional parameters, ``backend_params`` is ignored
-            - "loky", "multiprocessing" and "threading": default ``joblib`` backends
-              any valid keys for ``joblib.Parallel`` can be passed here, e.g.,
-              ``n_jobs``, with the exception of ``backend`` which is directly
-              controlled by ``backend``.
-              If ``n_jobs`` is not passed, it will default to ``-1``, other parameters
-              will default to ``joblib`` defaults.
-            - "joblib": custom and 3rd party ``joblib`` backends,
-              e.g., ``spark``. Any valid keys for ``joblib.Parallel``
-              can be passed here, e.g., ``n_jobs``,
-            ``backend`` must be passed as a key of ``backend_params`` in this case.
-              If ``n_jobs`` is not passed, it will default to ``-1``, other parameters
-              will default to ``joblib`` defaults.
-            - "dask": any valid keys for ``dask.compute`` can be passed,
-              e.g., ``scheduler``
-        """,
-    }
-
     def __init__(self):
         self._is_fitted = False
 
diff --git a/sktime/forecasting/base/_fh.py b/sktime/forecasting/base/_fh.py
index ff5341c96..7a71c6e9a 100644
--- a/sktime/forecasting/base/_fh.py
+++ b/sktime/forecasting/base/_fh.py
@@ -806,7 +806,7 @@ def _to_relative(fh: ForecastingHorizon, cutoff=None) -> ForecastingHorizon:
             absolute = _coerce_to_period(absolute, freq=fh.freq)
             cutoff = _coerce_to_period(cutoff, freq=fh.freq)
 
-        # TODO: 0.25.0:
+        # TODO: 0.26.0:
         # Check at every minor release whether lower pandas bound >=0.15.0
         # if yes, can remove the workaround in the "else" condition and the check
         #
diff --git a/sktime/forecasting/base/adapters/_statsmodels.py b/sktime/forecasting/base/adapters/_statsmodels.py
index 0e3d91ecb..2d11bc4b0 100644
--- a/sktime/forecasting/base/adapters/_statsmodels.py
+++ b/sktime/forecasting/base/adapters/_statsmodels.py
@@ -2,7 +2,7 @@
 # copyright: sktime developers, BSD-3-Clause License (see LICENSE file)
 """Implements adapter for statsmodels forecasters to be used in sktime framework."""
 
-__author__ = ["mloning"]
+__author__ = ["mloning", "ciaran-g"]
 __all__ = ["_StatsModelsAdapter"]
 
 import inspect
@@ -99,7 +99,7 @@ class _StatsModelsAdapter(BaseForecaster):
         # statsmodels requires zero-based indexing starting at the
         # beginning of the training series when passing integers
         start, end = fh.to_absolute_int(self._y.index[0], self.cutoff)[[0, -1]]
-        fh_abs = fh.to_absolute_index(self.cutoff)
+        fh_int = fh.to_absolute_int(self._y.index[0], self.cutoff) - len(self._y)
 
         # bug fix for evaluate function as test_plus_train indices are passed
         # statsmodels exog must contain test indices only.
@@ -108,7 +108,7 @@ class _StatsModelsAdapter(BaseForecaster):
             ind_drop = self._X.index
             X = X.loc[~X.index.isin(ind_drop)]
             # Entire range of the forecast horizon is required
-            X = X[: fh_abs[-1]]
+            X = X.iloc[: (fh_int[-1] + 1)]  # include end point
 
         if "exog" in inspect.signature(self._forecaster.__init__).parameters.keys():
             y_pred = self._fitted_forecaster.predict(start=start, end=end, exog=X)
@@ -117,7 +117,9 @@ class _StatsModelsAdapter(BaseForecaster):
 
         # statsmodels forecasts all periods from start to end of forecasting
         # horizon, but only return given time points in forecasting horizon
-        y_pred = y_pred.loc[fh_abs]
+        # if fh[0] > 1 steps ahead of cutoff then make relative to `start`
+        fh_int = fh_int - fh_int[0]
+        y_pred = y_pred.iloc[fh_int]
         # ensure that name is not added nor removed
         # otherwise this may upset conversion to pd.DataFrame
         y_pred.name = self._y.name
@@ -186,7 +188,9 @@ class _StatsModelsAdapter(BaseForecaster):
             return BaseForecaster._predict_interval(self, fh, X=X, coverage=coverage)
 
         start, end = fh.to_absolute_int(self._y.index[0], self.cutoff)[[0, -1]]
-        valid_indices = fh.to_absolute(self.cutoff).to_pandas()
+        fh_int = fh.to_absolute_int(self._y.index[0], self.cutoff) - len(self._y)
+        # if fh > 1 steps ahead of cutoff
+        fh_int = fh_int - fh_int[0]
 
         get_prediction_arguments = {"start": start, "end": end}
 
@@ -205,17 +209,15 @@ class _StatsModelsAdapter(BaseForecaster):
         var_names = self._get_varnames()
         var_name = var_names[0]
         columns = pd.MultiIndex.from_product([var_names, coverage, ["lower", "upper"]])
-        pred_int = pd.DataFrame(index=valid_indices, columns=columns)
+        preds_index = self._extract_conf_int(prediction_results, (1 - coverage[0]))
+        preds_index = preds_index.iloc[fh_int].index
+        pred_int = pd.DataFrame(index=preds_index, columns=columns)
 
         for c in coverage:
             pred_statsmodels = self._extract_conf_int(prediction_results, (1 - c))
 
-            pred_int[(var_name, c, "lower")] = pred_statsmodels.loc[
-                valid_indices, "lower"
-            ]
-            pred_int[(var_name, c, "upper")] = pred_statsmodels.loc[
-                valid_indices, "upper"
-            ]
+            pred_int[(var_name, c, "lower")] = pred_statsmodels.iloc[fh_int]["lower"]
+            pred_int[(var_name, c, "upper")] = pred_statsmodels.iloc[fh_int]["upper"]
 
         return pred_int
 
diff --git a/sktime/forecasting/compose/_pipeline.py b/sktime/forecasting/compose/_pipeline.py
index b802dfa38..cae9a4c1e 100644
--- a/sktime/forecasting/compose/_pipeline.py
+++ b/sktime/forecasting/compose/_pipeline.py
@@ -170,7 +170,7 @@ class _Pipeline(_HeterogenousMetaEstimator, BaseForecaster):
                         if len(levels) == 1:
                             levels = levels[0]
                         yt[ix] = y.xs(ix, level=levels, axis=1)
-                        # todo 0.25.0 - check why this cannot be easily removed
+                        # todo 0.26.0 - check why this cannot be easily removed
                         # in theory, we should get rid of the "Coverage" case treatment
                         # (the legacy naming convention was removed in 0.23.0)
                         # deal with the "Coverage" case, we need to get rid of this
diff --git a/sktime/forecasting/model_evaluation/_functions.py b/sktime/forecasting/model_evaluation/_functions.py
index 613360937..a0a62c01c 100644
--- a/sktime/forecasting/model_evaluation/_functions.py
+++ b/sktime/forecasting/model_evaluation/_functions.py
@@ -318,8 +318,6 @@ def _evaluate_window(x, meta):
         return result
 
 
-# todo 0.25.0: remove compute argument and docstring
-# todo 0.25.0: remove kwargs and docstring
 def evaluate(
     forecaster,
     cv,
@@ -330,7 +328,6 @@ def evaluate(
     return_data: bool = False,
     error_score: Union[str, int, float] = np.nan,
     backend: Optional[str] = None,
-    compute: bool = None,
     cv_X=None,
     backend_params: Optional[dict] = None,
     **kwargs,
@@ -348,13 +345,13 @@ def evaluate(
     the train/test folds produced by the generator ``cv_X.split_series(X)``
     (if ``X`` is ``None``, consider these to be ``None`` as well).
 
-    1. Set ``i = 1``
+    1. Initialize the counter to ``i = 1``
     2. Fit the ``forecaster`` to :math:`y_{train, 1}`, :math:`X_{train, 1}`,
-       with a ``fh`` to forecast :math:`y_{test, 1}`
-    3. The ``forecaster`` predict with exogeneous data :math:`X_{test, i}`
-       ``y_pred = forecaster.predict`` (or ``predict_proba`` or ``predict_quantiles``,
-       depending on ``scoring``)
-    4. Compute ``scoring`` on ``y_pred`` versus :math:`y_{test, 1}`
+       with ``fh`` set to the absolute indices of :math:`y_{test, 1}`.
+    3. Use the ``forecaster`` to make a prediction ``y_pred`` with the exogeneous
+      data :math:`X_{test, i}`. Predictions are made using either ``predict``,
+      ``predict_proba`` or ``predict_quantiles``, depending on ``scoring``.
+    4. Compute the ``scoring`` function on ``y_pred`` versus :math:`y_{test, i}`
     5. If ``i == K``, terminate, otherwise
     6. Set ``i = i + 1``
     7. Ingest more data :math:`y_{train, i}`, :math:`X_{train, i}`,
@@ -422,9 +419,6 @@ def evaluate(
         "threading" is unlikely to see speed ups due to the GIL and the serialization
         backend (``cloudpickle``) for "dask" and "loky" is generally more robust
         than the standard ``pickle`` library used in "multiprocessing".
-    compute : bool, default=True, deprecated and will be removed in 0.25.0.
-        If backend="dask", whether returned DataFrame is computed.
-        If set to True, returns `pd.DataFrame`, otherwise `dask.dataframe.DataFrame`.
     cv_X : sktime BaseSplitter descendant, optional
         determines split of ``X`` into test and train folds
         default is ``X`` being split to identical ``loc`` indices as ``y``
@@ -525,10 +519,10 @@ def evaluate(
                 "installed, but dask is not present in the python environment"
             )
 
-    # todo 0.25.0: remove kwargs and this warning
+    # todo 0.26.0: remove kwargs and this warning
     if kwargs != {}:
         warnings.warn(
-            "in evaluate, kwargs will no longer be supported from sktime 0.25.0. "
+            "in evaluate, kwargs are no longer supported. "
             "to pass configuration arguments to the parallelization backend, "
             "use backend_params instead. "
             f"The following kwargs were found: {kwargs.keys()}, pass these as "
@@ -537,28 +531,8 @@ def evaluate(
             stacklevel=2,
         )
 
-    # todo 0.25.0: remove compute argument and logic, and remove this warning
-    if compute is not None:
-        warnings.warn(
-            "the compute argument of evaluate is deprecated and will be removed "
-            "in sktime 0.25.0. For the same behaviour in the future, "
-            'use backend="dask_lazy"',
-            DeprecationWarning,
-            stacklevel=2,
-        )
-    if compute is None:
-        compute = True
-    if backend == "dask" and not compute:
-        backend = "dask_lazy"
-
     _check_strategy(strategy)
     cv = check_cv(cv, enforce_start_with_window=True)
-    # TODO: remove lines(four lines below) and 599-612 in v0.25.0
-    if isinstance(scoring, list):
-        raise_warn, num = True, len(scoring)
-    else:
-        raise_warn, num = False, 1
-    # removal until here
     scoring = _check_scores(scoring)
 
     ALLOWED_SCITYPES = ["Series", "Panel", "Hierarchical"]
@@ -669,21 +643,4 @@ def evaluate(
     # final formatting of results DataFrame
     results = results.reset_index(drop=True)
 
-    # TODO: remove 16 lines below and 451-455 in v0.25.0
-    if raise_warn:
-        warnings.warn(
-            "Starting v0.25.0 model_evaluation.evaluate module will rearrange "
-            "all metric columns to the left of its output result DataFrame. "
-            "Please use loc references when addressing the columns. You can "
-            "safely ignore this warning if you don't use evaluate function directly.",
-            DeprecationWarning,
-            stacklevel=2,
-        )
-        columns = results.columns.to_list()
-        non_first_metrics = []
-        for _ in range(1, num):
-            metric = columns.pop(1)
-            non_first_metrics.append(metric)
-        results = results.reindex(columns=columns + non_first_metrics)
-    #  removal until here
     return results
diff --git a/sktime/forecasting/model_selection/__init__.py b/sktime/forecasting/model_selection/__init__.py
index c55f73960..da08a4853 100644
--- a/sktime/forecasting/model_selection/__init__.py
+++ b/sktime/forecasting/model_selection/__init__.py
@@ -2,33 +2,98 @@
 # copyright: sktime developers, BSD-3-Clause License (see LICENSE file)
 """Implements functionality for selecting forecasting models."""
 
-__author__ = ["mloning", "kkoralturk"]
 __all__ = [
-    "CutoffSplitter",
-    "SameLocSplitter",
-    "SingleWindowSplitter",
-    "SlidingWindowSplitter",
-    "temporal_train_test_split",
-    "ExpandingGreedySplitter",
-    "ExpandingWindowSplitter",
-    "TestPlusTrainSplitter",
     "ForecastingGridSearchCV",
     "ForecastingRandomizedSearchCV",
     "ForecastingSkoptSearchCV",
+    "ExpandingWindowSplitter",
+    "SlidingWindowSplitter",
+    "temporal_train_test_split",
 ]
 
-from sktime.forecasting.model_selection._split import (
-    CutoffSplitter,
-    ExpandingGreedySplitter,
-    ExpandingWindowSplitter,
-    SameLocSplitter,
-    SingleWindowSplitter,
-    SlidingWindowSplitter,
-    TestPlusTrainSplitter,
-    temporal_train_test_split,
-)
 from sktime.forecasting.model_selection._tune import (
     ForecastingGridSearchCV,
     ForecastingRandomizedSearchCV,
     ForecastingSkoptSearchCV,
 )
+
+
+# todo 0.26.0 - check whether we should remove, otherwise bump
+# still used in blog posts and old tutorials
+def temporal_train_test_split(
+    y, X=None, test_size=None, train_size=None, fh=None, anchor="start"
+):
+    """Split time series data into temporal train and test sets.
+
+    DEPRECATED - use sktime.split.temporal_train_test_split instead.
+    """
+    from warnings import warn
+
+    from sktime.split import temporal_train_test_split as _tts
+
+    warn(
+        "WARNING - the old location of temporal_train_test_split in "
+        "sktime.forecasting.model_selection is deprecated and is scheduled for "
+        "imminent removal in a MINOR version. "
+        "Please update any import statements to "
+        "from sktime.split import temporal_train_test_split.",
+        DeprecationWarning,
+    )
+
+    return _tts(
+        y=y, X=X, test_size=test_size, train_size=train_size, fh=fh, anchor=anchor
+    )
+
+
+# todo 0.26.0 - check whether we should remove, otherwise bump
+# still used in blog posts and old tutorials
+def ExpandingWindowSplitter(fh=1, initial_window=10, step_length=1):
+    """Legacy export of Expanding window splitter.
+
+    DEPRECATED - use sktime.split.ExpandingWindowSplitter instead.
+    """
+    from warnings import warn
+
+    from sktime.split import ExpandingWindowSplitter as _EWSplitter
+
+    warn(
+        "WARNING - the old location of ExpandingWindowSplitter in "
+        "sktime.forecasting.model_selection is deprecated and is scheduled for "
+        "imminent removal in a MINOR version. "
+        "Please update any import statements to "
+        "from sktime.split import ExpandingWindowSplitter.",
+        DeprecationWarning,
+    )
+
+    return _EWSplitter(fh=fh, initial_window=initial_window, step_length=step_length)
+
+
+# todo 0.26.0 - check whether we should remove, otherwise bump
+# still used in blog posts and old tutorials
+def SlidingWindowSplitter(
+    fh=1, window_length=10, step_length=1, initial_window=None, start_with_window=True
+):
+    """Legacy export of Sliding window splitter.
+
+    DEPRECATED - use sktime.split.ExpandingWindowSplitter instead.
+    """
+    from warnings import warn
+
+    from sktime.split import SlidingWindowSplitter as _SWSplitter
+
+    warn(
+        "WARNING - the old location of SlidingWindowSplitter in "
+        "sktime.forecasting.model_selection is deprecated and is scheduled for "
+        "imminent removal in a MINOR version. "
+        "Please update any import statements to "
+        "from sktime.split import SlidingWindowSplitter.",
+        DeprecationWarning,
+    )
+
+    return _SWSplitter(
+        fh=fh,
+        window_length=window_length,
+        step_length=step_length,
+        initial_window=initial_window,
+        start_with_window=start_with_window,
+    )
diff --git a/sktime/forecasting/model_selection/_split.py b/sktime/forecasting/model_selection/_split.py
deleted file mode 100644
index 3d77d9d77..000000000
--- a/sktime/forecasting/model_selection/_split.py
+++ /dev/null
@@ -1,42 +0,0 @@
-#!/usr/bin/env python3 -u
-# copyright: sktime developers, BSD-3-Clause License (see LICENSE file)
-"""Implement dataset splitting for model evaluation and selection."""
-
-# This module is deprecated. Below the deprecation cycle :
-# - Until the 0.25.0 release, imports are allowed from this location but raise a
-# DeprecationWarning.
-# - In the 0.25.0 release, this file (sktime/forecasting/model_selection/_split.py)
-# will be removed.
-# todo 0.25.0 : Please ensure this file is deleted.
-
-__all__ = [
-    "ExpandingGreedySplitter",
-    "ExpandingWindowSplitter",
-    "SlidingWindowSplitter",
-    "CutoffSplitter",
-    "SingleWindowSplitter",
-    "SameLocSplitter",
-    "temporal_train_test_split",
-    "TestPlusTrainSplitter",
-]
-
-import warnings
-
-from sktime.split import (
-    CutoffSplitter,
-    ExpandingGreedySplitter,
-    ExpandingWindowSplitter,
-    SameLocSplitter,
-    SingleWindowSplitter,
-    SlidingWindowSplitter,
-    TestPlusTrainSplitter,
-    temporal_train_test_split,
-)
-
-warnings.warn(
-    "Please import sktime splitters from `sktime.split`. Importing splitters from "
-    "`sktime.forecasting.model_selection` is deprecated and will be removed in release "
-    "0.25.0",
-    DeprecationWarning,
-    stacklevel=2,
-)
diff --git a/sktime/forecasting/model_selection/_tune.py b/sktime/forecasting/model_selection/_tune.py
index fd0fcd585..5c11825bf 100644
--- a/sktime/forecasting/model_selection/_tune.py
+++ b/sktime/forecasting/model_selection/_tune.py
@@ -544,9 +544,7 @@ class ForecastingGridSearchCV(BaseGridSearch):
     --------
     >>> from sktime.datasets import load_shampoo_sales
     >>> from sktime.forecasting.model_selection import ForecastingGridSearchCV
-    >>> from sktime.split import (
-    ...     ExpandingWindowSplitter,
-    ...     ExpandingWindowSplitter)
+    >>> from sktime.split import ExpandingWindowSplitter
     >>> from sktime.forecasting.naive import NaiveForecaster
     >>> y = load_shampoo_sales()
     >>> fh = [1,2,3]
diff --git a/sktime/forecasting/trend/_pwl_trend_forecaster.py b/sktime/forecasting/trend/_pwl_trend_forecaster.py
index c3d947660..e59bc3fb9 100644
--- a/sktime/forecasting/trend/_pwl_trend_forecaster.py
+++ b/sktime/forecasting/trend/_pwl_trend_forecaster.py
@@ -60,8 +60,8 @@ class ProphetPiecewiseLinearTrendForecaster(_ProphetAdapter):
     --------
     >>> from sktime.datasets import load_airline
     >>> from sktime.forecasting.trend import ProphetPiecewiseLinearTrendForecaster
-    >>> from sktime.forecasting.model_selection import temporal_train_test_split
     >>> from sktime.forecasting.base import ForecastingHorizon
+    >>> from sktime.split import temporal_train_test_split
     >>> y =load_airline().to_timestamp(freq='M')
     >>> y_train, y_test = temporal_train_test_split(y)
     >>> fh = ForecastingHorizon(y.index, is_relative=False)
diff --git a/sktime/registry/_tags.py b/sktime/registry/_tags.py
index 1290a351d..b40bb47b3 100644
--- a/sktime/registry/_tags.py
+++ b/sktime/registry/_tags.py
@@ -492,6 +492,12 @@ ESTIMATOR_TAG_REGISTER = [
         "int",
         "sample size used in approximating other statistics if not available",
     ),
+    (
+        "capability:multioutput",
+        ["classifier", "regressor"],  # might need to add "early_classifier" here
+        "bool",
+        "can the estimator handle multioutput data?",
+    ),
 ]
 
 ESTIMATOR_TAG_TABLE = pd.DataFrame(ESTIMATOR_TAG_REGISTER)
diff --git a/sktime/regression/base.py b/sktime/regression/base.py
index 903cc0cf7..f89e1aa86 100644
--- a/sktime/regression/base.py
+++ b/sktime/regression/base.py
@@ -23,24 +23,16 @@ __all__ = [
 __author__ = ["mloning", "fkiraly"]
 
 import time
-from abc import ABC, abstractmethod
 
 import numpy as np
-import pandas as pd
-
-from sktime.base import BaseEstimator
-from sktime.datatypes import (
-    MTYPE_LIST_PANEL,
-    check_is_error_msg,
-    check_is_scitype,
-    convert_to,
-)
+
+from sktime.base import BasePanelMixin
+from sktime.datatypes import VectorizedDF
 from sktime.utils.sklearn import is_sklearn_transformer
 from sktime.utils.validation import check_n_jobs
-from sktime.utils.warnings import warn
 
 
-class BaseRegressor(BaseEstimator, ABC):
+class BaseRegressor(BasePanelMixin):
     """Abstract base class for time series regressors.
 
     The base regressor specifies the methods and method signatures that all
@@ -57,7 +49,9 @@ class BaseRegressor(BaseEstimator, ABC):
     _tags = {
         "object_type": "regressor",  # type of object
         "X_inner_mtype": "numpy3D",  # which type do _fit/_predict, support for X?
+        "y_inner_mtype": "numpy1D",  # which type do _fit/_predict, support for y?
         #    it should be either "numpy3D" or "nested_univ" (nested pd.DataFrame)
+        "capability:multioutput": False,  # whether regressor supports multioutput
         "capability:multivariate": False,
         "capability:unequal_length": False,
         "capability:missing_values": False,
@@ -75,6 +69,14 @@ class BaseRegressor(BaseEstimator, ABC):
         "is_equal_length",
     ]
 
+    # attribute name where vectorized estimators are stored
+    VECTORIZATION_ATTR = "regressors_"  # e.g., classifiers_, regressors_
+
+    # used in error messages
+    TASK = "regression"  # e.g., classification, regression
+    EST_TYPE = "regressor"  # e.g., classifier, regressor
+    EST_TYPE_PLURAL = "regressors"  # e.g., classifiers, regressors
+
     def __init__(self):
         self.fit_time_ = 0
         self._class_dictionary = {}
@@ -84,6 +86,9 @@ class BaseRegressor(BaseEstimator, ABC):
         # required for compatibility with some sklearn interfaces
         # i.e. CalibratedRegressorCV
         self._estimator_type = "regressor"
+        self._is_vectorized = False
+        self._is_timed = False
+        self._converter_store_y = {}
 
         super().__init__()
 
@@ -128,19 +133,31 @@ class BaseRegressor(BaseEstimator, ABC):
     def fit(self, X, y):
         """Fit time series regressor to training data.
 
+        State change:
+            Changes state to "fitted".
+
+        Writes to self:
+            Sets self.is_fitted to True.
+            Sets fitted model attributes ending in "_".
+
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of float, of shape [n_instances] - regression labels for fitting
-            indices correspond to instance indices in X
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
+        y : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            supported sktime types: np.ndarray (1D, 2D), pd.Series, pd.DataFrame
 
         Returns
         -------
@@ -151,21 +168,46 @@ class BaseRegressor(BaseEstimator, ABC):
         Changes state by creating a fitted model that updates attributes
         ending in "_" and sets is_fitted flag to True.
         """
+        self.reset()
+
+        # fit timer start
         start = int(round(time.time() * 1000))
+
+        # check and convert y for multioutput vectorization
+        y, y_metadata = self._check_y(y)
+        self._y_metadata = y_metadata
+        self._is_vectorized = isinstance(y, VectorizedDF)
+
+        if self._is_vectorized:
+            self._vectorize("fit", X=X, y=y)
+            # fit timer end
+            self.fit_time_ = int(round(time.time() * 1000)) - start
+            # this should happen last: fitted state is set to True
+            self._is_fitted = True
+            return self
+
+        # no vectorization needed, proceed with normal fit
+
         # convenience conversions to allow user flexibility:
         # if X is 2D array, convert to 3D, if y is Series, convert to numpy
-        X, y = _internal_convert(X, y)
-        X_metadata = _check_regressor_input(
+        X, y = self._internal_convert(X, y)
+
+        # y float coercion
+        if y is not None and isinstance(y, np.ndarray):
+            y = y.astype("float")
+
+        # input checks
+        X_metadata = self._check_input(
             X, y, return_metadata=self.METADATA_REQ_IN_CHECKS
         )
         self._X_metadata = X_metadata
-        missing = X_metadata["has_nans"]
-        multivariate = not X_metadata["is_univariate"]
-        unequal = not X_metadata["is_equal_length"]
+        X_mtype = X_metadata["mtype"]
+
         # Check this regressor can handle characteristics
-        self._check_capabilities(missing, multivariate, unequal)
+        self._check_capabilities(X_metadata)
+
         # Convert data as dictated by the regressor tags
-        X = self._convert_X(X)
+        X = self._convert_X(X, X_mtype)
         multithread = self.get_tag("capability:multithreading")
         if multithread:
             try:
@@ -177,7 +219,8 @@ class BaseRegressor(BaseEstimator, ABC):
 
         self._fit(X, y)
         self.fit_time_ = int(round(time.time() * 1000)) - start
-        # this should happen last
+
+        # this should happen last: fitted state is set to True
         self._is_fitted = True
         return self
 
@@ -186,44 +229,64 @@ class BaseRegressor(BaseEstimator, ABC):
 
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
 
         Returns
         -------
-        y : 1D np.array of float, of shape [n_instances] - predicted regression labels
-            indices correspond to instance indices in X
+        y_pred : sktime compatible tabular data container, Table scitype
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            predicted class labels
+            0-th indices correspond to instance indices in X
+            1-st indices (if applicable) correspond to multioutput vector indices in X
+            1D np.npdarray, if y univariate (one dimension)
+            otherwise, same type as y passed in fit
         """
         self.check_is_fitted()
 
+        # enter vectorized mode if needed
+        if self._is_vectorized:
+            return self._vectorize("predict", X=X)
+
         # boilerplate input checks for predict-like methods
         X = self._check_convert_X_for_predict(X)
 
-        return self._predict(X)
+        # call internal _predict, convert output
+        y_pred_inner = self._predict(X)
+        y_pred = self._convert_output_y(y_pred_inner)
+        return y_pred
 
-    def score(self, X, y) -> float:
+    def score(self, X, y, multioutput="uniform_average") -> float:
         """Scores predicted labels against ground truth labels on X.
 
         Parameters
         ----------
-        X : 3D np.array (any number of dimensions, equal length series)
-                of shape [n_instances, n_dimensions, series_length]
-            or 2D np.array (univariate, equal length series)
-                of shape [n_instances, series_length]
-            or pd.DataFrame with each column a dimension, each cell a pd.Series
-                (any number of dimensions, equal or unequal length series)
+        X : sktime compatible time series panel data container, Panel scitype, e.g.,
+            pd-multiindex: pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
+            numpy3D: 3D np.array (any number of dimensions, equal length series)
+            of shape [n_instances, n_dimensions, series_length]
             or of any other supported Panel mtype
-                for list of mtypes, see datatypes.SCITYPE_REGISTER
-                for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of float, of shape [n_instances] - regression labels (gnr truth)
-            indices correspond to instance indices in X
+            for list of mtypes, see datatypes.SCITYPE_REGISTER
+            for specifications, see examples/AA_datatypes_and_datasets.ipynb
+        y : 2D np.array of int, of shape [n_instances, n_dimensions] - regression labels
+            for fitting indices correspond to instance indices in X
+            or 1D np.array of int, of shape [n_instances] - regression labels for
+            fitting indices correspond to instance indices in X
+        multioutput : str, optional (default="uniform_average")
+            {"raw_values", "uniform_average", "variance_weighted"}, array-like of shape
+            (n_outputs,) or None, default="uniform_average".
+            Defines aggregating of multiple output scores. Array-like value defines
+            weights used to average scores.
 
         Returns
         -------
@@ -233,9 +296,8 @@ class BaseRegressor(BaseEstimator, ABC):
 
         self.check_is_fitted()
 
-        return r2_score(y, self.predict(X), normalize=True)
+        return r2_score(y, self.predict(X), normalize=True, multioutput=multioutput)
 
-    @abstractmethod
     def _fit(self, X, y):
         """Fit time series regressor to training data.
 
@@ -245,13 +307,21 @@ class BaseRegressor(BaseEstimator, ABC):
         ----------
         X : guaranteed to be of a type in self.get_tag("X_inner_mtype")
             if self.get_tag("X_inner_mtype") = "numpy3D":
-                3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            3D np.ndarray of shape = [n_instances, n_dimensions, series_length]
+            if self.get_tag("X_inner_mtype") = "pd-multiindex":
+            pd.DataFrame with columns = variables,
+            index = pd.MultiIndex with first level = instance indices,
+            second level = time indices
             if self.get_tag("X_inner_mtype") = "nested_univ":
-                pd.DataFrame with each column a dimension, each cell a pd.Series
+            pd.DataFrame with each column a dimension, each cell a pd.Series
             for list of other mtypes, see datatypes.SCITYPE_REGISTER
             for specifications, see examples/AA_datatypes_and_datasets.ipynb
-        y : 1D np.array of float, of shape [n_instances] - regression labels for fitting
-            indices correspond to instance indices in X
+        y : guaranteed to be of a type in self.get_tag("y_inner_mtype")
+            1D iterable, of shape [n_instances]
+            or 2D iterable, of shape [n_instances, n_dimensions]
+            class labels for fitting
+            if self.get_tag("capaility:multioutput") = False, guaranteed to be 1D
+            if self.get_tag("capaility:multioutput") = True, guaranteed to be 2D
 
         Returns
         -------
@@ -261,9 +331,8 @@ class BaseRegressor(BaseEstimator, ABC):
         -----
         Changes state by creating a fitted model that updates attributes ending in "_"
         """
-        ...
+        raise NotImplementedError("abstract method")
 
-    @abstractmethod
     def _predict(self, X) -> np.ndarray:
         """Predicts labels for sequences in X.
 
@@ -284,213 +353,4 @@ class BaseRegressor(BaseEstimator, ABC):
         y : 1D np.array of float, of shape [n_instances] - predicted regression labels
             indices correspond to instance indices in X
         """
-        ...
-
-    def _check_convert_X_for_predict(self, X):
-        """Input checks, capability checks, repeated in all predict/score methods.
-
-        Parameters
-        ----------
-        X : any object (to check/convert)
-            should be of a supported Panel mtype or 2D numpy.ndarray
-
-        Returns
-        -------
-        X: an object of a supported Panel mtype, numpy3D if X was a 2D numpy.ndarray
-
-        Raises
-        ------
-        ValueError if X is of invalid input data type, or there is not enough data
-        ValueError if the capabilities in self._tags do not handle the data.
-        """
-        X = _internal_convert(X)
-        X_metadata = _check_regressor_input(
-            X, return_metadata=self.METADATA_REQ_IN_CHECKS
-        )
-        missing = X_metadata["has_nans"]
-        multivariate = not X_metadata["is_univariate"]
-        unequal = not X_metadata["is_equal_length"]
-        # Check this regressor can handle characteristics
-        self._check_capabilities(missing, multivariate, unequal)
-        # Convert data as dictated by the regressor tags
-        X = self._convert_X(X)
-
-        return X
-
-    def _check_capabilities(self, missing, multivariate, unequal):
-        """Check whether this regressor can handle the data characteristics.
-
-        Parameters
-        ----------
-        missing : boolean, does the data passed to fit contain missing values?
-        multivariate : boolean, does the data passed to fit contain missing values?
-        unequal : boolea, do the time series passed to fit have variable lengths?
-
-        Raises
-        ------
-        ValueError if the capabilities in self._tags do not handle the data.
-        """
-        allow_multivariate = self.get_tag("capability:multivariate")
-        allow_missing = self.get_tag("capability:missing_values")
-        allow_unequal = self.get_tag("capability:unequal_length")
-
-        self_name = type(self).__name__
-
-        # identify problems, mismatch of capability and inputs
-        problems = []
-        if missing and not allow_missing:
-            problems += ["missing values"]
-        if multivariate and not allow_multivariate:
-            problems += ["multivariate series"]
-        if unequal and not allow_unequal:
-            problems += ["unequal length series"]
-
-        # construct error message
-        problems_and = " and ".join(problems)
-        problems_or = " or ".join(problems)
-        msg = (
-            f"Data seen by {self_name} instance has {problems_and}, "
-            f"but this {self_name} instance cannot handle {problems_or}. "
-            f"Calls with {problems_or} may result in error or unreliable results."
-        )
-
-        # raise exception or warning with message
-        # if self is composite, raise a warning, since passing could be fine
-        #   see discussion in PR 2366 why
-        if len(problems) > 0:
-            if self.is_composite():
-                warn(msg, obj=self)
-            else:
-                raise ValueError(msg)
-
-    def _convert_X(self, X):
-        """Convert equal length series from DataFrame to numpy array or vice versa.
-
-        Parameters
-        ----------
-        self : this regressor
-        X : pd.DataFrame or np.ndarray. Input attribute data
-
-        Returns
-        -------
-        X : input X converted to type in "X_inner_mtype" tag
-                usually a pd.DataFrame (nested) or 3D np.ndarray
-            Checked and possibly converted input data
-        """
-        inner_type = self.get_tag("X_inner_mtype")
-        # convert pd.DataFrame
-        X = convert_to(
-            X,
-            to_type=inner_type,
-            as_scitype="Panel",
-        )
-        return X
-
-
-def _check_regressor_input(
-    X,
-    y=None,
-    enforce_min_instances=1,
-    return_metadata=True,
-):
-    """Check whether input X and y are valid formats with minimum data.
-
-    Raises a ValueError if the input is not valid.
-
-    Parameters
-    ----------
-    X : check whether conformant with any sktime Panel mtype specification
-    y : check whether a pd.Series or np.array
-    enforce_min_instances : int, optional (default=1)
-        check there are a minimum number of instances.
-    return_metadata : bool, str, or list of str
-        metadata fields to return with X_metadata, input to check_is_scitype
-
-    Returns
-    -------
-    metadata : dict with metadata for X returned by datatypes.check_is_scitype
-
-    Raises
-    ------
-    ValueError
-        If y or X is invalid input data type, or there is not enough data
-    """
-    # Check X is valid input type and recover the data characteristics
-    X_valid, msg, X_metadata = check_is_scitype(
-        X, scitype="Panel", return_metadata=return_metadata
-    )
-    # raise informative error message if X is in wrong format
-    allowed_msg = (
-        f"Allowed scitypes for classifiers are Panel mtypes, "
-        f"for instance a pandas.DataFrame with MultiIndex and last(-1) "
-        f"level an sktime compatible time index. "
-        f"Allowed compatible mtype format specifications are: {MTYPE_LIST_PANEL} ."
-    )
-    if not X_valid:
-        check_is_error_msg(
-            msg, var_name="X", allowed_msg=allowed_msg, raise_exception=True
-        )
-
-    n_cases = X_metadata["n_instances"]
-    if n_cases < enforce_min_instances:
-        raise ValueError(
-            f"Minimum number of cases required is {enforce_min_instances} but X "
-            f"has : {n_cases}"
-        )
-
-    # Check y if passed
-    if y is not None:
-        # Check y valid input
-        if not isinstance(y, (pd.Series, np.ndarray)) and not (
-            isinstance(y, pd.DataFrame) and y.shape[1] == 1
-        ):
-            raise ValueError(
-                f"y must be a np.array, pd.Series or a 1-D pd.DataFrame,"
-                f"but found type: {type(y)}"
-            )
-        # Check matching number of labels
-        n_labels = y.shape[0]
-        if n_cases != n_labels:
-            raise ValueError(
-                f"Mismatch in number of cases. Number in X = {n_cases} nos in y = "
-                f"{n_labels}"
-            )
-        if isinstance(y, np.ndarray):
-            if y.ndim > 1:
-                raise ValueError(
-                    f"np.ndarray y must be 1-dimensional, "
-                    f"but found {y.ndim} dimensions"
-                )
-    return X_metadata
-
-
-def _internal_convert(X, y=None):
-    """Convert X and y if necessary as a user convenience.
-
-    Convert X to a 3D numpy array if already a 2D and convert y into an 1D numpy
-    array if passed as a Series.
-
-    Parameters
-    ----------
-    X : an object of a supported Panel mtype, or 2D numpy.ndarray
-    y : np.ndarray or pd.Series
-
-    Returns
-    -------
-    X: an object of a supported Panel mtype, numpy3D if X was a 2D numpy.ndarray
-    y: np.ndarray
-    """
-    if isinstance(X, np.ndarray):
-        # Temporary fix to insist on 3D numpy. For univariate problems,
-        # most regressors simply convert back to 2D. This squeezing should be
-        # done here, but touches a lot of files, so will get this to work first.
-        if X.ndim == 2:
-            X = X.reshape(X.shape[0], 1, X.shape[1])
-    if y is not None and isinstance(y, pd.Series):
-        # y should be a numpy array, although we allow Series for user convenience
-        y = pd.Series.to_numpy(y)
-    if y is not None and isinstance(y, np.ndarray):
-        y = y.astype("float")
-    if y is None:
-        return X
-    return X, y
+        raise NotImplementedError("abstract method")
diff --git a/sktime/regression/deep_learning/rnn.py b/sktime/regression/deep_learning/rnn.py
index c32bd8199..068fe1084 100644
--- a/sktime/regression/deep_learning/rnn.py
+++ b/sktime/regression/deep_learning/rnn.py
@@ -11,8 +11,13 @@ from sklearn.utils import check_random_state
 from sktime.networks.rnn import RNNNetwork
 from sktime.regression.deep_learning.base import BaseDeepRegressor
 from sktime.utils.validation._dependencies import _check_dl_dependencies
+from sktime.utils.warnings import warn
 
 
+# todo: 0.26.0 - please remove num_epochs parameter and related logic,
+# because parameter num_epochs usage is deprecated and will be renamed
+# to n_epochs.  Also n_epochs should be moved to become the first
+# argument in the __init__ and super().__init__.
 class SimpleRNNRegressor(BaseDeepRegressor):
     """Simple recurrent neural network.
 
@@ -51,7 +56,7 @@ class SimpleRNNRegressor(BaseDeepRegressor):
 
     def __init__(
         self,
-        num_epochs=100,
+        num_epochs=None,
         batch_size=1,
         units=6,
         callbacks=None,
@@ -63,10 +68,26 @@ class SimpleRNNRegressor(BaseDeepRegressor):
         activation="linear",
         use_bias=True,
         optimizer=None,
+        n_epochs=100,
     ):
         _check_dl_dependencies(severity="error")
         super().__init__()
+        # todo: 0.26.0 - remove this, replace by
+        # self._n_epochs = n_epochs
+        # Deprecated Parameter Handling
+        if num_epochs is not None:
+            warn(
+                "In SimpleRNNRegressor, the parameter 'num_epochs' is deprecated "
+                "and will be removed in v0.26.0. It will be renamed to n_epochs. "
+                "To avoid this warning, update your code to use 'n_epochs'. ",
+                DeprecationWarning,
+                stacklevel=2,
+            )
+            self._n_epochs = num_epochs
+        else:
+            self._n_epochs = n_epochs
         self.num_epochs = num_epochs
+        # end remove
         self.batch_size = batch_size
         self.verbose = verbose
         self.units = units
@@ -80,6 +101,7 @@ class SimpleRNNRegressor(BaseDeepRegressor):
         self.optimizer = optimizer
         self.history = None
         self._network = RNNNetwork(random_state=random_state, units=units)
+        self.n_epochs = n_epochs
 
     def build_model(self, input_shape, **kwargs):
         """Construct a compiled, un-trained, keras model that is ready for training.
@@ -183,7 +205,7 @@ class SimpleRNNRegressor(BaseDeepRegressor):
             X,
             y,
             batch_size=self.batch_size,
-            epochs=self.num_epochs,
+            epochs=self._n_epochs,
             verbose=self.verbose,
             callbacks=self.callbacks_,
         )
@@ -212,7 +234,7 @@ class SimpleRNNRegressor(BaseDeepRegressor):
         """
         params1 = {}
         params2 = {
-            "num_epochs": 50,
+            "n_epochs": 50,
             "batch_size": 2,
             "units": 5,
             "use_bias": False,
diff --git a/sktime/split/sameloc.py b/sktime/split/sameloc.py
index 040c748d3..c94e1bb6c 100644
--- a/sktime/split/sameloc.py
+++ b/sktime/split/sameloc.py
@@ -141,10 +141,7 @@ class SameLocSplitter(BaseSplitter):
             `create_test_instance` uses the first (or only) dictionary in `params`
         """
         from sktime.datasets import load_airline
-        from sktime.forecasting.model_selection import (
-            ExpandingWindowSplitter,
-            SingleWindowSplitter,
-        )
+        from sktime.split import ExpandingWindowSplitter, SingleWindowSplitter
 
         y = load_airline()
         y_temp = y[:60]
diff --git a/sktime/split/slidingwindow.py b/sktime/split/slidingwindow.py
index 031ffe6e2..ed55667e4 100644
--- a/sktime/split/slidingwindow.py
+++ b/sktime/split/slidingwindow.py
@@ -51,9 +51,9 @@ class SlidingWindowSplitter(BaseWindowSplitter):
 
     Parameters
     ----------
-    fh : int, list or np.array
+    fh : int, list or np.array, optional (default=1)
         Forecasting horizon
-    window_length : int or timedelta or pd.DateOffset
+    window_length : int or timedelta or pd.DateOffset, optional (default=10)
         Window length
     step_length : int or timedelta or pd.DateOffset, optional (default=1)
         Step length between windows
diff --git a/sktime/transformations/base.py b/sktime/transformations/base.py
index 0dff77c85..23005244c 100644
--- a/sktime/transformations/base.py
+++ b/sktime/transformations/base.py
@@ -138,6 +138,7 @@ class BaseTransformer(BaseEstimator):
     }
 
     # default config values
+    # see set_config documentation for details
     _config = {
         "input_conversion": "on",
         # controls input checks and conversions,
@@ -156,41 +157,12 @@ class BaseTransformer(BaseEstimator):
         #  {None, "dask", "loky", "multiprocessing", "threading"}
         #  None: no parallelization
         #  "loky", "multiprocessing" and "threading": uses `joblib` Parallel loops
+        #  "joblib": uses custom joblib backend, set via `joblib_backend` tag
         #  "dask": uses `dask`, requires `dask` package in environment
         "backend:parallel:params": None,  # params for parallelization backend
     }
 
     _config_doc = {
-        "backend:parallel": """
-        backend:parallel : str, optional, default="None"
-            backend to use for parallelization when broadcasting/vectorizing, one of
-
-            - "None": executes loop sequentally, simple list comprehension
-            - "loky", "multiprocessing" and "threading": uses ``joblib.Parallel``
-            - "joblib": custom and 3rd party ``joblib`` backends, e.g., ``spark``
-            - "dask": uses ``dask``, requires ``dask`` package in environment
-        """,
-        "backend:parallel:params": """
-        backend:parallel:params : dict, optional, default={} (no parameters passed)
-            additional parameters passed to the parallelization backend as config.
-            Valid keys depend on the value of ``backend:parallel``:
-
-            - "None": no additional parameters, ``backend_params`` is ignored
-            - "loky", "multiprocessing" and "threading": default ``joblib`` backends
-              any valid keys for ``joblib.Parallel`` can be passed here, e.g.,
-              ``n_jobs``, with the exception of ``backend`` which is directly
-              controlled by ``backend``.
-              If ``n_jobs`` is not passed, it will default to ``-1``, other parameters
-              will default to ``joblib`` defaults.
-            - "joblib": custom and 3rd party ``joblib`` backends,
-              e.g., ``spark``. Any valid keys for ``joblib.Parallel``
-              can be passed here, e.g., ``n_jobs``,
-            ``backend`` must be passed as a key of ``backend_params`` in this case.
-              If ``n_jobs`` is not passed, it will default to ``-1``, other parameters
-              will default to ``joblib`` defaults.
-            - "dask": any valid keys for ``dask.compute``
-              can be passed, e.g., ``scheduler``
-        """,
         "input_conversion": """
         input_conversion : str, one of "on", "off", valid mtype string
             controls input checks and conversions,
@@ -1186,10 +1158,20 @@ class BaseTransformer(BaseEstimator):
         #   skipped for output_scitype = "Primitives"
         #       since then the output always is a pd.DataFrame
         if case == "case 2: higher scitype supported" and output_scitype == "Series":
-            Xt = convert_to(
-                Xt,
-                to_type=["pd-multiindex", "numpy3D", "df-list", "pd_multiindex_hier"],
-            )
+            if self.get_tags()["scitype:transform-input"] == "Panel":
+                # Conversion from Series to Panel done for being compatible with
+                # algorithm. Thus, the returned Series should stay a Series.
+                pass
+            else:
+                Xt = convert_to(
+                    Xt,
+                    to_type=[
+                        "pd-multiindex",
+                        "numpy3D",
+                        "df-list",
+                        "pd_multiindex_hier",
+                    ],
+                )
             Xt = convert_to_scitype(Xt, to_scitype=X_input_scitype)
 
         # now, in all cases, Xt is in the right scitype,
@@ -1203,7 +1185,6 @@ class BaseTransformer(BaseEstimator):
         if output_scitype == "Series":
             # output mtype is input mtype
             X_output_mtype = X_input_mtype
-
             # exception to this: if the transformer outputs multivariate series,
             #   we cannot convert back to pd.Series, do pd.DataFrame instead then
             #   this happens only for Series, not Panel
@@ -1229,6 +1210,13 @@ class BaseTransformer(BaseEstimator):
                     )
                 if X_input_mtype == "pd.Series" and not metadata["is_univariate"]:
                     X_output_mtype = "pd.DataFrame"
+            elif self.get_tags()["scitype:transform-input"] == "Panel":
+                # Input has always to be Panel
+                X_output_mtype = "pd.DataFrame"
+            else:
+                # Input can be Panel or Hierachical, since it is supported
+                # by the used mtype
+                output_scitype = X_input_scitype
                 # Xt_mtype = metadata["mtype"]
             # else:
             #     Xt_mtype = X_input_mtype
@@ -1241,10 +1229,10 @@ class BaseTransformer(BaseEstimator):
             #     store=_converter_store_X,
             #     store_behaviour="freeze",
             # )
-            Xt = convert_to(
+            return convert_to(
                 Xt,
                 to_type=X_output_mtype,
-                as_scitype=X_input_scitype,
+                as_scitype=output_scitype,
                 store=_converter_store_X,
                 store_behaviour="freeze",
             )
@@ -1260,7 +1248,7 @@ class BaseTransformer(BaseEstimator):
                 # else this is only zeros and should be reset to RangeIndex
                 else:
                     Xt = Xt.reset_index(drop=True)
-            Xt = convert_to(
+            return convert_to(
                 Xt,
                 to_type="pd_DataFrame_Table",
                 as_scitype="Table",
diff --git a/sktime/transformations/merger.py b/sktime/transformations/merger.py
new file mode 100644
index 000000000..10c7024fb
--- /dev/null
+++ b/sktime/transformations/merger.py
@@ -0,0 +1,134 @@
+# copyright: sktime developers, BSD-3-Clause License (see LICENSE file)
+"""Implements a merger for panel data."""
+import numpy as np
+
+__author__ = ["benHeid"]
+
+
+from sktime.transformations.base import BaseTransformer
+
+
+class Merger(BaseTransformer):
+    """Aggregates Panel data containing overlapping windows of one time series.
+
+    The input data contains multiple overlapping time series elements that could
+    arranged as follows:
+    xxxx.....
+    .xxxx....
+    ..xxxx...
+    ...xxxx..
+    ....xxxx.
+    .....xxxx
+    ......xxxx
+    .......xxxx
+    ........xxxx
+    .........xxxx
+    The merger aggregates the data by aligning the time series windows as shown above
+    and applying a aggregation function to the overlapping data points.
+    The aggregation function can be one of "mean" or "median". I.e., the `mean` or
+    `median` of each column is calculated, resulting in a univariate time series.
+
+    Parameters
+    ----------
+    method : {`median`, `mean`}, default="median"
+        The method to use for aggregation. Can be one of "mean" or "median".
+    stride : int, default=0
+        The stride to use for the aggregation. The stride determines the number of
+        shifts between consecutive instances. A stride of 0 means no shift. A
+        stride of 1 means that the time series is aggregated as above.
+
+    Examples
+    --------
+    >>> from sktime.transformations.merger import Merger
+    >>> from sktime.utils._testing.panel import _make_panel
+    >>> y = _make_panel(n_instances=10, n_columns=3, n_timepoints=5)
+    >>> result = Merger(method="median").fit_transform(y)
+    >>> result.shape
+    (3, 5)
+
+    >>> from sktime.transformations.merger import Merger
+    >>> from sktime.utils._testing.panel import _make_panel
+    >>> y = _make_panel(n_instances=10, n_columns=3, n_timepoints=5)
+    >>> result = Merger(method="median", stride=1).fit_transform(y)
+    >>> result.shape
+    (14, 3)
+    """
+
+    _tags = {
+        "scitype:transform-input": "Panel",
+        "scitype:transform-output": "Series",
+        "X_inner_mtype": "numpy3D",
+        "fit_is_empty": True,
+    }
+
+    def __init__(self, method="median", stride=0):
+        self.stride = stride
+        if method not in ["median", "mean"]:
+            raise ValueError(f"{method} must be 'mean' or 'median'.")
+        self.method = method
+        super().__init__()
+
+    def _transform(self, X=None, y=None):
+        """Merge the Panel data by aligning them temporally.
+
+        Parameters
+        ----------
+        X : pd.DataFrame
+            The input panel data.
+        y : pd.Series
+            ignored
+
+        Returns
+        -------
+        returns a single time series
+        """
+        horizon = X.shape[-1]
+
+        if self.method == "mean":
+            result = np.nanmean(self._align_temporal(horizon, X), axis=0)
+        elif self.method == "median":
+            result = np.nanmedian(self._align_temporal(horizon, X), axis=0)
+        else:
+            raise ValueError(f"{self.method} must be 'mean' or 'median'.")
+        return result
+
+    def _align_temporal(self, horizon, x):
+        x = x.astype(float)
+        if self.stride > 0:
+            x = np.insert(
+                x, np.arange(1, x.shape[0]).repeat(self.stride - 1), np.nan, axis=0
+            )
+        elif self.stride == 0:
+            return x
+        r = []
+        for i in range(horizon):
+            _res = np.concatenate(
+                [
+                    np.full(fill_value=np.nan, shape=(i, x.shape[1])),
+                    x[:, :, i],
+                    np.full(fill_value=np.nan, shape=((horizon - 1 - i), x.shape[1])),
+                ]
+            )
+            r.append(_res)
+        return np.stack(r)
+
+    @classmethod
+    def get_test_params(cls, parameter_set="default"):
+        """Return testing parameter settings for the estimator.
+
+        Parameters
+        ----------
+        parameter_set : str, default="default"
+            Name of the set of test parameters to return, for use in tests. If no
+            special parameters are defined for a value, will return `"default"` set.
+            There are currently no reserved values for forecasters.
+
+        Returns
+        -------
+        params : dict or list of dict, default = {}
+            Parameters to create testing instances of the class
+            Each dict are parameters to construct an "interesting" test instance, i.e.,
+            `MyClass(**params)` or `MyClass(**params[i])` creates a valid test instance.
+            `create_test_instance` uses the first (or only) dictionary in `params`
+        """
+        return [{"method": "mean"}, {"method": "median"}]


Test patch:
diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml
index 282bb2233..3a285e30b 100644
--- a/.github/workflows/test.yml
+++ b/.github/workflows/test.yml
@@ -44,6 +44,23 @@ jobs:
         run: build_tools/run_examples.sh
         shell: bash
 
+  run-blogpost-examples:
+    needs: code-quality
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v4
+      - name: Set up Python
+        uses: actions/setup-python@v5
+        with:
+          python-version: 3.9
+      - name: Install dependencies
+        run: |
+          python -m pip install --upgrade pip
+          python -m pip install .[all_extras,binder,dev,mlflow]
+      - name: Run example notebooks
+        run: build_tools/run_blogposts.sh
+        shell: bash
+
   test-nodevdeps:
     needs: code-quality
     runs-on: ubuntu-latest
diff --git a/sktime/classification/tests/test_all_classifiers.py b/sktime/classification/tests/test_all_classifiers.py
index d84533eac..cda9a0245 100644
--- a/sktime/classification/tests/test_all_classifiers.py
+++ b/sktime/classification/tests/test_all_classifiers.py
@@ -4,6 +4,7 @@ __author__ = ["mloning", "TonyBagnall", "fkiraly"]
 
 
 import numpy as np
+import pandas as pd
 import pytest
 
 from sktime.classification.tests._expected_outputs import (
@@ -207,7 +208,29 @@ class TestAllClassifiers(ClassifierFixtureGenerator, QuickTester):
         X, y = make_classification_problem()
         y[:] = 42
 
-        error_msg = "single class label"
+        error_msg = "single label"
 
         with pytest.warns(UserWarning, match=error_msg):
             estimator_instance.fit(X, y)
+
+    def test_multioutput(self, estimator_instance):
+        """Test multioutput classification for all classifiers.
+
+        All classifiers should follow the same interface,
+        those that do not genuinely should vectorize/broadcast over y.
+        """
+        n_instances = 20
+        X, y = make_classification_problem(n_instances=n_instances)
+        y_mult = pd.DataFrame({"a": y, "b": y})
+
+        estimator_instance.fit(X, y_mult)
+        y_pred = estimator_instance.predict(X)
+
+        assert isinstance(y_pred, pd.DataFrame)
+        assert y_pred.shape == y_mult.shape
+
+        vectorized = estimator_instance.get_tag("capability:multioutput")
+        if vectorized:
+            assert hasattr(estimator_instance, "classifiers_")
+            assert isinstance(estimator_instance.classifiers_, pd.DataFrame)
+            assert estimator_instance.classifiers_.shape == (1, 2)
diff --git a/sktime/classification/tests/test_base.py b/sktime/classification/tests/test_base.py
index 9713be6f0..33d0409f1 100644
--- a/sktime/classification/tests/test_base.py
+++ b/sktime/classification/tests/test_base.py
@@ -120,8 +120,6 @@ class _DummyConvertPandas(BaseClassifier):
 multivariate_message = r"multivariate series"
 missing_message = r"missing values"
 unequal_message = r"unequal length series"
-incorrect_X_data_structure = r"must be a np.array or a pd.Series"
-incorrect_y_data_structure = r"must be 1-dimensional"
 
 
 def test_base_classifier_fit():
@@ -144,8 +142,11 @@ def test_base_classifier_fit():
     test_X3 = _create_example_dataframe(cases=cases, dimensions=1, length=length)
     test_X4 = _create_example_dataframe(cases=cases, dimensions=3, length=length)
     test_y1 = np.random.randint(0, 2, size=(cases))
+    test_y2 = pd.DataFrame({"0": [1] * cases, "1": [0] * cases})
     result = dummy.fit(test_X1, test_y1)
     assert result is dummy
+    result = dummy.fit(test_X3, test_y2)
+    assert result is dummy
     with pytest.raises(ValueError, match=multivariate_message):
         result = dummy.fit(test_X2, test_y1)
     assert result is dummy
@@ -154,15 +155,6 @@ def test_base_classifier_fit():
     with pytest.raises(ValueError, match=multivariate_message):
         result = dummy.fit(test_X4, test_y1)
     assert result is dummy
-    # Raise a specific error if y is in a 2D matrix (1,cases)
-    test_y2 = np.array([test_y1])
-    # What if y is in a 2D matrix (cases,1)?
-    test_y2 = np.array([test_y1]).transpose()
-    with pytest.raises(ValueError, match=incorrect_y_data_structure):
-        result = dummy.fit(test_X1, test_y2)
-    # Pass a data fram
-    with pytest.raises(ValueError, match=incorrect_X_data_structure):
-        result = dummy.fit(test_X1, test_X3)
 
 
 TF = [True, False]
@@ -178,36 +170,41 @@ def test_check_capabilities(missing, multivariate, unequal):
     handle it and that cannot. Obvs could loop, but I think its clearer to just
     explicitly test;
     """
+    X_metadata = {
+        "has_nans": missing,
+        "is_univariate": not multivariate,
+        "is_equal_length": not unequal,
+    }
     handles_none = _DummyClassifier()
     handles_none_composite = _DummyComposite(_DummyClassifier())
 
     # checks that errors are raised
     if missing:
         with pytest.raises(ValueError, match=missing_message):
-            handles_none._check_capabilities(missing, multivariate, unequal)
+            handles_none._check_capabilities(X_metadata)
     if multivariate:
         with pytest.raises(ValueError, match=multivariate_message):
-            handles_none._check_capabilities(missing, multivariate, unequal)
+            handles_none._check_capabilities(X_metadata)
     if unequal:
         with pytest.raises(ValueError, match=unequal_message):
-            handles_none._check_capabilities(missing, multivariate, unequal)
+            handles_none._check_capabilities(X_metadata)
     if not missing and not multivariate and not unequal:
-        handles_none._check_capabilities(missing, multivariate, unequal)
+        handles_none._check_capabilities(X_metadata)
 
     if missing:
         with pytest.warns(UserWarning, match=missing_message):
-            handles_none_composite._check_capabilities(missing, multivariate, unequal)
+            handles_none_composite._check_capabilities(X_metadata)
     if multivariate:
         with pytest.warns(UserWarning, match=multivariate_message):
-            handles_none_composite._check_capabilities(missing, multivariate, unequal)
+            handles_none_composite._check_capabilities(X_metadata)
     if unequal:
         with pytest.warns(UserWarning, match=unequal_message):
-            handles_none_composite._check_capabilities(missing, multivariate, unequal)
+            handles_none_composite._check_capabilities(X_metadata)
     if not missing and not multivariate and not unequal:
-        handles_none_composite._check_capabilities(missing, multivariate, unequal)
+        handles_none_composite._check_capabilities(X_metadata)
 
     handles_all = _DummyHandlesAllInput()
-    handles_all._check_capabilities(missing, multivariate, unequal)
+    handles_all._check_capabilities(X_metadata)
 
 
 def test_convert_input():
@@ -228,26 +225,26 @@ def test_convert_input():
     test_X1 = np.random.uniform(-1, 1, size=(cases, length))
     test_X2 = np.random.uniform(-1, 1, size=(cases, 2, length))
     tester = _DummyClassifier()
-    tempX = tester._convert_X(test_X2)
+    tempX = tester._convert_X(test_X2, "numpy3D")
     assert tempX.shape[0] == cases and tempX.shape[1] == 2 and tempX.shape[2] == length
     instance_list = []
     for _ in range(0, cases):
         instance_list.append(pd.Series(np.random.randn(10)))
     test_X3 = _create_example_dataframe(cases=cases, dimensions=1, length=length)
     test_X4 = _create_example_dataframe(cases=cases, dimensions=3, length=length)
-    tempX = tester._convert_X(test_X3)
+    tempX = tester._convert_X(test_X3, "nested_univ")
     assert tempX.shape[0] == cases and tempX.shape[1] == 1 and tempX.shape[2] == length
-    tempX = tester._convert_X(test_X4)
+    tempX = tester._convert_X(test_X4, "nested_univ")
     assert tempX.shape[0] == cases and tempX.shape[1] == 3 and tempX.shape[2] == length
     tester = _DummyConvertPandas()
-    tempX = tester._convert_X(test_X2)
+    tempX = tester._convert_X(test_X2, "numpy3D")
     assert isinstance(tempX, pd.DataFrame)
     assert tempX.shape[0] == cases
     assert tempX.shape[1] == 2
     test_y1 = np.random.randint(0, 1, size=(cases))
     test_y1 = pd.Series(test_y1)
     tempX, tempY = _internal_convert(test_X1, test_y1)
-    assert isinstance(tempY, np.ndarray)
+    assert isinstance(tempY, pd.Series)
     assert isinstance(tempX, np.ndarray)
     assert tempX.ndim == 3
 
@@ -261,9 +258,10 @@ def test__check_classifier_input():
     4. Test incorrect: y as a list
     5. Test incorrect: too few cases or too short a series
     """
+    clf = _DummyClassifier()
 
     def _check_classifier_input(X, y=None, enforce_min_instances=1):
-        return BaseClassifier._check_classifier_input(None, X, y, enforce_min_instances)
+        return clf._check_input(X, y, enforce_min_instances)
 
     # 1. Test correct: X: np.array of 2 and 3 dimensions vs y:np.array and np.Series
     test_X1 = np.random.uniform(-1, 1, size=(5, 10))
diff --git a/sktime/forecasting/compose/tests/test_pipeline.py b/sktime/forecasting/compose/tests/test_pipeline.py
index d3b650dca..299cd7229 100644
--- a/sktime/forecasting/compose/tests/test_pipeline.py
+++ b/sktime/forecasting/compose/tests/test_pipeline.py
@@ -485,7 +485,7 @@ def test_forecastx_logic():
 )
 def test_forecastx_fit_behavior():
     from sktime.forecasting.compose import ForecastX
-    from sktime.forecasting.model_selection import temporal_train_test_split
+    from sktime.split import temporal_train_test_split
 
     y, X = load_longley()
     y_train, y_test, X_train, X_test = temporal_train_test_split(y, X)
diff --git a/sktime/forecasting/tests/test_exp_smoothing.py b/sktime/forecasting/tests/test_exp_smoothing.py
index fe1a2cf21..2d2633faa 100644
--- a/sktime/forecasting/tests/test_exp_smoothing.py
+++ b/sktime/forecasting/tests/test_exp_smoothing.py
@@ -1,11 +1,13 @@
 """Test exponential smoothing forecasters."""
 
-__author__ = ["mloning", "big-o"]
+__author__ = ["mloning", "big-o", "ciaran-g"]
 __all__ = ["test_set_params"]
 
+import pandas as pd
 import pytest
 from numpy.testing import assert_array_equal
 
+from sktime.datasets import load_airline
 from sktime.forecasting.exp_smoothing import ExponentialSmoothing
 from sktime.split import temporal_train_test_split
 from sktime.utils._testing.forecasting import make_forecasting_problem
@@ -35,3 +37,23 @@ def test_set_params():
     y_pred = f.predict()
 
     assert_array_equal(y_pred, expected)
+
+
+@pytest.mark.skipif(
+    not _check_soft_dependencies("statsmodels", severity="none"),
+    reason="skip test if required soft dependency not available",
+)
+def check_panel_expsmooth():
+    """Test exponential smoothing on panel data with datetime index."""
+    # make panel with hour of day panel and datetime index
+    y = load_airline()
+    y.index = pd.date_range(start="1960-01-01", periods=len(y.index), freq="H")
+    y.index.names = ["datetime"]
+    y.name = "passengers"
+    y = y.to_frame()
+    y["hour_of_day"] = y.index.hour
+    y = y.reset_index().set_index(["hour_of_day", "datetime"]).sort_index()
+
+    forecaster = ExponentialSmoothing(trend="add", sp=1)
+    forecaster.fit(y)
+    forecaster.predict(fh=[1, 3])
diff --git a/sktime/forecasting/tests/test_theta.py b/sktime/forecasting/tests/test_theta.py
index f6b6685b1..a31a22fd3 100644
--- a/sktime/forecasting/tests/test_theta.py
+++ b/sktime/forecasting/tests/test_theta.py
@@ -1,9 +1,10 @@
 """Tests for ThetaForecaster."""
 # copyright: sktime developers, BSD-3-Clause License (see LICENSE file)
 
-__author__ = ["big-o", "kejsitake"]
+__author__ = ["big-o", "kejsitake", "ciaran-g"]
 
 import numpy as np
+import pandas as pd
 import pytest
 
 from sktime.datasets import load_airline
@@ -121,3 +122,24 @@ def test_theta_and_thetamodular():
     y_pred_thetamodular = f1.predict(fh=fh)
 
     np.testing.assert_allclose(y_pred_theta, y_pred_thetamodular, rtol=0.06)
+
+
+@pytest.mark.skipif(
+    not _check_soft_dependencies("statsmodels", severity="none"),
+    reason="skip test if required soft dependency not available",
+)
+def check_panel_theta_quantiles():
+    """Test predict quantiles with theta on panel data with datetime index."""
+    # make panel with hour of day panel and datetime index
+    y = load_airline()
+    y.index = pd.date_range(start="1960-01-01", periods=len(y.index), freq="H")
+    y.index.names = ["datetime"]
+    y.name = "passengers"
+    y = y.to_frame()
+    y["hour_of_day"] = y.index.hour
+    y = y.reset_index().set_index(["hour_of_day", "datetime"]).sort_index()
+
+    forecaster = ThetaForecaster(sp=1)
+    forecaster.fit(y)
+    forecaster.predict(fh=[1, 3])
+    forecaster.predict_quantiles(fh=[1, 3], alpha=[0.1, 0.5, 0.9])
diff --git a/sktime/forecasting/trend/tests/test_pwl_trend.py b/sktime/forecasting/trend/tests/test_pwl_trend.py
index a1e4d938c..b5ab77109 100644
--- a/sktime/forecasting/trend/tests/test_pwl_trend.py
+++ b/sktime/forecasting/trend/tests/test_pwl_trend.py
@@ -11,11 +11,11 @@ import pytest
 
 from sktime.datasets import load_airline
 from sktime.forecasting.base import ForecastingHorizon
-from sktime.forecasting.model_selection import temporal_train_test_split
 from sktime.forecasting.trend import (
     PolynomialTrendForecaster,
     ProphetPiecewiseLinearTrendForecaster,
 )
+from sktime.split import temporal_train_test_split
 from sktime.tests.test_switch import run_test_for_class
 
 
diff --git a/sktime/regression/tests/test_all_regressors.py b/sktime/regression/tests/test_all_regressors.py
index 3dfea1768..c3e8b1c4a 100644
--- a/sktime/regression/tests/test_all_regressors.py
+++ b/sktime/regression/tests/test_all_regressors.py
@@ -4,10 +4,12 @@ __author__ = ["mloning", "TonyBagnall", "fkiraly"]
 
 
 import numpy as np
+import pandas as pd
 import pytest
 
 from sktime.datatypes import check_is_scitype
 from sktime.tests.test_all_estimators import BaseFixtureGenerator, QuickTester
+from sktime.utils._testing.panel import make_regression_problem
 from sktime.utils._testing.scenarios_classification import (
     ClassifierFitPredictMultivariate,
 )
@@ -80,3 +82,25 @@ class TestAllRegressors(RegressorFixtureGenerator, QuickTester):
         assert isinstance(y_pred, np.ndarray)
         assert y_pred.shape == (X_new_instances,)
         assert np.issubdtype(y_pred.dtype, np.floating)
+
+    def test_multioutput(self, estimator_instance):
+        """Test multioutput regression for all classifiers.
+
+        All classifiers should follow the same interface,
+        those that do not genuinely should vectorize/broadcast over y.
+        """
+        n_instances = 20
+        X, y = make_regression_problem(n_instances=n_instances)
+        y_mult = pd.DataFrame({"a": y, "b": y})
+
+        estimator_instance.fit(X, y_mult)
+        y_pred = estimator_instance.predict(X)
+
+        assert isinstance(y_pred, pd.DataFrame)
+        assert y_pred.shape == y_mult.shape
+
+        vectorized = estimator_instance.get_tag("capability:multioutput")
+        if vectorized:
+            assert hasattr(estimator_instance, "regressors_")
+            assert isinstance(estimator_instance.regressors_, pd.DataFrame)
+            assert estimator_instance.regressors_.shape == (1, 2)
diff --git a/sktime/split/temporal_train_test_split.py b/sktime/split/temporal_train_test_split.py
index 2509b0ec9..55e7e9195 100644
--- a/sktime/split/temporal_train_test_split.py
+++ b/sktime/split/temporal_train_test_split.py
@@ -91,8 +91,8 @@ def temporal_train_test_split(
 
     Examples
     --------
-    >>> from sktime.forecasting.model_selection import temporal_train_test_split
     >>> from sktime.datasets import load_airline, load_osuleaf
+    >>> from sktime.split import temporal_train_test_split
     >>> from sktime.utils._testing.panel import _make_panel
     >>> # univariate time series
     >>> y = load_airline()
diff --git a/sktime/split/testplustrain.py b/sktime/split/testplustrain.py
index 73981a9c6..a803a52c9 100644
--- a/sktime/split/testplustrain.py
+++ b/sktime/split/testplustrain.py
@@ -131,10 +131,7 @@ class TestPlusTrainSplitter(BaseSplitter):
             `MyClass(**params)` or `MyClass(**params[i])` creates a valid test instance.
             `create_test_instance` uses the first (or only) dictionary in `params`
         """
-        from sktime.forecasting.model_selection import (
-            ExpandingWindowSplitter,
-            SingleWindowSplitter,
-        )
+        from sktime.split import ExpandingWindowSplitter, SingleWindowSplitter
 
         cv_1 = ExpandingWindowSplitter(fh=[2, 4], initial_window=24, step_length=12)
         cv_2 = SingleWindowSplitter(fh=[2, 4], window_length=24)
diff --git a/sktime/tests/_config.py b/sktime/tests/_config.py
index 14435f3c4..7b7783d7d 100644
--- a/sktime/tests/_config.py
+++ b/sktime/tests/_config.py
@@ -77,6 +77,7 @@ EXCLUDED_TESTS = {
         "test_persistence_via_pickle",
         "test_fit_does_not_overwrite_hyper_params",
         "test_save_estimators_to_file",
+        "test_multiprocessing_idempotent",  # see 5658
     ],
     "ProximityForest": [
         "test_persistence_via_pickle",
@@ -182,6 +183,7 @@ EXCLUDED_TESTS = {
     "LTSFLinearForecaster": ["test_predict_time_index_in_sample_full"],
     "LTSFDLinearForecaster": ["test_predict_time_index_in_sample_full"],
     "LTSFNLinearForecaster": ["test_predict_time_index_in_sample_full"],
+    "WEASEL": ["test_multiprocessing_idempotent"],  # see 5658
 }
 
 # We use estimator tags in addition to class hierarchies to further distinguish
diff --git a/sktime/transformations/tests/test_all_transformers.py b/sktime/transformations/tests/test_all_transformers.py
index 493fb63d3..b6a9699ac 100644
--- a/sktime/transformations/tests/test_all_transformers.py
+++ b/sktime/transformations/tests/test_all_transformers.py
@@ -78,6 +78,8 @@ class TestAllTransformers(TransformerFixtureGenerator, QuickTester):
                 return "Panel"
             if X_scitype in ["Panel", "Hierarchical"]:
                 return "Hierarchical"
+        if trafo_input == "Panel" and trafo_output == "Series":
+            return "Series"
 
     def test_fit_transform_output(self, estimator_instance, scenario):
         """Test that transform output is of expected scitype."""
diff --git a/sktime/utils/_testing/estimator_checks.py b/sktime/utils/_testing/estimator_checks.py
index a275d8d64..ff56eb21f 100644
--- a/sktime/utils/_testing/estimator_checks.py
+++ b/sktime/utils/_testing/estimator_checks.py
@@ -20,6 +20,7 @@ from sktime.clustering.base import BaseClusterer
 from sktime.datatypes._panel._check import is_nested_dataframe
 from sktime.dists_kernels import BasePairwiseTransformer, BasePairwiseTransformerPanel
 from sktime.forecasting.base import BaseForecaster
+from sktime.registry import scitype
 from sktime.regression.base import BaseRegressor
 from sktime.tests._config import VALID_ESTIMATOR_TYPES
 from sktime.transformations.base import BaseTransformer
@@ -216,6 +217,11 @@ def _has_capability(est, method: str) -> bool:
             return True
         return get_tag(est, "capability:pred_int", False)
     # skip transform for forecasters that have it - pipelines
-    if method == "transform" and isinstance(est, BaseForecaster):
+    if method == "transform" and scitype(est) in (
+        "classifier",
+        "forecaster",
+    ):
+        return False
+    if method == "predict" and scitype(est) == "transformer":
         return False
     return True
diff --git a/sktime/utils/_testing/scenarios_classification.py b/sktime/utils/_testing/scenarios_classification.py
index ae12056bc..dc34b2173 100644
--- a/sktime/utils/_testing/scenarios_classification.py
+++ b/sktime/utils/_testing/scenarios_classification.py
@@ -14,9 +14,7 @@ __all__ = [
 from inspect import isclass
 
 from sktime.base import BaseObject
-from sktime.classification.base import BaseClassifier
-from sktime.classification.early_classification import BaseEarlyClassifier
-from sktime.regression.base import BaseRegressor
+from sktime.registry import scitype
 from sktime.utils._testing.hierarchical import _make_hierarchical
 from sktime.utils._testing.panel import _make_classification_y, _make_panel_X
 from sktime.utils._testing.scenarios import TestScenario
@@ -71,12 +69,14 @@ class ClassifierTestScenario(TestScenario, BaseObject):
             else:
                 return obj.get_tag(tag_name)
 
-        regr_or_classf = (BaseClassifier, BaseEarlyClassifier, BaseRegressor)
-
         # applicable only if obj inherits from BaseClassifier, BaseEarlyClassifier or
         #   BaseRegressor. currently we test both classifiers and regressors using these
         #   scenarios
-        if not isinstance(obj, regr_or_classf) and not issubclass(obj, regr_or_classf):
+        if not scitype(obj) in (
+            "classifier",
+            "early_classifier",
+            "regressor",
+        ):
             return False
 
         # if X is multivariate, applicable only if can handle multivariate
diff --git a/sktime/utils/_testing/scenarios_forecasting.py b/sktime/utils/_testing/scenarios_forecasting.py
index 2dc06aecd..9281153fd 100644
--- a/sktime/utils/_testing/scenarios_forecasting.py
+++ b/sktime/utils/_testing/scenarios_forecasting.py
@@ -18,7 +18,7 @@ import pandas as pd
 
 from sktime.base import BaseObject
 from sktime.datatypes import mtype_to_scitype
-from sktime.forecasting.base import BaseForecaster
+from sktime.registry import scitype
 from sktime.utils._testing.hierarchical import _make_hierarchical
 from sktime.utils._testing.panel import _make_panel_X
 from sktime.utils._testing.scenarios import TestScenario
@@ -49,7 +49,7 @@ class ForecasterTestScenario(TestScenario, BaseObject):
                 return obj.get_tag(tag_name)
 
         # applicable only if obj inherits from BaseForecaster
-        if not isinstance(obj, BaseForecaster) and not issubclass(obj, BaseForecaster):
+        if scitype(obj) != "forecaster":
             return False
 
         # applicable only if number of variables in y complies with scitype:y