API reference

andrey.rcd#

andrey.rcd(data, *, max_explanatory_vars=2, cor_alpha=0.01, ind_alpha=0.01, shapiro_alpha=0.01)[source]#

Learn a linear non-Gaussian graph with hidden common causes by repetitive causal discovery.

RCD [1] repeatedly regresses variables on candidate ancestors and tests the residuals for independence, to find each variable’s ancestors, then its parents, and finally the pairs that share a hidden common cause.

Parameters:
  • data (ArrayLike) – Observations, one variable per column. A DataFrame’s column names become the structure’s labels.

  • max_explanatory_vars (int) – Most variables regressed out together when testing whether one variable is an ancestor of another, at least 1. Larger values find more ancestors and take longer.

  • cor_alpha (float) – Significance level of the Pearson correlation tests, in (0, 1).

  • ind_alpha (float) – Significance level of the HSIC independence tests, in (0, 1).

  • shapiro_alpha (float) – Significance level of the Shapiro-Wilk tests that decide whether a residual is non-Gaussian, in (0, 1).

Returns:

A GraphStructure of kind pag: i -> j for a direct cause and i <-> j for a pair with a hidden common cause. metadata["algorithm"] is "RCD"; metadata["confounded_pairs"] lists the <-> pairs as sorted [i, j] lists; weighted_adjacency holds the edge coefficients.

Return type:

StructureOutput

Raises:
  • TypeError – If max_explanatory_vars is not an int.

  • ValueError – If data has fewer than 6 rows, if max_explanatory_vars is less than 1, or if a significance level is not in (0, 1).

Notes

RCD assumes linear relations and non-Gaussian noise. It needs more rows than DirectLiNGAM: its tests compare residuals, and a weak hidden cause leaves only a small trace in them.

Limitations. Measured on linear data with uniform noise, one node in five hidden, and 2,000 rows, from 4 to 24 observed variables: its median directed-edge errors were below the empty graph’s at every size but 8 variables, where they equaled it. A fit took about 2 s at 4 variables and 18 s at 24, on two cores. Larger problems are untested; the work grows with the number of variables and with max_explanatory_vars.

References

Examples

>>> import numpy as np
>>> import andrey
>>> rng = np.random.default_rng(0)
>>> u = lambda: rng.uniform(-1, 1, 2000)
>>> h = u()  # hidden: a common cause of x0 and x1
>>> x0, x1 = 1.5 * h + u(), 1.5 * h + u()
>>> data = np.column_stack([x0, x1, x1 + u(), x0 + u()])
>>> out = andrey.rcd(data)
>>> out.metadata["confounded_pairs"]
[[0, 1]]