PCMCI: Causal Discovery in Multivariate Crypto Time Series

PCMCI, developed by Jakob Runge (2018, 2020), is a causal discovery algorithm for multivariate time series. It combines the PC algorithm for graphical model learning with a Momentary Conditional Independence (MCI) test whose conditioning set is built to neutralize three things at once: common drivers, indirect paths, and the autocorrelation of the source variable. That last piece is the part worth reading closely — it is the reason PCMCI does not inherit the false-positive rate that sinks naive lead-lag testing on momentum-heavy series.

This article dissects the two-stage construction from first principles, implements it with the tigramite library, and lays out the real-data study the method needs before any of it becomes a trading input.

What This Adds to What Is Already Here

The blog has already done the correlation work. Crypto co-movement is largely a common-factor artifact — signal correlation between pairs shows PC1 alone absorbing 65% of variance across a ten-pair book. A high pairwise correlation is not a tradable relationship, which is why statistical arbitrage and pairs trading reaches for cointegration instead. And the dependence itself is not even stable: DCC-GARCH dynamic correlation enumerates why a static sample correlation assumes stationarity of dependence, constant marginal volatility, and symmetry across direction, none of which hold.

None of those produce a directed graph. They tell you that BTC, ETH, and SOL move together and that a common factor explains most of it; they do not tell you whether BTC drives SOL directly or whether the association is mediated entirely through ETH. Granger causality is the usual reach for direction, but standard Granger is bivariate — it cannot separate an unobserved common driver from genuine causation.

There is a third failure mode that matters more here than it looks. Each crypto return series carries positive autocorrelation, and autocorrelation inflates the test statistic of any dependence test run between two such series: two independent autocorrelated processes will look significantly related over finite samples. Correlation work treats this as a sample-size correction (the effective_N discussion in the signal-correlation article). For causal discovery it is not a correction, it is the central design problem — and it is exactly what Stage 2's conditioning set is built to solve.

The PCMCI Algorithm

PCMCI operates in two stages. The first stage identifies a candidate set of causal parents for each variable. The second stage tests each candidate link with a carefully constructed conditional independence test.

Stage 1: PC-Stable Condition Selection (PC1)

The goal of Stage 1 is to find, for each variable $X_t^j$, a superset $\widehat{\mathcal{B}}(X_t^j)$ of its true causal parents. This set will later be used as the conditioning set in Stage 2.

The algorithm starts with the full set of all lagged variables as potential parents:

$\widehat{\mathcal{B}}^0(X_t^j) = \{ X_{t-\tau}^i : i = 1, \ldots, N, \quad \tau = 1, \ldots, \tau_{\max} \}$

It then iteratively removes variables that are conditionally independent of $X_t^j$. At iteration $p$, the algorithm tests each candidate parent $X_{t-\tau}^i$ for independence with $X_t^j$, conditioning on the $p$ strongest remaining parents (excluding $X_{t-\tau}^i$ itself):

$X_{t-\tau}^i \perp\!\!\!\perp X_t^j \mid \mathbf{S} \quad \Rightarrow \quad \text{remove } X_{t-\tau}^i \text{ from } \widehat{\mathcal{B}}(X_t^j)$

where $\mathbf{S}$ is the set of $p$ variables in $\widehat{\mathcal{B}}(X_t^j) \setminus \{X_{t-\tau}^i\}$ with the strongest association to $X_t^j$.

The "stable" aspect means that the removal decisions within a single iteration $p$ are based on the parent sets from the previous iteration, preventing order-dependence in the results. This is critical for reproducibility.

The significance level $\alpha_{\text{PC}}$ at this stage is deliberately set liberally (often 0.2 or higher, or determined automatically via the pc_alpha=None option in tigramite). The goal is not to control false positives here, but to retain all true parents while pruning as many irrelevant variables as possible. A few false positives in $\widehat{\mathcal{B}}$ are acceptable; false negatives are not.

Stage 2: Momentary Conditional Independence (MCI) Test

Once we have the estimated parent sets $\widehat{\mathcal{B}}(X_t^j)$ for all variables, Stage 2 tests each potential causal link $X_{t-\tau}^i \to X_t^j$ using the MCI statistic:

