Skip to content

chc.dynamics_id

Stability tier: stable.

Causal identification of the residual's control channel, not just of a scalar effect.

chc.train fits a residual by prediction error. Under a confounded logging policy that is the wrong object: if the historical action was chosen from a covariate that also drives the state, then d(dx/dt)/du fitted on (x, u, x_next) is the observational response, and a controller that optimises against it moves the state the wrong way. That is the same failure the scalar headline benchmark shows -- and until now nothing protected the residual from it.

This module closes that gap with Robinson partialling-out lifted from a scalar effect to a state-dependent matrix. With nuisances g(x,z) = E[y | x,z] and m(x,z) = E[u | x,z], where y = (x_next - x)/dt - f_known(t, x, u) is the part the residual must explain::

y - g(x,z)  =  B_θ(x) (u - m(x,z))  +  eps,     E[eps | x, z, u - m] = 0

so the estimator regresses the residualised state rate on the residualised action, with K-fold cross-fitting of both nuisances. Neyman orthogonality is what buys the guarantee: the channel error is second order in nuisance error, which is exactly Results 18 (order transfer p -> 2p) and 19 (debias every channel via cross-fit Robinson DML). This module is those results' missing consumer, not a new one.

HONEST SCOPE, four limits worth stating before the code:

  • Only the channel B_θ is interventional. a_θ is fitted on the remainder and therefore absorbs whatever the omitted confounder contributes to the drift in-sample, so it is an observational-conditional drift. Planning is unbiased in the direction the optimiser moves; the predicted trajectory level still shifts if the confounder's distribution does. On real data this bites harder than it sounds: a confounder that trends with the state gets charged to positive feedback in x, so the fitted drift can come back unstable even where the channel is clean -- and an MPC uses the drift too. Check the drift's spectrum, and give measured exogenous drivers a place in the drift as regressors instead of leaving them only in adjust_for.
  • The residual must be control-affine (chc.residual.ControlAffineResidual). A general r_θ(x, u) has no partialling-out moment and gets no guarantee here.
  • With no adjustment set and no instrument nothing in the log identifies the channel, and an instrument identifies it only along the directions it moves the action. The estimator reports identified=False in both cases rather than returning a confident wrong answer -- that case belongs to chc.sensitivity, which prices the radius instead of pretending it away.
  • The premise is Markov: given x, z and u, a transition's noise does not depend on the past. Where the noise persists within a unit and a lever persists too, the state is a common effect of the past lever and the past noise, so E[eps | x, z, u - m] is not 0 and the channel is biased by an amount channel_error does not cover. persistence_check reads the noise's persistence off the moment's residual; adjusting for the state, the levers and z a period earlier removes the bias where the noise is AR(1) (ADR 0064).

Integrator = Literal['euler', 'rk4'] module-attribute

Which one-step map the fit is asked to be consistent with. See fit_causal_residual.

ConfoundedControlAffineSystem dataclass

Euler-logged dx/dt = f_known(x,u) + A x + B_true u + C z under u = K z + eta.

The vector-state analogue of chc.causal.ConfoundedLinearSystem and the DGP fit_causal_residual is scored against. z is the confounder: it drives the logged action through K and the state rate through C, so with z withheld the response of the rate to the action is B_true + C sigma_z^2 K^T / var(u) rather than B_true, where var(u) = K K^T sigma_z^2 + G G^T + sigma_eta^2 collects every source of action variance.

A scalar channel with B=1, C=2, K=-1.5, sigma_z=1, sigma_eta=0.5 and no instrument gives 1 - 2*1.5/2.5 = -0.2: the sign flip of the headline benchmark, reproduced per channel. Turning the instrument on adds G G^T to that denominator and dilutes the bias toward zero -- the exogenous shifter makes the log more informative, so the observational fit looks less wrong while being no more identified. Worth knowing before reading a number off it.

The Euler discretisation is deliberately the one fit_causal_residual inverts, so a failing recovery test means the identification is wrong, not that the integrator disagreed.

sample(n, key, known)

