The permutation test is not a conditional independence test

[kslote1]

shuffle_test (discovery.py:903) builds its null by permuting the rows of X while leaving Z fixed:

X_perm = X[rng.permutation(len(X)), :]

This destroys the X–Z dependence as well as the X–Y dependence, so the null is "X independent of (Y, Z)", not H0: X, Y | Z. When the candidate parent is correlated with the conditioning set — the normal situation in an autocorrelated network time series — this inflates false positives.

A restricted/local permutation scheme (permuting within Z-neighborhoods, as in Runge's CMIknn) is the standard remedy. This is consistent with the false-positive rates I measured in the Gaussian path below:

Related: p_value = np.mean(null_cmi >= observed_cmi) can be exactly 0, and the package does report P_Value = 0.0 to users. The standard estimator (1 + #{null ≥ obs}) / (1 + n_shuffles) should be used.

[chaewon-research]

I reproduced the concern in this issue on current upstream main (ae9868f) and ran a controlled validation without modifying any library source files.

The clearest result is the effect of the current global row permutation on the conditioning structure. Using a synthetic conditional-null model

X = 3Z + eps_x
Y = -2Z + eps_y

with independent residuals, so X ⟂ Y | Z, I measured corr(X,Z) = +0.952000.

Across 500 surrogate draws:

• current global permutation: mean corr(X_perm,Z) = -0.000144
• local Z-neighbor permutation (k=5, Chebyshev distance): mean corr(X_local,Z) = +0.950885

In this setup, the current global permutation removes essentially all of the observed X-Z dependence, while the local permutation largely preserves it. This is consistent with the motivation for the local Z-neighbor permutation used in CMIknn-style conditional-independence testing.

I also ran a 2×2 calibration experiment separating the permutation strategy from the finite-permutation p-value correction. Each condition used N=200, 200 shuffles, 500 Monte Carlo replicates, and nominal alpha=0.05.

coupling	global + percentile	global + corrected p	local + percentile	local + corrected p
c=0	7.0%	6.6%	5.6%	5.2%
c=1	6.2%	5.6%	6.0%	6.0%
c=3	4.2%	3.8%	4.4%	4.0%

I don’t think these calibration results by themselves show a uniform Type-I-error improvement from the local permutation since the differences are fairly small and vary by coupling. The clearer result from this experiment is that the global permutation does not preserve the X-Z dependence structure, whereas the local surrogate largely does.

Separately, I tested the finite-permutation correction

p = (1 + #{null >= observed}) / (B + 1)

which slightly reduced rejection rates in these experiments and prevents reported p=0.

As a power sanity check, I also tested a strong conditional alternative:

X = 3Z + eps_x
Y = -3Z + 0.5X + eps_y

All four variants rejected 500/500 datasets, so I did not observe a power loss from the local shuffle in this particular strong-alternative setting.

I have the scripts and raw results preserved in an isolated clean environment, and I haven’t changed production code yet. If this direction looks right, I’d be happy to work on a focused implementation of the local Z-neighbor permutation and the finite-permutation p-value correction.