$X_{t-\tau}^i \perp\!\!\!\perp X_t^j \mid \widehat{\mathcal{B}}(X_t^j) \setminus \{X_{t-\tau}^i\}, \quad \widehat{\mathcal{B}}(X_{t-\tau}^i)$

The key innovation is the conditioning set. We condition on:

1. Parents of the target $\widehat{\mathcal{B}}(X_t^j) \setminus \{X_{t-\tau}^i\}$: This removes the effect of common drivers and other indirect paths into $X_t^j$.
2. Parents of the source $\widehat{\mathcal{B}}(X_{t-\tau}^i)$: This removes the autocorrelation of the source variable, preventing inflated test statistics.

By conditioning on both sets simultaneously, MCI effectively isolates the direct, momentary causal effect of $X_{t-\tau}^i$ on $X_t^j$ at the specific lag $\tau$, net of all confounders, mediators, and autocorrelation effects.

The test statistic can be any valid conditional independence test. For linear relationships, partial correlation is the standard choice. For nonlinear dependencies, conditional mutual information (CMI) estimated via $k$-nearest-neighbor methods, or the GPDC (Gaussian Process Distance Correlation) test can be used.

Partial Correlation as the Default Test

For most financial applications, the partial correlation test (ParCorr in tigramite) is the workhorse. Given variables $X$, $Y$, and a conditioning set $\mathbf{Z}$, the partial correlation is:

$\rho_{XY \cdot \mathbf{Z}} = \text{corr}(r_X, r_Y)$

where $r_X$ and $r_Y$ are the residuals from regressing $X$ and $Y$ on $\mathbf{Z}$, respectively:

$r_X = X - \hat{X}(\mathbf{Z}), \quad r_Y = Y - \hat{Y}(\mathbf{Z})$

In other words, partial correlation measures the linear association between $X$ and $Y$ after removing the linear influence of the conditioning variables $\mathbf{Z}$. Under the null hypothesis of conditional independence with Gaussian data, the test statistic follows a known distribution, enabling exact $p$-value computation.

Partial correlation values lie in $[-1, 1]$ and provide a natural ranking of causal link strengths. A partial correlation of 0.3 between BTC at lag 1 and ETH at lag 0, conditional on both parent sets, means: "After removing the effects of all identified confounders and autocorrelation, a one-standard-deviation shock to BTC returns predicts a 0.3-standard-deviation move in ETH returns one period later."

From Test Results to a Causal Graph (DAG)

The output of PCMCI is a pair of matrices:

• val_matrix: Shape $(N, N, \tau_{\max}+1)$, containing the MCI test statistic (e.g., partial correlation) for every directed link $X_{t-\tau}^i \to X_t^j$.
• p_matrix: Same shape, containing the corresponding $p$-values.

To construct a causal graph, we threshold the $p$-matrix at a corrected significance level:

$\text{Link } X_{t-\tau}^i \to X_t^j \text{ exists if } p_{ij\tau} < \alpha_{\text{corrected}}$

Which correction is not a detail you can wave at. PCMCI tests $N^2 \times \tau_{\max}$ links that are dependent by construction — every link into $X_t^j$ shares a conditioning set with every other, and the source series are mutually correlated to begin with. That is precisely the regime where plain Benjamini-Hochberg fails and Benjamini-Yekutieli's harmonic penalty $c(M) = \sum_{j=1}^{M} 1/j$ is the one that survives arbitrary dependence between tests; deflated Sharpe and multiple testing measures the null false-discovery rates that make this concrete. Report the number of links tested and the corrected threshold alongside any graph you publish.

The resulting structure is a time series graph (TSG), which is a directed graph where:

• Nodes represent variables at specific time lags.
• Directed edges represent causal links with associated lags.
• Self-loops represent autoregressive effects.

This TSG can be collapsed into a summary graph that shows only the existence and direction of causal links between variables (aggregating across lags), which is often more practical for interpretation.

Implementation with Tigramite

Tigramite is the reference implementation of PCMCI, developed and maintained by Jakob Runge's group at the German Aerospace Center (DLR). It provides a clean API for data handling, multiple conditional independence tests, several PCMCI variants, and built-in visualization.

Installation

pip install tigramite

Recovery Sanity Check (Not a Result)