Draw n transitions: x (n,d_x), z, u, x_next, w (instrument).

CausalDynamicsFit dataclass

A fitted causal residual next to what is and is not known about it.

identified is the load-bearing field: False means nothing the fit was given pins the channel down from the log -- no adjustment set and no instrument, or an instrument whose moment falls short of rank (instrument_rank) -- and the residual is then the observational fit, kept only so the caller can compare.

error_cosines property

nu, the cosines of the periods the channel's error was read with where the fit was given periods, counted over the distinct periods, and the degrees of freedom of the t a test read off it is sized against; None where it was given none.

OmittedConfounderBound dataclass

How far a confounder the adjustment set left out could move a linear functional of the channel, and how strong it would have to be to move it to null. Experimental.

bias_scale is identified; strength is the analyst's, as Gamma is (docs/concepts/gamma.md), in partial R^2s rather than odds.

ChannelRegion dataclass

A region that holds the fitted channel's coefficients at a stated level, from one log.

Experimental: it may change or be withdrawn in any release. State s's coefficients c_s (the channel's row s, raveled lever by lever, each lever's features in turn) lie in {c : ||c_hat_s - c||_V <= radius[s]} for every state at once with probability at least level. V = D'D + ridge P is the matrix the fit's least squares solved, D the channel's design on the action residuals. channel_region says where it holds.

contains(channel)

Whether every state's row of channel, shaped as channel, lies in its ellipsoid.

operator_radius(state)

A bound, at the region's level, on ||B_hat(x) - B(x)||_2 at the state x: the largest move of the rate that the channel's error makes per unit action norm. It holds at every state at once, since each is read off the same ellipsoids.

rate_error(states, reach)

