andrey.core.ci_kernel.RCoT

class andrey.core.ci_kernel.RCoT(data, *, n_features=5, n_features_cond=100, seed=None)

Randomized conditional correlation test (RCoT), a fast approximation of KCI.

Each query maps x, y, and the conditioning set through random Fourier features of a Gaussian kernel, removes from the x and y features what a linear ridge regression on the conditioning-set features explains, and tests the cross-covariance of what is left. The cost of a query grows linearly with the number of rows.

Parameters:

• data (ndarray) – Observed data, one row per sample. Each column is standardized to unit variance.

• n_features (int) – Random Fourier features sampled for x and for y.

• n_features_cond (int) – Random Fourier features sampled for the conditioning set. More features remove more of the conditioning set’s effect, at a cost that grows with their square. The default is four times the 25 of [1]: with 25, two conditioning variables that drive both x and y leave enough of their effect to reject most true independences at 500 rows.

• seed (int | None) – Seed for the random features. None uses ANDREY_SEED, or 0 when it is unset.

Raises:

ValueError – If data is not 2-D, has fewer than 2 rows, contains NaN or infinity, or has values too large or too small to standardize, or a feature count is below 1.

Warns:

SampleSizeWarning – On the first query with a conditioning set, when data has fewer than 500 rows, the smallest sample size [1] evaluates.

Notes

How it works. A variable v with median pairwise distance s (over its first 500 rows) gets features sqrt(2) cos(v W + b), with frequencies W from the kernel’s spectral distribution, normal of standard deviation 1 / s, and phases b uniform on [0, 2 pi]: the random Fourier features [5] of the Gaussian kernel exp(-|v - v'|^2 / (2 s^2)). The median distance is that kernel’s bandwidth, the median heuristic [4] that [1] uses; [1] prints the kernel as exp(-|v - v'|^2 / s). The frequencies and phases are sampled by randomized quasi-Monte Carlo, the quasi-Monte Carlo feature maps of [6]: the points of a scrambled Halton sequence, one coordinate per variable and one for the phase, give the frequencies through the normal quantile function and the phases as 2 pi times the last coordinate. Each point is uniform on the unit cube, so the features stay unbiased for the kernel, and the points cover the cube more evenly than independent samples, so a few features approximate the kernel more closely [6]. Every feature is standardized [1]. With A and B the features of x and y, and C those of the conditioning set, A and B are replaced by their residuals on C (ridge 1e-10). The statistic is n times the squared Frobenius norm of the residuals’ cross-covariance matrix, as [1] defines it. Under the null it is a weighted sum of chi-squared variables whose weights are the eigenvalues of the covariance of the products of residuals; the p-value matches the first three cumulants of that sum with a shifted and scaled chi-squared law (Hall-Buckley-Eagleson) [2] [3]. An empty conditioning set uses A and B, centered.

The features of x and of y are sampled for each conditioning set, from generators seeded by seed, the variable, and the set, and those of the conditioning set from one seeded by seed and the set; each column enters by its tie_ranks() rank, computed from the column’s own data. The queries of one pair with different conditioning sets therefore use independent features, as in [1], where every test samples its own. PC removes an edge when any conditioning set gives a p-value above alpha, and features shared by those queries would make them reject together, which keeps edges between conditionally independent variables. Queries on the same variable and conditioning set share its features, and the test object keeps them, and the residuals, for later queries. A p-value does not depend on the order of the queries, of the arguments, or of the columns, and the same seed reproduces it.

Limitations. Measured in PC on simulated data from sparse random graphs (mean degree 2), each variable a weighted sum of z + sin(z), or of cos(z), over its parents plus Gaussian noise, with 20 datasets per size on 4 cores. From 5 to 50 variables with 500 to 5,000 rows, its median errors (SHD) stayed below the empty graph’s at every size, such as 9.5 against 50 at 50 variables and 5,000 rows. The median fit took 0.2 to 0.5 seconds at 5 variables, 0.4 to 1.2 seconds at 10, 0.8 to 3.4 seconds at 20, and 3 to 15 seconds at 50, and the slowest took 70 seconds. Measured on columns independent of each other, its rate of rejecting a true conditional independence with two conditioning variables was about twice alpha at 100 rows and near alpha from 500 rows. With three, four, or six conditioning variables the rate averaged about 1.3 times alpha from 500 to 2,000 rows, with no row count from which it held at alpha. The p-value is random through the features: a different seed moves it, most near the threshold. The regression on the conditioning-set features removes a strong nonlinear effect of that set only approximately, and what remains reads as dependence: with two or more conditioning variables that drive both x and y, the test rejects a true independence more often than alpha, increasingly so as the rows grow (on simulated data with three such variables and 2,000 rows, nearly always). RCoT tests conditional uncorrelatedness of n_features features per variable, which can miss a dependence the full kernel test finds, such as one between an odd and an even function of the same cause on a few hundred rows. Because the features are sampled for each conditioning set, PC tests a weak dependence against every conditioning set with a new chance to miss it. On the datasets above, PC with RCoT missed 14 to 20% of the true edges at 500 rows and 1.5 to 4% at 5,000 rows, and 2 to 6% of the edges it kept were not true edges.

References

[1] (1,2,3,4,5,6,7)

Strobl, E. V., Zhang, K., and Visweswaran, S. (2019). Approximate kernel-based conditional independence tests for fast non-parametric causal discovery. Journal of Causal Inference 7(1), 20180017.

[2]

Hall, P. (1983). Chi squared approximations to the distribution of a sum of independent random variables. Annals of Probability 11(4), 1028-1036.

[3]

Buckley, M. J., and Eagleson, G. K. (1988). An approximation to the distribution of quadratic forms in normal random variables. Australian Journal of Statistics 30A, 150-159.

[4]

Gretton, A., Fukumizu, K., Teo, C. H., Song, L., Scholkopf, B., and Smola, A. J. (2007). A kernel statistical test of independence. Advances in Neural Information Processing Systems 20, 585-592.

[5]

Rahimi, A., and Recht, B. (2007). Random features for large-scale kernel machines. Advances in Neural Information Processing Systems 20, 1177-1184.

[6] (1,2)

Yang, J., Sindhwani, V., Avron, H., and Mahoney, M. W. (2014). Quasi-Monte Carlo feature maps for shift-invariant kernels. Proceedings of the 31st International Conference on Machine Learning, PMLR 32(1), 485-493.