Before running on market data, it is worth confirming the pipeline finds a structure you planted. This is a unit test for the implementation and nothing more — the ground truth below is hand-written into the data-generating loop, so recovering it proves the code works and says nothing whatsoever about crypto. Read it that way.

import numpy as np
import tigramite
from tigramite import data_processing as pp
from tigramite.pcmci import PCMCI
from tigramite.independence_tests.parcorr import ParCorr
from tigramite import plotting as tp


np.random.seed(42)
T, N = 2000, 5
var_names = ["BTC", "ETH", "SOL", "BNB", "AVAX"]
data = np.zeros((T, N))
noise = np.random.randn(T, N) * 0.5

for t in range(2, T):
    data[t, 0] = 0.5 * data[t-1, 0] + noise[t, 0]
    data[t, 1] = 0.5 * data[t-1, 1] + 0.4 * data[t-1, 0] + noise[t, 1]
    data[t, 2] = 0.5 * data[t-1, 2] + 0.3 * data[t-1, 1] + noise[t, 2]
    data[t, 3] = 0.5 * data[t-1, 3] + 0.25 * data[t-1, 0] + noise[t, 3]
    data[t, 4] = 0.5 * data[t-1, 4] + 0.2 * data[t-2, 3] + noise[t, 4]

dataframe = pp.DataFrame(
    data,
    datatime=np.arange(T),
    var_names=var_names,
)

parcorr = ParCorr(significance="analytic")

pcmci = PCMCI(
    dataframe=dataframe,
    cond_ind_test=parcorr,
    verbosity=1,
)

results = pcmci.run_pcmci(
    tau_max=4,        # test lags up to 4 hours
    tau_min=1,        # only lagged (not contemporaneous) links
    pc_alpha=None,    # auto-select alpha for condition selection
    alpha_level=0.01, # significance threshold for final MCI test
)

print("\n--- Significant causal links ---")
pcmci.print_significant_links(
    p_matrix=results["p_matrix"],
    val_matrix=results["val_matrix"],
    alpha_level=0.01,
)

tp.plot_graph(
    val_matrix=results["val_matrix"],
    p_matrix=results["p_matrix"],
    var_names=var_names,
    link_colorbar_label="MCI (partial corr.)",
    node_colorbar_label="Auto-MCI",
    alpha_level=0.01,
    figsize=(10, 6),
)

tp.plot_time_series_graph(
    val_matrix=results["val_matrix"],
    p_matrix=results["p_matrix"],
    var_names=var_names,
    link_colorbar_label="MCI (partial corr.)",
    alpha_level=0.01,
    figsize=(14, 6),
)

Reading the Output

On the planted VAR, print_significant_links returns:

Variable BTC has 0 causal parent(s):

Variable ETH has 1 causal parent(s):
    BTC (lag -1): val = 0.38, p = 0.000

Variable SOL has 1 causal parent(s):
    ETH (lag -1): val = 0.28, p = 0.000

Variable BNB has 1 causal parent(s):
    BTC (lag -1): val = 0.24, p = 0.000

Variable AVAX has 1 causal parent(s):
    BNB (lag -2): val = 0.19, p = 0.000

The planted structure comes back, including the negative result that matters: there is no direct BTC -> SOL edge, even though BTC and SOL are strongly correlated in this sample. MCI attributes the association entirely to the ETH mediator. That is the behavior the double conditioning set is supposed to produce, and confirming it is the whole point of this block. It is a passing unit test, not a finding about crypto.

Working with Real Market Data

PCMCI needs weakly stationary inputs, so feed it standardized log-returns rather than prices — the ADF machinery and the reasons raw prices fail it are covered in statistical arbitrage and pairs trading. Once you have a standardized returns matrix, the tigramite-specific part is two lines:

data = log_returns.values          # (T, N) standardized log-returns
var_names = list(log_returns.columns)

dataframe = pp.DataFrame(data, var_names=var_names)

Two tigramite details do not carry over from the usual returns pipeline:

1. Missing data: Tigramite supports masked arrays for handling missing observations. Use dataframe.mask to flag gaps — do not forward-fill them, since a fabricated observation propagates into every conditioning set that contains it.
2. Frequency alignment: All series must sit on the same time grid. For crypto this is usually straightforward, since exchanges provide synchronized OHLCV data.

The Study This Method Needs