A bound, at the region's level, on ||(B_hat(x) - B(x)) a|| over the (N, n) states and every action deviation a with |a_j| <= reach[j]: the rate error that the plan's tube prices as one standard error (chc.decision's _model_error), read here at a level. The largest of a convex form over a box is at a corner, so it reads the 2^m corners.

stability_level(rho, state)

The largest level at which operator_radius at state is within rho: the probability the region gives that the channel's error there is at most rho, 0 where no level reaches it.

A certainty-equivalent gain whose stability holds within a radius rho of the channel (chc.regret.ce_explicit_constant_certificate) then stabilises the plant with at least this probability. That certificate's rho is the true plant's, and its channel is the discrete step's, dt times this one under Euler: pass its radius over dt. Read at the estimate, rho is the estimate's, which the region does not check.

Solves sum_s (R_s A + b_s)^2 = (rho / g)^2 for A, g the operator radius per unit of the radii, and A^2 = 2 log(n / delta) + log_det_ratio for delta (validation/tube_levels.mac).

PersistenceCheck dataclass

Whether a fit's transition noise persists within units, read off its channel moment's residual by persistence_check. Experimental.

The fit's premise is Markov: given the state, the levers and the adjustment set, a transition's noise does not depend on the past. Where the noise persists within a unit and a lever persists too, the state is a common effect of the past lever and the past noise, so given the state this period's lever and noise move together, and the channel is biased by an amount its error does not cover (ADR 0064).

ClosedLoopAttribution dataclass

Which coefficients of a control-affine fit the log identifies, on a tracked plant.

ClosedLoopAttributionCertificate dataclass

Two arms: an interaction the plant has, and one the tracking loop manufactures.

solve_channel_moment(state_residual, action_residual, states, *, instrument_action=None, degree=1, ridge=1e-06, weights=None)

Solve E[(y_res - B_θ(x) u_res) (x) (z (x) phi(x))] = 0 for the channel.

The Robinson moment on its own, separated from the nuisance estimation fit_causal_residual wraps around it. Public for two reasons: cross-fitted ridge-polynomial nuisances are a default rather than a commitment, and the point of a debiased score is that g and m may come from any learner -- gradient boosting, a neural net, a model the caller already had -- so the caller needs a way in that does not go through ours.

It solves the moment it is handed and checks no rank. Along a direction of the channel that instrument_action does not move, the moment holds at every channel, and the ridge picks one. Whether an instrument identifies the channel is the caller's to settle before the solve, which sees neither the raw shifter nor the covariates the nuisances took out, and a rank that decides it reads both.

Parameters:

Name Type Description Default
state_residual Array

y - g(x, z), shape (N, n).

required
action_residual Array

u - m(x, z), shape (N, m) -- the regressor.

required
states Array

x, shape (N, n) -- the channel's own feature argument, which is not residualised: B_θ is allowed to depend on the state.

required
instrument_action Array | None

the action-shaped variable z that enters the moment, when it differs from the regressor: the projection of the action on an exogenous shifter, for the case where the confounder is latent. None means the action residual instruments itself, which is the orthogonal (adjusted) case.

None
degree int

monomial degree of B_θ's dependence on the state.

1
ridge float

Tikhonov term on the moment's Gram matrix, on each coefficient scaled by its regressor column's mean square, weighted as the rows are, so that it reads the same in any units of the actions and the state. Added as a constant, it outweighed actions logged in units a millionth of their own, and set their channel near zero.

1e-06
weights Array | None

one weight per transition, shape (N,), on its moment; None weighs all alike. A weight that is a function of the state alone keeps the moment orthogonal (see fit_causal_residual). Experimental: it may change or be withdrawn in any release.

None

Returns:

Type Description
Array

The channel coefficients, shape (n, m, n_features), consumable directly as

Array

attr:chc.residual.ControlAffineResidual.channel.

fit_causal_residual(known, data, dt, *, adjust_for=(), instrument=None, degree=1, channel_degree=None, nuisance_degree=2, folds=2, ridge=1e-06, seed=0, integrator='euler', drivers=(), weights=None, influence=False, clusters=None, periods=None)

Fit a ControlAffineResidual whose channel is the interventional control response.

Parameters:

Name Type Description Default
known Dynamics

the physics kept fixed; its own control dependence is known, not estimated.

required
data dict[str, Array]

columns x (N,n), u (N,m), x_next (N,n), plus any named in adjust_for and instrument.

required
adjust_for tuple[str, ...]

observed confounders. Empty and no instrument => identified=False.

()
instrument str | None

name of an exogenous action shifter, for when the confounder is latent. It enters as the moment's instrument, not as the regressor, so this is real 2SLS. Expect a variance premium, not a free lunch: identification here rides on however much of the action the shifter explains, and on the reference DGP that is 18%, which costs roughly an order of magnitude in channel error against adjusting for a logged confounder (0.10 vs 0.002 at N=4000). Still ~10x better than not identifying at all.

It identifies the channel only where its moment has rank: where the first stage's push on the actions moves with the actions along every direction of the channel the actions move (all but unmoved's), each less what the covariates' features predict, read without the ridge and in any units (CausalDynamicsFit.instrument_relevance), and the actions move one at least. Short of that rank, the moment holds at every channel along a direction it misses, and the fit reads as one with no adjustment and no instrument: identified=False, method="observational", the channel of the action residuals as their own instrument, kept to compare, and no error. It does so beside adjust_for too: naming an instrument says the covariates leave a confounder, so the fit does not fall back on them. Without the check an instrument of zeros read identified, the channel set by the moment's ridge and noise: on the reference plant [0.81, -0.06] with an error of 6.0, where the truth is [1.0, 0.5]. At nuisance_degree=0 the first stage is a constant, which no instrument enters, so no instrument has rank there.

The rank does not grade an instrument: a column of noise drawn apart from the action keeps it, and reads identified with a channel that misses. Its relevance grades it: on the reference plant such columns read 0.002 to 0.035, the medians of eight at 500, 2000 and 8000 rows, of the order of 1 / sqrt(N), where the plant's instrument read 0.39 to 0.48. No inference here is robust to a weak instrument.

None
degree int

the feature degree of the drift, and of the channel unless channel_degree says otherwise. 1 fits an affine drift and, by default, a channel affine in the state, which contains the constant channel §18/§19 cover without being restricted to it.

1
channel_degree int | None

