A statistical concept

The “Markov blanket” began as a purely statistical idea. In his book on Bayesian networks, Judea Pearl defined the Markov blanket of a variable as a set of variables such that, once they are known, the rest of the network provides no further information about that variable (Pearl 1988). In a directed graph it consists of the variable’s parents, its children and its children’s other parents. In symbols, if $b$ is the Markov blanket of the internal state $\mu$ and $\eta$ is the remaining external state, then

$$p(\mu \mid b, \eta) = p(\mu \mid b).$$

Friston applied the idea to living systems: a cell, an organism, even a person can be described as internal states separated from their environment by a Markov blanket (Friston 2013). The blanket has two kinds of states: sensory states (the environment affects them and they affect the inside) and active states (the inside affects them and they affect the environment). Internal and external states are not directly connected and interact only through the blanket. Kirchhoff and colleagues discussed how such a boundary relates to biological autonomy (Kirchhoff et al. 2018).

A toy model that can be tested

The program below builds a linear-Gaussian loop: the external state $\eta$ produces a sensory state $s$; the internal state $\mu$ depends on $s$ and on its own past; the active state $a$ depends on $\mu$ and pushes back on $\eta$. If the blanket idea holds here, the correlation between $\mu$ and $\eta$ should vanish once the blanket is controlled for.

# A Markov blanket in a linear-Gaussian loop (a toy model). External state eta, sensory state s,
# internal state mu, active state a:
#   s_t = eta_t + e1;   mu_t = 0.9 mu_{t-1} + 0.3 s_t + e2;   a_t = mu_t + e3;
#   eta_{t+1} = 0.95 eta_t - 0.5 a_t + e4   (action pushes back on the world)
# Test: is mu_t independent of eta_t once the blanket (s_t and mu's own past) is known?
import random, math
rng = random.Random(7)
T = 60000
eta, mu = 0.0, 0.0
rows = []
for t in range(T):
    s = eta + rng.gauss(0, 0.5)
    mu_prev = mu
    mu = 0.9 * mu + 0.3 * s + rng.gauss(0, 0.5)
    a = mu + rng.gauss(0, 0.5)
    eta_next = 0.95 * eta - 0.5 * a + rng.gauss(0, 0.5)
    rows.append((eta, s, mu_prev, mu, a, eta_next))
    eta = eta_next
rows = rows[1000:]
col = {k: [r[i] for r in rows] for i, k in enumerate(("eta", "s", "mu_prev", "mu", "a", "eta_next"))}

def residual(y, xs):
    """Ordinary least squares of y on xs (with intercept); returns the residuals."""
    X = [[1.0] + [x[i] for x in xs] for i in range(len(y))]
    p = len(X[0])
    XtX = [[sum(r[i] * r[j] for r in X) for j in range(p)] for i in range(p)]
    Xty = [sum(r[i] * yy for r, yy in zip(X, y)) for i in range(p)]
    for i in range(p):                                   # Gauss-Jordan elimination
        piv = XtX[i][i]
        XtX[i] = [v / piv for v in XtX[i]]; Xty[i] /= piv
        for j in range(p):
            if j != i:
                f = XtX[j][i]
                XtX[j] = [vj - f * vi for vj, vi in zip(XtX[j], XtX[i])]; Xty[j] -= f * Xty[i]
    return [yy - sum(b * v for b, v in zip(Xty, r)) for r, yy in zip(X, y)]

def corr(x, y):
    mx, my = sum(x) / len(x), sum(y) / len(y)
    sxy = sum((a - mx) * (b - my) for a, b in zip(x, y))
    return sxy / math.sqrt(sum((a - mx) ** 2 for a in x) * sum((b - my) ** 2 for b in y))

def pcorr(x, y, given):
    g = [col[k] for k in given]
    return corr(residual(col[x], g), residual(col[y], g)) if g else corr(col[x], col[y])

print(f"samples: {len(rows)}")
for x, y, given in (("mu", "eta", []), ("mu", "eta", ["s"]), ("mu", "eta", ["s", "mu_prev"]),
                    ("mu", "eta_next", []), ("mu", "eta_next", ["a", "eta"])):
    label = f"corr({x}, {y}" + (f" | {', '.join(given)})" if given else ")")
    print(f"{label:34s} = {pcorr(x, y, given):+.3f}")
