[ENH] Add multirocket multivariate cython implementation

[leandrepi]

Reference Issues/PRs

See also sktime/sktime#11292 regarding the original_implementation flag.

What does this implement/fix? Explain your changes.

This PR adds a cython estimator for multirocket (multivariate). In doing so it tries to reuse identical code from the minirocket implementation by @sssilvar (fitting the dilations, and computing the per-combination quantiles, which I've rewritten in a vectorized fashion, along with a non-regression test).

This cython implem has the advantage of accepting either one of float32/float64 inputs, without enforcing a cast, relying on the cython floating dtype.

Does your contribution introduce a new dependency? If yes, which one?

No

What should a reviewer concentrate their feedback on?

Probably the conjunction of this PR and the one related to the multirocket bug in the current sktime implementation. I'm not quite sure whether merging this PR before the sktime one previously referenced is the correct way to go; upon merge of this PR, we'd also need to wrap the estimator in the main repo, as was done for _minirocket_cython.

Any other comments?

PR checklist

For all contributions

• I've added unit tests and made sure they pass locally.

For new estimators

• I've added the estimator to the online documentation.

[sssilvar]

It's more of a side note after testing aaround (@fkiraly):

The output array is sized from this parameter, and the differenced block starts at half that width (.pyx L272-276):

num_total_features = (num_features + num_features1) * num_features_per_kernel
num_features_per_transform = num_total_features >> 1

_one ignores the parameter and always writes 4 features per kernel, one block of num_features columns each (.pyx L176-183). With k = 2, the array is 4 nf wide, and the differenced pass writes up to 6 nf. Each row writes into the next one, and the last row writes past the end of the array. Bounds checking is off, so nothing stops it. With k = 8, the array is
16 nf wide but only 8 nf columns are written.

sktime has the same code (_multirocket_multi_numba.py L703-707: the same allocation, the same features.shape[1] / 2 offset, 4 writes per kernel):

X = rng.normal(size=(64, 3, 500))

# this PR
p = multirocket_fit(X, num_kernels=10_000, random_state=0)
for _ in range(50):
    multirocket_transform(X, p, n_features_per_kernel=2)
# Out: Segmentation fault (exit 139)

# sktime
t = MultiRocketMultivariate(num_kernels=10_000, n_features_per_kernel=2, random_state=0).fit(X)
for _ in range(50):
    t.transform(X)
# Out: Segmentation fault (exit 139)

[leandrepi]

Yeah that's a bit silly, this shouldn't be a parametrizable value at all... In the original implem it's set to 4 within the __init__ block but it's not an arg of the method. There are always 4 features per kernel.

[sssilvar]

LGTM @fkiraly Let me know if we can merge.

[fkiraly]

quick question, is this stacked on #6?

[fkiraly]

High-level review: can we leave setup.py to accept explicit modules? There might be other modules than just rocket ones, e.g., another contributor could add cython code pertinent for another sktime module.

[sssilvar]

quick question, is this stacked on #6?

No, both #6 and #7 are rebased on main. I'll rebase it if we merge #6 first. Order should not matter here.

High-level review: can we leave setup.py to accept explicit modules? There might be other modules than just rocket ones, e.g., another contributor could add cython code pertinent for another sktime module.

Do you mean something like this?

_MODULES = [
    "sktime_cython.transformations.rocket._minirocket_multivariate_cython",
    "sktime_cython.transformations.rocket._multirocket_multivariate_cython",
]

Cause otherwise Cython accepts something like:

extensions = [
    Extension(
        "*",
        ["sktime_cython/**/*.pyx"],
        include_dirs=[np.get_include()],
        define_macros=[("NPY_NO_DEPRECATED_API", "NPY_1_7_API_VERSION")],
        extra_compile_args=_fast,
    )
]

And it should do it automatically.