the feature degree of the channel; None is degree. 0 fits the constant channel §18/§19 cover beside a drift of degree degree: n m coefficients where the affine channel has n m (n + 1), and as many fewer tangents in the rk4 fixed point. An affine channel is as accurate as a constant one only where the log is centred. On eight logs of chc.zones.ZoneMarketSystem, against one RK4 period's response (scripts/bench_channel_degree.py), the constant channel was 3.6% (linear matching) and 4.6% (harmonic) off in relative RMS wherever it was read; the affine one was as close at the logs' mean state, 16% off at the do-nothing point and 159% and 123% at x = 0, where chc.decision.Prescription.reach read it up to 0.12.

None
nuisance_degree int

flexibility of g and m. Richer nuisances are the whole point of cross-fitting -- orthogonality is what makes their error enter only at second order.

2
folds int

cross-fitting folds. 1 fits the nuisances on the same rows it residualises. That is not biased for the default nuisances and should not be sold as such: residualising y and u by the same linear projection whose span contains the truth is Frisch-Waugh-Lovell, so own-sample partialling-out is exactly unbiased and out-of-fold prediction only adds variance -- measurably worse at small N. Cross-fitting earns its keep against learners whose fit is adaptive to the sample (feature selection, trees, early stopping) or saturated enough to memorise it, where own-sample residuals collapse; folds>=2 is the safe default for that reason.

2
ridge float

the Tikhonov term of the nuisances, the channel's moment, the drift regression and the instrument's first stage. The nuisances read standardised covariates. The drift regression and the first stage scale each slope's term by its column's variance, beside a free bias; the moment, which has no bias, scales each coefficient's by its column's mean square, on the log's raw actions. So the ridge reads the same in any units of the state, the actions, the drivers and the instrument. Added as a constant, it outweighed a column logged in units a millionth of its own: such actions read a channel near zero, and such a state a drift slope near zero.

1e-06
integrator Integrator

the one-step map the fitted field is made consistent with.

"euler" (default, and what every release so far did) reads the rate off the log as (x_next - x)/dt, which is exactly right if the caller treats the result as a discrete-time model and steps it the same way. It is not right if the caller then integrates with chc.integrate.rk4_step, as chc.plan.causal_plan and chc.integrate.rollout do: the two disagree by the RK4 amplification factor 1 + z + z^2/2 + z^3/6 + z^4/24 at z = A dt. On a linear plant at theta*dt = 0.7 that is a decay fitted at -0.502 against a true -0.700 and, more to the point, a control channel fitted at 0.574 against a true 0.800.

"rk4" makes the field consistent with RK4 instead, without inverting RK4 in closed form, which does not exist for a nonlinear field. Coefficients read a rate off the log: their own fitted rate plus what their RK4 step leaves of x_next. The rk4 fit is the one the estimator returns when handed its own reading, found by Newton's method from the Euler fit in three to six steps between |A|dt = 0.05 and 1.2. Measured on the same plant, it recovers -0.699 and 0.800. Its standard errors are the fixed point's own (CausalDynamicsFit.channel_error).

Up to 0.6.0 the fixed point was iterated instead -- add the leftover to the target, refit -- and each state stopped once its defect's RMS stopped falling. That is not where the fixed point is: the RMS flattens under the noise first, and the iteration contracts the drift only by 1 - (1 + z + z^2/2 + z^3/6) a pass, 0.51 at z = -0.7. Even run until the RMS stopped falling at all, it stopped 0.9 standard errors of the channel and 1.4 of the decay short at theta*dt = 0.7 and N = 1000, and 0.09 of the decay short at 1.25 on a log solved exactly.

A log that falls by more in one step than RK4 can follow has no rk4 reading at all. RK4's growth factor bottoms out at 0.2704, at z = -1.596, and an exact linear mode falls below that past |A|dt = 1.308. The fit raises there rather than return a field that is neither reading. A channel that reads the state can run out sooner: on the plant in chc.mmm with its fastest channel decaying at theta*dt = 1.2, one seed in eight stalls in that channel's dependence on its own adstock, at a gap a damped Newton cannot close either, at both precisions. integrator_defect floors at the noise under both.

