因果发现：用 AI 从数据中学习因果图（PC、LiNGAM、NOTEARS）及其边界

import numpy as np
rng = np.random.default_rng(0)
n = 4000
X = rng.normal(size=n)
Y = 1.2 * X + rng.normal(size=n)
Z = 0.9 * Y + rng.normal(size=n)
def pcorr(x, y, z):
    rx = x - np.polyval(np.polyfit(z, x, 1), z)
    ry = y - np.polyval(np.polyfit(z, y, 1), z)
    return np.corrcoef(rx, ry)[0, 1]
print("corr(X,Z)      =", round(np.corrcoef(X, Z)[0,1], 3))
print("pcorr(X,Z | Y) =", round(pcorr(X, Z, Y), 3))

预期输出

corr(X,Z)      = 0.73
pcorr(X,Z | Y) = 0.00

import numpy as np
rng = np.random.default_rng(1)
n = 5000
X = rng.normal(size=n); Y = 0.8 * X + rng.normal(size=n)  # X -> Y
# Fit both directions; both explain the data symmetrically.
b_xy = np.polyfit(X, Y, 1)[0]
b_yx = np.polyfit(Y, X, 1)[0]
print("slope Y~X =", round(b_xy, 3), " slope X~Y =", round(b_yx, 3))
print("corr identical both ways:", round(np.corrcoef(X, Y)[0,1], 3))

预期输出

slope Y~X = 0.80  slope X~Y = 0.49
corr identical both ways: 0.62

import numpy as np
from scipy import stats
rng = np.random.default_rng(2)
n = 5000
X = rng.laplace(size=n)              # non-Gaussian cause
Y = 0.9 * X + rng.laplace(size=n)    # X -> Y
res_fwd = Y - np.polyval(np.polyfit(X, Y, 1), X)   # regress Y on X
res_bwd = X - np.polyval(np.polyfit(Y, X, 1), Y)   # regress X on Y
# Lower dependence (corr with regressor) indicates the causal direction.
print("forward |corr(res, X)| =", round(abs(np.corrcoef(res_fwd, X)[0,1]), 3))
print("backward|corr(res, Y)| =", round(abs(np.corrcoef(res_bwd, Y)[0,1]), 3))

预期输出

forward |corr(res, X)| = 0.00
backward|corr(res, Y)| = 0.21