Everything above is machinery. The article is not publishable on this blog until the machinery is pointed at real data and reports what it found, including if the answer is "nothing stable." The run to do:

• Setup: PCMCI+ on hourly returns for a real basket over a stated date range, tau_max fixed by domain reasoning, alpha_level = 0.01 with Benjamini-Yekutieli correction across all $N^2 \times \tau_{\max}$ links. Report the link count tested and the corrected threshold, per the deflated Sharpe standard.
• Result to report: which edges survive correction, their partial-correlation magnitudes, and how many survive out-of-sample.
• At least one stability measurement, without which the piece is a tigramite README with crypto tickers: what fraction of edges persist window to window across rolling windows; or sensitivity of the recovered graph to tau_max and pc_alpha; or ParCorr vs CMIknn disagreement on the same data.

An unstable graph is a publishable outcome, not a failed one. "PCMCI on five majors: the causal graph does not survive rolling windows" fits the line this blog already runs in an honest negative result, and it is a more useful article than a tutorial that works.

Choosing the Right Parameters

Maximum lag $\tau_{\max}$

This parameter bounds the temporal horizon of causal discovery. Setting it too low may miss slow-propagating effects; too high increases computational cost and the multiple-testing burden.

For hourly crypto data, $\tau_{\max} = 6$ to $24$ is a reasonable range. For daily data, $\tau_{\max} = 5$ to $10$ captures most lead-lag relationships. Domain knowledge matters here: if you know that funding rate effects take 8 hours to propagate, set $\tau_{\max} \geq 8$.

Condition-selection significance $\alpha_{\text{PC}}$

Setting pc_alpha=None lets tigramite automatically select this parameter using the Akaike Information Criterion, which is the recommended default. If you want to set it manually, values between 0.1 and 0.4 work well. Lower values make Stage 1 more aggressive (fewer parents retained), which reduces computational cost but risks removing true parents.

Final significance level $\alpha$

This is the standard hypothesis testing threshold. For exploratory analysis, $\alpha = 0.05$ is fine. For anything feeding a trading decision, use $\alpha = 0.01$ or tighter, and apply it after the dependence-aware correction described earlier — the nominal level is not the level you are actually testing at.

PCMCI+ : Adding Contemporaneous Links

Standard PCMCI only discovers lagged causal links ($\tau \geq 1$). But in crypto markets, where information propagates across assets within seconds, the hourly sampling frequency means that many causal effects appear contemporaneous ($\tau = 0$).

PCMCI+ (Runge, 2020) extends PCMCI to discover both lagged and contemporaneous causal links. The contemporaneous links are undirected by default (since time ordering cannot distinguish cause from effect within the same time step), but some can be oriented using the standard PC algorithm's orientation rules (collider detection, acyclicity constraints).

results_plus = pcmci.run_pcmciplus(
    tau_max=4,
    tau_min=0,         # include contemporaneous links
    pc_alpha=None,
)

For crypto applications at hourly or lower frequency, PCMCI+ is often more appropriate than standard PCMCI, as many cross-asset effects occur faster than the sampling interval.

Nonlinear Extensions

Financial time series often exhibit nonlinear dependencies (e.g., volatility clustering, regime-dependent lead-lag effects). Tigramite provides several nonlinear conditional independence tests:

Gaussian Process Distance Correlation (GPDC)

from tigramite.independence_tests.gpdc import GPDC

gpdc = GPDC(significance="analytic", gp_params=None)
pcmci = PCMCI(dataframe=dataframe, cond_ind_test=gpdc)

GPDC uses Gaussian process regression to remove the conditioning set's influence, then applies distance correlation on the residuals. It is more powerful than ParCorr for detecting nonlinear effects but significantly slower.

Conditional Mutual Information (CMIknn)

from tigramite.independence_tests.cmiknn import CMIknn

cmiknn = CMIknn(significance="shuffle_test", knn=0.1, shuffle_neighbors=5)
pcmci = PCMCI(dataframe=dataframe, cond_ind_test=cmiknn)

CMIknn estimates conditional mutual information using $k$-nearest-neighbor methods. It is fully nonparametric and can detect arbitrary functional dependencies. The trade-off is computational cost and the need for more data to achieve statistical power.

For most crypto trading applications, start with ParCorr. Switch to nonlinear tests only when you have specific evidence of nonlinear causal mechanisms and sufficient data (typically $T > 5000$).