The default stays "euler" so no shipped fit changes meaning. chc.decision.prescribe defaults the other way, because it knows the field goes to the planner.

'euler'
drivers tuple[str, ...]

exogenous columns that push the state and that nothing the plan does can move -- weather, a demand forecast. data carries each at the start of the transition under its own name and at the end under f"{name}_next", and between the two it moves linearly, which is how chc.dynamics.DrivenDynamics replays a forecast and what the rk4 fit steps. Each enters the drift regression jointly with the state features, so the drift stops carrying the driver's correlation with the state (this module's first scope note): at its mean over the step under rk4, the push a linear path puts into a step, and at the start under euler, all the Euler map reads. Both ends join the channel's nuisance covariates. A logging policy that reads the forecast sets the action on where the driver is going, and with the start alone that is confounding: on a random-walk driver logged that way the euler channel came back 0.656 against the Euler map's own 0.780. Partialling out is safe for any driver, since a column nothing downstream of the action moves can be neither a mediator nor a collider of its effect. The fitted gain is CausalDynamicsFit.driver_gain; its standard error is inside drift_error.

()
weights Callable[[Array], Array] | None

a weight on each transition's channel moment, a function of its state alone. Experimental: it may change or be withdrawn in any release. A callable is called with x (N, n) and nothing else, since a weight that reads the action or the next state biases the channel at first order, while any weight on the state keeps the moment orthogonal. The weights are scaled to mean 1. The channel's standard error is robust, as for every fit, and here it has to be: a weight that loads a few rows loads the nuisances' error at them too. effective_sample_size says how few.

Under a channel class that contains the truth every weight estimates the same channel, and only the variance moves. The inverse of the rate's noise variance is the efficient weight in the limit, but it loads the states where the noise is small, which can be where the nuisances extrapolate: on a log whose rate noise ran exp(-x/2) against the action's exp(x/2), the channel scattered 0.87 as far as the unweighted fit's at N = 64000, its limit, and 1.15 as far at N = 4000.

Under a class that misses the truth, each weight estimates a different channel: the projection of the true one under w s^2 times the log's law of states, where s^2(x) is the variance the adjustment set leaves in the action at x. The weight then chooses where the fit is right. For a one-shot decision on a scalar action, made at states drawn from Q with regret curvature kappa(x) in the channel, w = kappa (dQ/dP) / s^2 makes the fit the best of its class for that decision. A weight that has to be estimated first -- s^2, or the inverse noise -- carries its own error into the channel at first order there, since the estimand moves with the weight. No weight on the state alone serves a dynamic plan, which reads the channel's slope along its path as well as its level: there the costate weight left 36 to 178 times the regret of the best fit in the class.

None
influence bool

whether to keep each transition's influence on the fitted parameters, CausalDynamicsFit.influence, which chc.misspecification.misspecification_cost reads to compare two fits of one log, and beside it the channel's representer, which omitted_confounder_bound reads with the moment residual. Off by default: under euler it costs a reverse pass per parameter rather than per channel coefficient, and it is N n p numbers the fit then carries, and N n q more.

False
clusters ndarray | Array | None

each transition's cluster, N labels of any one sortable kind, at least two distinct; one cluster's transitions are read across their periods by periods instead. The channel's covariance then sums each cluster's scores over its transitions and states before squaring them, CR1 with the factor G / (G - 1) (N - 1) / (N - k) over G clusters and k coefficients a state, so it holds whatever the dependence within a cluster; and a reader of influence sums its rows the same way. None takes every transition's every state as independent of the rest, which a log of units followed over time breaks wherever the model leaves out something that persists: the transitions of one unit then share it. The estimate does not move; only its error does. On panels of 5 to 80 units whose noise persists within each unit, a 5 % test of a zero channel rejects 12.5 % to 35.5 % of the time row by row and 3.75 % to 14.75 % by unit: CR1 still over-rejects with few clusters (ADR 0053). Two labels a transition, shape (N, 2), cluster it two ways, by its unit and by its period say. The error is then the largest of three (MacKinnon, Nielsen and Webb 2024): the two dimensions' sums less their intersection's (Cameron, Gelbach and Miller 2011), at CR1's factor for the smaller dimension's G and less any part that reads below nothing; and each dimension's sums alone, at its own. A shock every unit shares in a period, met by levers every unit moves together, makes one period's scores move together across units, which a sum within units leaves out (ADR 0061).

