I need to apply the coint function from the statsmodels library to 207 times series with 1397 points each, two by two.

Currently, it takes between 35-40 minutes on my computer with an Intel 24 Cores CPU, last generation.

I tried to use Cython, or Data processing hacks from this article but I get the exact same processing time.

Here is the code to reproduce it:

# Data generation (no improvement needed)
df_timeseries = pd.DataFrame(np.random.uniform(low=2.25, high=2784.07,size=(1397, 207)))
# Downcasting to float 16 (could be unsigned too)
df_timeseries[df_timeseries.columns] = df_timeseries[df_timeseries.columns].astype(np.float16)

The cointegration function requires two time series, so I build a permutation of time series data frame column names:

from more_itertools import distinct_permutations as idp
# Shape: 42642 rows, 2 columns
df_permut = pd.DataFrame(idp(df_timeseries.columns, 2), columns=['ts1', 'ts2'])

Then, I apply the coint function to the permutation dataframe and extract only the p-value from the return (coint function returns coint_t, pvalue, crit_value):

import statsmodels.tsa.stattools as st
df_permut["pvalue"] = df_permut.apply(
                        lambda x: 
                            )][1], axis=1)

This last part takes about 40 mins to run. I know that statsmodels is quite optimized and there is no chance I can imporved the coint method (I did try by extracting the code and removing the numerous unwanted check, but the linear regression runs from statsmodels take most of the time) I don't believe there is many places for improvements on my code as the coint method is the most resource greedy.

How to drastically improve the speed? If there is no way, is that a path to move the coint method to GPU?


1 Answer 1


Your code is short, clear, and time consuming.

How to drastically improve the speed?

You need to compute fewer figures. For N time series you perform O(N^2) cointegration tests.

The OP does not describe the use case nor the observed data pattern. Fortunately you do know what pattern the results showed, since you patiently waited 40 minutes for the exhaustive comparison.

break symmetry: 2x

Currently for each TS pair (a, b) you also test (b, a). Consult your corpus of results to see if that's really necessary. Many pairs, of course, will show no cointegration, in either direction, showing no "hit". Is there anything you know about each time series that lets you predict in which direction we might see a hit?

For example, given a time series of Close prices you're likely to also know its corresponding volume, volatility, vertical segment, and which market it trades on. Given hits from your results corpus, there's an opportunity to do K-means clustering to identify predictive features.

autolag: 10x

By default we optimize an information criterion: AIC. Running coint with autolag=None skips that, so it runs an order of magnitude quicker. Consult your results corpus to see if AIC produces an important difference in the filtering result.

common sense

Your results corpus might show for example that a high volume stock will only forecast the price of a low volume stock, or vice versa. Or that only stocks within the same vertical have predictive power for one another. Prune your search space accordingly.

simple features

The quadratic search is killing performance. You want to be able to process each time series exactly once to extract relevant features, and then perform simple all-pairs comparisons. One candidate is Change Point Detection.

Sometimes an exogenous news event will affect a vertical market, inducing immediate or delayed effects in small or large players. Set parameter K to a conveniently small number, like 1. Identify the K most likely Change Points in each time series, sort them, and scan the sorted list of timestamps to identify candidate pairs of time series to analyze.


Prune the search space.

There is knowledge, from the real world and from your results corpus, that you can bring into this problem so you don't perform a quadratic number of expensive analyses.

  • \$\begingroup\$ Many thanks for the detailed answer. \$\endgroup\$
    – Begoodpy
    Jun 2, 2023 at 19:07
  • \$\begingroup\$ break symmetry: 2x from the following q&a I understood that cointegration is non-directional but there might be cases where it can be, thus I test both (a, b) and (b, a): stats.stackexchange.com/questions/300209/… \$\endgroup\$
    – Begoodpy
    Jun 2, 2023 at 19:09
  • \$\begingroup\$ autolag Just tested it with None and indeed it is 10 times faster, but the results are far too different \$\endgroup\$
    – Begoodpy
    Jun 2, 2023 at 19:10
  • \$\begingroup\$ Oh, that PCA remark offers an interesting possibility! Yeah, I had a feeling that lag details would show up differently in your source data than in my random data. The comment about "searching for optimal lag is expensive" still stands, though. Perhaps AIC finds just one or a handful of lags that matter to you, such as 7 and 30 ? You might be able to do a cheaper autolag search. \$\endgroup\$
    – J_H
    Jun 2, 2023 at 19:37
  • \$\begingroup\$ Thanks. After checking, the autolag is indeed the most greedy function. I improved it by computing just the aic value, splitting the processing time by two. \$\endgroup\$
    – Begoodpy
    Jun 3, 2023 at 8:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.