# samples: 59000
# corr(mu, eta)                      = +0.228
# corr(mu, eta | s)                  = -0.045
# corr(mu, eta | s, mu_prev)         = -0.001
# corr(mu, eta_next)                 = -0.174
# corr(mu, eta_next | a, eta)        = -0.002

The internal and external states are correlated at 0.228: the inside “knows” something about the outside. Controlling for the current sensory state alone is not enough (partial correlation $-0.045$), because the internal state’s past carries information about the external state’s past; once the internal state’s own past is counted as part of the blanket, the partial correlation falls to $-0.001$. The same holds on the other side: the internal state and the next external state are correlated at $-0.174$ (action pushes back on the world), and almost zero once action and the current external state are controlled for. Two things follow: in this model the blanket does separate inside from outside, and in a dynamic system the “blanket” must include time; it cannot be read off the sensory input of a single instant.

Changing the input channel

If the internal states can learn about the world only through the sensory states, then changing the source of the sensory states changes the internal states. The program below takes the same loop and replaces part of the sensory input with a low-quality source: biased (mean $-1$), noisy, and unrelated to the actual world.

# Same loop as before, but part of the sensory channel now comes from a low-quality source n
# (biased and noisy: mean -1, sd 1.5) instead of the world: s_t = (1 - w) eta_t + w n_t + e1.
# The internal state's long-run distribution depends on what crosses the blanket (a toy model).
import random, math, statistics as st

def run(w, T=40000, seed=11):
    rng = random.Random(seed)
    eta = mu = 0.0
    mus, etas = [], []
    for t in range(T):
        n = rng.gauss(-1.0, 1.5)
        s = (1 - w) * eta + w * n + rng.gauss(0, 0.5)
        mu = 0.9 * mu + 0.3 * s + rng.gauss(0, 0.5)
        a = mu + rng.gauss(0, 0.5)
        if t > 1000:
            mus.append(mu); etas.append(eta)
        eta = 0.95 * eta - 0.5 * a + rng.gauss(0, 0.5)
    return st.fmean(mus), st.pstdev(mus), st.correlation(mus, etas), st.fmean(etas)

print("share of input from the low-quality source | internal mean  sd  corr(internal, world) | world mean")
for w in (0.0, 0.2, 0.4, 0.6, 0.8):
    m, sd, r, me = run(w)
    print(f"{w:42.0%} | {m:+13.2f} {sd:4.2f} {r:+21.2f} | {me:+10.2f}")
# share of input from the low-quality source | internal mean  sd  corr(internal, world) | world mean
#                                         0% |         +0.01 1.31                 +0.22 |      +0.00
#                                        20% |         -0.01 1.27                 +0.18 |      +0.24
#                                        40% |         -0.05 1.24                 +0.13 |      +0.63
#                                        60% |         -0.13 1.24                 +0.04 |      +1.39
#                                        80% |         -0.33 1.27                 -0.09 |      +3.43

The result has two sides, and the second matters more. First, as the share of low-quality input rises, the correlation between the internal state and the real world falls from 0.22 to below zero: the internal model gradually decouples from the world. Second, the mean of the internal state moves only a little (from 0.01 to $-0.33$), because the system uses action to keep its internal state stable; the price is that the world itself is pushed away, the mean external state drifting from 0 to $+3.43$. In a closed loop, biased input does not only produce a wrong view; through the actions taken on that view it also changes the environment systematically.

This is a simplified analogy, but it points to a practical question: how much of the information a person receives each day comes from sources connected with the real environment they live in? The audit in Chapter 8 makes this its third step.

Two uses of the boundary

Taken with the previous chapter, this one shows two kinds of entry point for changing a self-sustaining loop. One lies outside the blanket: changing the environment, for example moving to surroundings that no longer reinforce the old pattern. The other lies in the blanket itself: changing what comes in (the sensory states) or what goes out (the active states). In Chapter 2’s model, “keep opening in the sandbox” is an intervention on the active states; in this chapter’s model, “change the information source” is an intervention on the sensory states. Neither rewrites the internal belief directly; both change the evidence from which the belief is updated.