None
periods ArrayLike | None

each transition's period, a whole number of steps, on a log whose transitions fall in one group, one unit's or the units' of one declared cluster, which leaves no second cluster to sum within: one label a transition, two periods at least. The channel's covariance then sums the scores of each period's transitions over them and their states, and reads the sums' dependence across the periods by nu cosines of them, the equal-weighted cosine estimate with nu the most whole cosines at or below 0.4 T^(2/3) for T distinct periods (Lazarus, Lewis, Stock and Watson 2018), at the rows' factor N / (N - k), and a reader of influence sums its rows the same way. Summed within each period first, it holds whatever the dependence between one period's transitions, as between units that share a shock (Driscoll and Kraay 1998). A test read off it is sized against t(nu) (CausalDynamicsFit.error_cosines). None takes the transitions as independent, which a unit's persistent noise met by a lever that persists breaks: the scores, each the lever's residual times the noise, then persist too. On one unit of 200 transitions whose noise is AR(0.9) and whose lever is AR(0.7), a 95 % interval covered the channel on 0.703 of 1000 logs row by row and on 0.925 by cosines; where the noise is drawn afresh, the cosines' interval is 11 % wider. On three units of one cluster over 200 periods, which share a shock AR(0.9) and half their lever AR(0.7), it covered 0.692 row by row and 0.934 by the period sums (ADR 0075).

None

Returns:

Name Type Description
A CausalDynamicsFit

class:CausalDynamicsFit. Read identified before residual.

Raises:

Type Description
ValueError

on a negative channel_degree; on clusters that do not label every transition once or once in each of two dimensions, or name fewer than two clusters in a dimension; on periods that do not label every transition once with a whole number, or name one period alone, or that come with clusters; on an rk4 fixed point that does not converge; on a covariate named u or x_next, an instrument or a driver named x, u or x_next, a driver named twice, or an instrument that is also a covariate or a driver, at the step's start or its end. Each would read a column the fit reads in another role: adjusted for the action, the channel came back as the confounded regression's.

omitted_confounder_bound(fit, functional, *, cf_y, cf_d, rho=1.0, level=0.95, null=0.0)

Bound a linear functional of the channel against a confounder adjust_for left out. Experimental.

Chernozhukov, Cinelli, Newey, Sharma and Syrgkanis (2022) bound the bias of a linear functional of a partially linear fit by |rho| sqrt(cf_y cf_d / (1 - cf_d)) sigma nu: sigma^2 the rate's residual variance after the fit, nu^2 the second moment of the functional's Riesz representer, both identified, and two shares that are not. cf_y is the share of that residual variance the latent explains; cf_d the share of the long representer's second moment the fit's short one misses, which for one lever is the latent's partial R^2 with it. rho is the correlation between the two gaps, 1 at the worst. On a linear Gaussian plant with one latent the bound is attained, for any functional of any number of levers, and with two latents it is strict unless they move the lever and the rate in proportion (validation/omitted_confounder_bound.mac). Over several states the bound sums sigma_s nu_s, taking the shares as common to them; under Euler a functional of one state's channel reads that state alone.

functional weighs the channel's coefficients, in its shape (n, m, n_features). One lever's effect on one state is a single 1. Moving spend from lever a to lever b is +1 at b and -1 at a, and since its value is linear in the channel, its robustness_value is the confounding at which the move stops paying.