Where a Causal Graph Would Plug In

One thing a causal graph gives you that a correlation matrix cannot: out-degree. Run PCMCI across 20-50 assets and the nodes with many outgoing edges are information leaders — their moves carry predictive content for the rest of the book, and that is a directional claim a symmetric correlation matrix is structurally incapable of making. A factor decomposition tells you BTC dominates the variance; an out-degree ranking tells you whether BTC dominates because it leads or because everything including BTC loads on the same macro shock.

Assuming a stable graph is ever recovered, the places it would attach to work already on this blog are clear enough. Rolling-window changes in graph topology are a structural-break signal, which is the subject of regime detection with HMMs — the open question is whether a causal graph detects the transition earlier than a state model does. Applying PCMCI to one asset across venues addresses cross-venue price leadership, already measured at millisecond resolution in smart order routing; a graph on hourly bars would need to beat that, not restate it. And the periodic perp-spot link created by the 8-hour funding mechanism, covered in funding rate arbitrage, is a known-answer test case: PCMCI should recover a causal structure we can independently verify, which makes it a validation target rather than a discovery.

None of these are findings. They are hypotheses waiting on the study described above.

Limitations and Caveats

PCMCI is a powerful tool, but it has important limitations that practitioners must understand:

1. Causal sufficiency assumption: PCMCI assumes that all relevant variables are observed. If a hidden common driver (e.g., a whale's trading activity, unreleased news) affects two observed assets, PCMCI may incorrectly report a direct causal link between them. The LPCMCI variant partially addresses this by allowing for latent confounders, at the cost of returning fewer oriented edges.

2. Stationarity assumption: The causal structure is assumed constant over the analysis window. In practice, crypto market dynamics change rapidly. Use rolling-window analysis to detect structural breaks.

3. Linear vs. nonlinear: With ParCorr, only linear causal effects are detected. A nonlinear causal mechanism (e.g., "BTC causes ETH to drop only when BTC drops more than 5%") would be invisible to the linear test.

4. Sampling frequency matters: Causal effects that occur faster than the sampling frequency appear as contemporaneous ($\tau = 0$) links in standard PCMCI, and their direction may be ambiguous. Use PCMCI+ and consider higher-frequency data.

5. Multiple testing: $N^2 \times \tau_{\max}$ dependent tests, corrected as described above — see deflated Sharpe and multiple testing for why the dependence structure rules out the obvious choice.

6. Sample size requirements: Reliable causal discovery requires sufficient data. As a rough guideline, aim for $T > 500$ for ParCorr with $N < 10$ variables, and $T > 2000$ for nonlinear tests or larger variable sets.

PCMCI vs. Other Methods

PCMCI's main advantage over Granger causality and transfer entropy is the MCI test, which correctly accounts for autocorrelation and common drivers in a multivariate setting. This is precisely the scenario encountered in crypto markets, where dozens of correlated, autocorrelated assets interact simultaneously.

Conclusion

The two-stage construction is the contribution worth taking away: sparse condition selection to bound the parent set, then an MCI test whose conditioning set includes the parents of the source as well as the target. That second half is what separates PCMCI from Granger and transfer entropy, and it is what makes the method viable on series as autocorrelated as crypto returns.

What the method does not come with is evidence. A recovered graph is a hypothesis about information flow, not a signal, and on this blog it does not count for anything until it has been run on real data with a stated basket, a stated date range, a dependence-aware correction, and at least one measurement of whether the graph holds still when the window moves. Until that run exists, treat everything here as tooling.

If the graph turns out not to be stable, that is the article. It would not be the first negative result here worth more than a working tutorial.

References

• Runge, J., Nowack, P., Kretschmer, M., Flaxman, S., and Sejdinovic, D. (2019). Detecting and quantifying causal associations in large nonlinear time series datasets. Science Advances, 5(11), eaau4996.
• Runge, J. (2020). Discovering contemporaneous and lagged causal relations in autocorrelated nonlinear time series datasets. Proceedings of the 36th Conference on Uncertainty in Artificial Intelligence (UAI), PMLR 124:1388-1397.
• Tigramite documentation: https://jakobrunge.github.io/tigramite/
• Tigramite GitHub repository: https://github.com/jakobrunge/tigramite