I'm currently exploring a couple of statistical forecasting methods. My problem rose when using VARMA(2,1) fixed window rolling forecast. The example code that I'm using is the following: Here I only used 2 series out of 24. import pandas as pd from statsmodels.tsa.statespace.varmax import VARMAX # Parameters window = 90 # Prepare the two series and train/test split series_df = dfp[['H0', 'H1']] n_test = int(len(series_df) * 0.1) n_train = len(series_df) - n_test test_index = series_df.index[n_train:] forecasts = pd.DataFrame(index=test_index, columns=series_df.columns, dtype=float) # Rolling one‐step‐ahead forecasts using VARMAX(2,1) for i, dt in enumerate(test_index): end = n_train + i start = end - window history_df = series_df.iloc[start:end] model = VARMAX( history_df, order=(2, 1), # AR=2, MA=1 trend='c', enforce_stationarity=False, enforce_invertibility=False ) res = model.fit(maxiter=10, disp=False) # forecast returns a DataFrame with one row and two columns fc = res.forecast(steps=1) forecasts.iloc = fc.iloc[0] The forecast window size is 90 days, and when applying it to the test set, we talk about 131 1-step-ahead forecasts. For my application it is necessary that I take the full 24-dimensional dataframe. My question now is, is there a way to speed this up. I'm running it on the 2 series in a colab environment and it takes a minute or so. Using the entire 24 dimensional dataframe, after hours the forecast hasnt been concluded yet. Not only that I tried in Python but also in R too, and encountered the same issue. I also tried to give the VARMAX function a warmstart, by feeding the prior estimates into each new window, also without success. Another question that I have is, I've investigated the internal construction of the VARMAX function of statsmodels, and by default, any VARMA model is internally transferred into a linear state space model. Would it be possible to talk for example the tensorflow time series libraries, that have a state space model available that runs on GPU and it would lead to the same (or very similar results) ? Has anyone an idea? As I said, when doing the very reduced VARMA(2,1) with the prior code, everything worked and finished soon-ish. I tried the full data VARMA with many different settings (low tolerance, low maxiter, low_memeory, different optimizer for estimation the MA-part of MLE, ...) but nothing seems to work to finish the rolling forecast in any manner of time. Continue reading...