The sampling error follows DoubleML's: each confidence bound is one-sided at level, and carries the estimate's influence and bias_scale's, whose nu^2 part is nu^2 - alpha^2 for a representer alpha. Its sign convention for the estimate differs from this one in the cross term of the two influences, which vanishes in expectation. A fit given clusters has both influences summed within each cluster, as its channel's error is, and one given periods reads them across its periods by cosines, each period's transitions summed first. The bounds take a t's quantile with G - 1 degrees of freedom, G the clusters, two-way the smaller dimension's count, or the rows where there are none; or nu, the cosines', under periods (ADR 0075). At 5 units, a 5 % test of the channel read off its error rejected 9.75 % and 14.75 % of panels against the normal's quantile, and 3.25 % and 5.75 % against t(4)'s (ADR 0062). DoubleML takes the normal's.

Parameters:

Name Type Description Default
fit CausalDynamicsFit

an identified fit of fit_causal_residual, by adjustment, unweighted, made with influence=True.

required
functional Array

the weight on each channel coefficient.

required
cf_y float

the share of the residual variance the latent explains, in [0, 1).

required
cf_d float

the share of the long representer's second moment the short one misses, in [0, 1).

required
rho float

the correlation of the two gaps, in [-1, 1].

1.0
level float

each confidence bound's one-sided level, in [0.5, 1).

0.95
null float

the value the robustness values measure the distance to.

0.0

Raises:

Type Description
ValueError

on a fit whose channel is not identified by adjustment -- an instrument's representer is not the moment's -- or that is weighted, or kept no representer; on a functional of another shape; on shares, rho or level outside their ranges.

channel_region(fit, level=0.95, *, bound=None)

The self-normalised confidence ellipsoid of the fitted channel, at level.

Experimental: it may change or be withdrawn in any release. Abbasi-Yadkori, Pal and Szepesvari (NeurIPS 2011, Theorems 1 and 2): for y_t = <X_t, c> + eta_t with X_t predictable and eta_t conditionally R-sub-Gaussian, the ridge estimate c_hat satisfies ||c_hat - c||_V <= R sqrt(2 log(det(V)^(1/2) det(V0)^(-1/2) / delta)) + ||c||_V0 with probability at least 1 - delta, V = V0 + X'X. The regressors need only be predictable, so a log the plant was controlled in qualifies: a closed-loop policy's actions move with the states and the noise before them. Here V0 = ridge P, the fit's own ridge, one ellipsoid a state at delta / n.

It gives the plan's rate error a level (ChannelRegion.rate_error) and the certainty-equivalent certificate a probability (ChannelRegion.stability_level), beside the one-standard-error scales the tube carries, which it does not move.

Premises, and where they fail:

  • No hidden confounder. The noise must be mean zero given the actions; a confounder the fit does not adjust for pushes the rate with the action, and the ellipsoid misses the channel, at any level.
  • No persistent noise. eta must be a martingale difference; AR noise is not, unless its lags are in the fit (persistence_check).
  • R is read off the log: each state's root mean square of the moment's residual, over N - q rows, a plug-in for the theorem's known R, which also counts the cross-fitted nuisances' error. The action residuals use nuisances fitted on other rows, so the regressors are predictable to the nuisance's first order only.
  • The prior term ||c||_V0 is bound where given, else the estimate's own, a plug-in; at the default ridge it is about 1e-3 of the coefficients' effect on the rate.

Under a hidden confounder an instrument's sign-perturbed sums (Szentpeteri and Csaji, 2023) stay valid; they are not built here.

Parameters:

Name Type Description Default
fit CausalDynamicsFit

an identified Euler fit by adjustment, unweighted, whose actions move every direction of the channel.

required
level float

the joint level over the states, in (0, 1).

0.95
bound float | None

an upper bound, for every state, on the true coefficients' norm in the ridge's metric, sqrt(sum_l ridge P_l c_l^2); None reads the estimate's.

None

Raises:

Type Description
ValueError

on a fit that is not identified, read through an instrument, fitted under rk4, weighted, with a direction its actions never move, or with no ridge; a level outside (0, 1); or a negative bound.

persistence_check(fit, units, periods)

The lag-1 autocorrelation of fit's moment residual within units. Experimental.

units and periods label the fit's transitions, N each: the unit, of any one sortable kind, and the period the transition starts in, a whole number of steps. Two transitions of one unit whose periods differ by one are a pair. Each state's residual is divided by its root mean square, so each state counts alike in any units. The pairs' products, summed over the states, are summed within each unit before the test, since the transitions of one unit need not be independent; the statistic is their total over its CR1 spread, each unit's sum less its pairs' share of the total, against t(G - 1) (ADR 0062). Pairs of one unit alone, or of a fit given periods, whose units share one group's noise, are read across their periods instead: the products summed within each period, and the total over the spread of their projections on nu cosines of the pairs' periods, against t(nu), as the fit reads its error (ADR 0075). Pairs that all fall in one period leave no spread to read.

A rejection says the premise fails, not by how much the channel is off: that depends on how far a lever persists as well (PersistenceCheck). Adjusting the fit for the state, the levers and the adjustment set a period earlier removes the bias where the noise is AR(1): on panels of 200 units over 100 periods whose noise and lever are AR(0.7), a channel of 0.8 read 0.0130 low, 0.0007 low with them, and 0.0250 low with the state's lag alone (ADR 0064).

Raises:

Type Description
ValueError

the fit kept no moment residual, as where its channel is not identified; the labels are not one a transition; a period is not a whole number; or a unit has two transitions starting in one period.

closed_loop_gain_attribution(states, actions, rates)

Attribute a control-affine fit's coefficients to the controller and to the plant.

Result 41 left one item open: why the interaction coefficient b1 of dx/dt = d + a*x + (b0 + b1*x)*u comes out large and negative on a setpoint-tracked zone. The answer is a property of the log, not of the plant. A proportional loop u = gain*(setpoint - x) + u0 puts every sample on an affine manifold x = c + m*u with m = -1/gain, and restricted to that manifold the four-term class collapses to a quadratic in the action::

d + a*x + (b0 + b1*x)*u  =  (d + a*c) + (a*m + b0 + b1*c)*u + (b1*m)*u^2

Only b1 reaches the u^2 term, so b1 is the identified coefficient and the pole is not -- the inverse of the usual reading, and the reason Result 41 (a) found the reported pole to be a units artefact while lambda(u) = a + b1*u held. Matching against an observed response A + B*u + C*u^2 gives

``b1 = C/m = -gain*C``,

so the sign is decided by the curvature of the response in the action and the magnitude by the controller: a tighter tracker reports a bigger interaction from identical physics. proofs/closed_loop_attribution.v proves the matching identity and exhibits the explicit one-parameter family that leaves d, a and b0 free.

This also refutes the standing guess that b1 = -1/gain. That is the manifold slope, not the interaction; equating them forces C = 1/gain^2, a constraint on the plant rather than an identity (guess_is_a_constraint_not_an_identity).

exploration_budget is what buys the rest back: variation off the manifold restores the design's rank, and with it the separation of a from b0. Reading a pole off a fit whose budget is ~0 is reading the regulariser.

closed_loop_attribution_certificate(gains=(0.5, 1.0, 2.6, 8.0), exploration=(0.0, 0.05, 0.25, 1.0), samples=4000, seed=0)

Separate the interaction a plant HAS from the one a tracking loop MANUFACTURES.

Arm A: the plant really is d + a*x + (b0 + b1*x)*u. Under exact tracking the design is singular, yet b1 comes back exactly while the drift a does not -- the non-identification is real and it is the pole that suffers, not the interaction.

Arm B: the plant has no interaction, but its response carries a u^2 term the fitted class cannot represent. The fit answers with b1 = -gain*curvature: a coefficient that is entirely an artefact of misspecification amplified by the loop, growing linearly with the controller's gain. That is the mechanism behind Result 41's large negative b1, and why the number cannot be read as authority-falls-with-temperature without checking the budget.

The exploration sweep is the remedy: variation off the manifold restores the design's rank and the drift with it. All three claims are measured here, and the middle one is the one that would embarrass the entry if b1 turned out to be as unstable as a.