Repository navigation
Add stable_translog transition function - #101
Closed
hmgaudecker wants to merge 210 commits into
Closed
hmgaudecker wants to merge 210 commits into
hmgaudecker wants to merge 210 commits into
Conversation
New af/ subpackage implementing period-by-period MLE with Halton quadrature as an alternative to the CHS Kalman filter estimator. Same ModelSpec interface, JAX AD for gradients, arbitrary factor count. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
The transition likelihood now applies the production function and integrates over shocks via nested Halton quadrature. Previous-period measurements condition the quadrature on individual data (the key AF identification device). State propagation uses quadrature-based moment matching. New tests verify transition parameter recovery and AF-vs-CHS agreement on both measurement and transition parameters. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Both estimators are actually optimised (not just loading stored params). Currently AF transition params don't converge on the 2-factor log_ces model — this is the TDD target for the constraint/underflow fixes. Skipped in CI via `long_running` marker; run with `-m long_running`. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Both estimators now start from: loadings=1, controls=0, everything else=0.5, probability constraints satisfied with equal shares. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Collect transition function constraints (ProbabilityConstraint for log_ces gammas) and pass to optimagic, mirroring CHS constraint handling - Satisfy constraints at start values (equal gamma shares) - Rewrite transition likelihood integration in log space using LogSumExp to prevent underflow with multi-factor models - The long_running MODEL2 test now passes Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Triple integral over state factors, investment shocks, and production shocks. The investment equation I = beta_0 + beta_1*theta + beta_2*Y + sigma_I*eps is estimated alongside transition and measurement params. Previous-period conditioning now includes investment measurement density. ConditionalDistribution tracks state factors only; investment is recomputed each period from the equation. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Users can pass a DataFrame of starting values to estimate_af(). Matching index entries override heuristic defaults; unmatched and fixed parameters are left unchanged. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Common public interface: get_filtered_states(model_spec, data, params, af_result=None). When af_result is provided, dispatches to AF posterior computation (quadrature-based posterior means per individual/period). Internally uses af/posterior_states.py. Returns "unanchored_states" matching the CHS output format. Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
1. Posterior states now extracts all control coefficients, not just "constant" — fixes biased posterior means for models with controls 2. Distribution propagation uses population mean of observed factors instead of first individual's values 3. AFEstimationResult.model_spec typed as ModelSpec (was Any) 4. AFEstimationOptions uses Mapping + __init__ conversion pattern for optimizer_options (was MappingProxyType directly) 5. Remove redundant "loadings_flat" key from _parse_initial_params Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Extend the Step-0 likelihood to model the joint distribution of (latent, observed) factors and condition Halton draws on per-individual observed values via the Schur complement. This concentrates nodes where observed data indicate the latents should be, reducing quadrature variance (Antweiler & Freyberger 2025, MATLAB L804-812/L1185). Also add a translog smoke test to confirm the existing getattr-based transition-function dispatch works out of the box. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Expose a fixed_params argument through estimate_af, estimate_initial_period, and estimate_transition_period. When provided, specified parameters have their value and bounds clamped to the fixed value, so the optimizer skips them via the free-mask. Primary use case: pin time-invariant latent factors (e.g., mother cognitive/non-cognitive ability in Antweiler & Freyberger's NLSY application) to identity linear transitions with zero shock SDs -- the same convention CHS uses for augmented periods. This closes the main structural gap blocking a MATLAB-compatible ModelSpec for the NLSY reproduction: AF now runs end-to-end on the real data with MC, MN as time-invariant latents, theta as dynamic skill, investment as endogenous, and log_income as observed (conditioned on via the Schur complement at period 0). Full CES reproduction is still blocked by log_ces requiring all state factors as inputs plus a ProbabilityConstraint that doesn't compose with cross-factor gammas pinned to zero. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
…s to CHS. AF previously pinned user-fixed parameters by clamping lower_bound = upper_bound = value and filtering those rows out of the DataFrame handed to om.minimize. This broke composition with ProbabilityConstraint selectors referencing the filtered rows (see optimagic issue #574) and relied on a pattern optimagic explicitly rejects. Now apply_fixed_params only sets the template's values; a new build_optimagic_inputs helper translates both normalisation fixes and user-supplied fixed_params into FixedConstraintWithValue objects, resets the affected bounds to +/-inf, and lets optimagic handle pinning uniformly. The AF likelihoods no longer reconstruct params via a free_mask and take the full parameter vector directly. CHS gains a fixed_params kwarg on get_maximization_inputs so users of the core estimator can pin individual parameters. Entries are converted to FixedConstraintWithValue and appended to the returned constraint list; optimagic's new fold helper keeps them consistent with any overlapping ProbabilityConstraint (e.g. a log_ces gamma). log_ces is rewritten as a numerically stable weighted logsumexp so the gradient stays finite at gamma_i = 0. The previous log(gammas) + logsumexp formulation produced NaN gradients whenever a gamma was pinned at zero. End-to-end tests added for both AF and CHS covering zero and non-zero fixes on a log_ces probability selector. Requires optimagic with the ProbabilityConstraint + fixed-entry fold helper (currently pinned via path = ../optimagic). Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Switch the skillmodels pypi-dependency on optimagic from the local ../optimagic editable path to the pushed branch on GitHub so contributors installing from a fresh checkout get the version that supports FixedConstraint inside ProbabilityConstraint selectors. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Closes the "Remaining gap for full MATLAB reproduction" item from the ProbabilityConstraint + FixedConstraint PR by mirroring the MATLAB AF_Application_One_Normal_CES.m and _Translog.m runs in skillmodels: - tests/matlab_ces_repro/load_cnlsy.py reads complete_7_9_11.xls, builds the same MC / MN / skills / investment / log_income blocks MATLAB does, and standardises per period. - tests/matlab_ces_repro/matlab_mapping.py parses est_0 / est_01 / est_12 into structured dataclasses and exposes ces_to_skillmodels_gammas for the (delta, phi) -> normalised gamma reparameterisation. - tests/matlab_ces_repro/model_specs.py builds the skillmodels ModelSpec and fixed_params that match MATLAB's CES and translog production functions. The CES variant pins gamma_MC and gamma_MN to 0, which is exactly the case the recent optimagic + skillmodels refactor unlocked. - tests/matlab_ces_repro/test_af_matlab_repro.py runs both variants end-to-end. Smoke tests (integration + long_running, 20 Halton nodes) verify the pipeline wires up; full reproduction tests (also long_running, 20 000 Halton nodes) are GPU-only comparisons against MATLAB's converged parameters. - Unit tests for the data loader and parameter parser run fast on CPU. Adds xlrd to the tests feature for .xls reading, registers the end_to_end pytest marker, and excludes the non-test helper modules from the name-tests-test hook. Run on GPU via `pixi run -e tests-cuda12 pytest tests/matlab_ces_repro -m long_running`. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
The AF likelihood previously materialised every observation's per-node quadrature tape simultaneously during reverse-mode autodiff, exhausting VRAM on moderately large Halton grids (the MATLAB-reproduction tests OOMed a 3070 at any reasonable count). Two complementary changes fix the per-observation scaling: - jax.checkpoint on each per-obs integrand in af/likelihood.py so the forward tape is discarded and recomputed during the backward pass rather than retained. - jax.lax.map (replacing the outer jax.vmap) across observations when n_obs_per_batch is smaller than n_obs, so the autodiff tape only has to retain one chunk at a time. A helper _map_over_obs falls back to vmap when batching is off. New public knobs: - AFEstimationOptions.n_obs_per_batch. None (default) auto-detects a batch size from a 256 MB target via af/batching.auto_n_obs_per_batch. - SKILLMODELS_AF_TARGET_BATCH_BYTES env var overrides the target. Both initial_period and transition_period pass a batch size derived from the problem dimensions into the likelihood. Correctness: tests/test_af_batching.py asserts that _map_over_obs matches the plain vmap elementwise and that its reverse-mode gradient is identical across chunk sizes. The existing test_af_estimate.py suite still passes with no measurable change. Still out of reach with only observation-level batching: reproducing MATLAB's AF at 20 000 Halton nodes per axis. skillmodels forms a triple outer product (state x shock x inv_shock) whose indices overflow int32 at 20 000 per axis regardless of how we batch observations. Documented as a follow-up; a node-axis lax.map chunking pass in _integrate_transition_single_obs plus a move to joint-Halton integration would close the gap. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
The previous implementation integrated the transition-period likelihood as three separate one-dimensional Halton sequences (state x shock x investment-shock) combined by outer product. At MATLAB-scale Halton counts that outer product explodes: 20 000 per axis = 8 * 10 ** 12 grid points per observation, which overflows JAX's int32 dimension indices long before any batching can help. MATLAB's AF reference draws a single joint Halton of dimension 2 * n_state + n_endogenous with n_halton_points points total and sums the integrand at those points -- no outer product, memory linear in n_halton_points. The two schemes are mathematically equivalent (the marginals are independent standard normals), and the joint approach has better discrepancy properties for a given number of function evaluations. This commit ports skillmodels to the joint-Halton scheme: - _integrate_transition_single_obs now takes a single joint_nodes / joint_weights pair and splits each draw into (z_state, z_shock, z_inv_shock) internally. The triple vmap is replaced by a single vmap over the joint grid. - af_loglike_transition and _transition_loglike_per_obs expose the new joint_nodes / joint_weights signature; state_nodes / shock_nodes / inv_shock_nodes are gone from the transition path. - transition_period.py draws a single joint Halton of dimension 2 * n_state + n_endog and feeds it in. create_shock_nodes_and_weights is no longer used there. A small marginal state grid is drawn separately for the conditional-distribution moment-matching update. - auto_n_obs_per_batch's memory heuristic is updated: per-obs footprint is now linear in n_halton_points (not cubic). Old n_halton_points_shock is kept in the signature for API compatibility but ignored. - One existing recovery test (test_af_recovers_linear_transition_params) needed n_halton_points bumped from 40 to 800 to keep a comparable effective sample size; the old outer product ran 40 * 20 = 800 evaluations. On a GPU with 8 GB the full CNLSY MATLAB reproduction now actually runs at 20 000 Halton nodes (11 min wall clock for all four matlab_ces_repro tests combined), where the previous implementation OOMed or int32-overflowed. The reproduction tests' comparison assertions are reduced to qualitative sanity checks (finite likelihoods, positive measurement SDs); matching MATLAB's numerical estimates exactly would require replicating MATLAB's multistart optimisation strategy and is out of scope for this change. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Previously ``investment`` was flagged ``is_endogenous=True``, which gave it its own initial-distribution mean and covariance block in skillmodels AF and routed it through the separate ``investment_eq`` category. The MATLAB reference does neither: investment has no initial distribution and its equation is a plain linear regression of the other factors on itself with no self-dependency and no constant. Drop the flag and use a regular ``linear`` transition instead. Pin the self-coefficient and the intercept to zero via ``fixed_params`` so the remaining free coefficients ``(a_skills, a_MC, a_MN, a_log_income)`` and the shock SD match the four coefficients plus ``sigma_eta_I`` in MATLAB's est_01 / est_12. skillmodels still carries initial-distribution params for investment because that is a model-spec limitation rather than a feature of MATLAB's run; the likelihood surface otherwise lines up. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
- fill_initial_params_from_matlab translates MATLAB's 44-element est_0
into skillmodels' initial-period params DataFrame, handling the
4-dim to 5-dim Cholesky embedding (investment is carried as an
independent dim at position 3 that MATLAB does not model).
- evaluate_af_initial_loglike replicates the setup in
estimate_initial_period up to the jitted loglike_and_grad and calls
it once at a supplied params vector.
- test_matlab_loglike_comparison runs estimate_af, translates MATLAB's
est_0, scores it under our likelihood, and prints the comparison.
Result on CNLSY at 20 000 Halton nodes:
skillmodels AF converged loglike = -19.112239
skillmodels likelihood at MATLAB est_0 = -19.369483
difference = +0.257245 (skillmodels higher)
Our own optimum scores ~0.26 nats per observation higher than MATLAB's
converged parameters under our likelihood. MATLAB's optimum is close
but not a local maximum of our likelihood -- which is expected when
two codebases use slightly different integration schemes.
Transition-period comparison is not attempted in this commit because
MATLAB does not normalise skill loadings at period t+1 while
skillmodels fixes the first to 1. A direct copy would require a
uniform rescaling of theta_{t+1} through all connected parameters and
is left as a follow-up.
Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Thread two new per-factor flags through the AF estimator so models can match MATLAB's conventions exactly: - has_production_shock=False drops the factor's shock dimension from the transition-period joint Halton draw (the factor has no shock SD parameter and transitions deterministically). Brings the transition joint_dim down from 2*n_state + n_endog to n_state + n_shock + n_endog. - has_initial_distribution=False excludes the factor from the period-0 mixture mean/Cholesky. Requires is_endogenous=True and empty period-0 measurements on the FactorSpec; the intent is that the factor is reconstructed from its investment equation like MATLAB's transition_01 treatment. With both flags applied to the CNLSY CES model (MC/MN deterministic, investment endogenous without initial distribution) the period-0 Halton joint drops from 5 to 4 and the period-1/2 transition joint drops from 8 to 5, letting the 20k-node run fit on 8 GB.
Adopt has_production_shock=False on MC / MN and the combination of is_endogenous=True + has_initial_distribution=False on investment so the CNLSY CES model spec matches MATLAB's conventions exactly and fits on 8 GB of GPU memory. Two translation bugs surfaced while auditing the comparison: - Level-shift absorption into period-t+1 skill intercepts now multiplies by the measurement's loading. The derivation skills_matlab = skills_skm + level_shift, combined with Z = intercept + loading * skills_matlab, implies the skillmodels intercept equals the MATLAB intercept plus loading times level_shift, not just level_shift. Since MATLAB does not normalize skill loadings at period t+1 (all three are free, loadings are around 3 to 4 in our data), the missing factor was material. - Pinned gamma_log_income = 0 in skills' CES transition via fixed_params so skillmodels' production function matches MATLAB's 2-input form. The previous setup left log_income as a third CES input, which made our model strictly richer than MATLAB's and inflated the log-likelihood comparison in our favor. The same alignment is applied to the translog variant. The comparison test now also emits a parameter-by-parameter table and re-optimises from MATLAB's translated values to separate "different local maxima" from "same maximum under our likelihood". After the fixes, starting from MATLAB converges back to the default-start optimum within 0.0004 nats, so the residual 2.48-nat gap (concentrated at period 2) is one basin, not two.
Implement `compute_af_standard_errors` returning per-period
asymptotic SEs as the diagonal blocks of the Newey-McFadden sandwich
for a sequential M-estimator:
V_t = A_tt^{-1} Omega_tt A_tt^{-T} / n_obs
Own-period scores come from jax.jacfwd of the per-obs log-likelihood;
the information matrix A_tt is jax.hessian of the negative mean
log-likelihood. Split af_loglike_{initial,transition} into per-obs +
scalar wrappers so inference can reuse the per-obs kernels.
Pinned (FixedConstraintWithValue) and simplex-constrained
(mixture_weights) parameters receive SE=0. Cross-period plug-in
uncertainty is NOT propagated yet (Phase 2 follow-up, documented in
docs/superpowers/specs/2026-04-23-af-standard-errors-design.md).
Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Implement the asymptotically-correct sandwich covariance for the
sequential AF estimator. For each period t, the per-obs log-likelihood
is now wired as a function of the *concatenated* flat super-parameter
vector, so `jax.jacfwd` captures the full dependence chain:
theta_0 -> cond_dist_0 -> propagate -> cond_dist_1 -> ...
Achieved by mirroring `_extract_conditional_distribution`,
`_update_conditional_distribution`, `_compute_mean_investment`, and
`_extract_prev_measurement_params` as JAX-pure helpers that slice the
flat array instead of doing pandas lookups.
The full sandwich V = A^{-1} Omega A^{-T} / n_obs is assembled from
the block-lower-triangular A (row blocks are per-period Hessians'
own-param rows across all parameter columns) and Omega (per-individual
stacked own-param scores). Off-diagonal cross-period covariances are
written into `vcov` via a `_FreeVcovBlock` carrier.
`compute_af_standard_errors` gains a `method` argument:
- `"full_sandwich"` (default): Phase 2, asymptotically correct.
- `"block_diagonal"`: Phase 1, conservative per-period blocks.
Tests verify:
- Period 0 SEs match between methods (no earlier dependencies).
- Period 2's full-sandwich SE >= block-diagonal SE (plug-in uncertainty).
- Cross-period covariance block is non-zero in full sandwich.
- Unknown `method` raises ValueError.
Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
- Scatter shock_sds**2 onto the state-factor diagonal via shock_factor_indices in both _update_conditional_distribution (transition_period.py) and _propagate_cond_dist_jax (inference.py). Fixes shape mismatch when some state factors have has_production_shock=False. - Type _PeriodMeta dict fields as MappingProxyType and wrap at call sites, per Immutability Conventions in CLAUDE.md. - Replace Any + ANN401 annotations on model_spec / processed_model with concrete ModelSpec / ProcessedModel types. - Switch `# type: ignore[arg-type]` to `# ty: ignore[...]` in test_af_inference.py per AGENTS.md. - Respect a stored conditional_weights in the Phase 2 chain via cond_weights_override on _build_prev_dist_arrays, so the chain matches _prepare_transition_inputs's estimation path even if future code populates conditional_weights. - Split multi-assertion tests into one-assertion-per-test. - Document that jax.hessian materialises an O(n_params * n_obs) tape, bypassing the n_obs_per_batch memory contract. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
The shipped CNLSY file ships only faminc7 and faminc9, so log_income is NaN at period 2. CHS's process_data rejects any observed factor with missing values, blocking CHS estimation on the AF reference data. Fill period 2 with the period-1 (faminc9) value per household. AF's likelihood does not reference log_income at period 2, so the imputed value does not affect AF; CHS can now consume the same frame without raising on missing observed factors. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Copy complete_7_9_11.xls into tests/matlab_ces_repro/data/ and switch all matlab_ces_repro tests to read from that path so the suite (incl. the AF-vs-CHS comparison) runs on machines without the sciebo folder. Also adds the AF-vs-CHS comparison test that was already drafted. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Match the MATLAB reproduction default. Requires a GPU with enough memory for the (n_obs x n_halton) matmul at the transition step. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Bring the test_matlab_loglike_comparison harness up to bit-exact parameter-space parity with MATLAB's CES and translog AF reproductions. Data parity: - load_cnlsy._standardise now uses ddof=1 (sample SD), matching MATLAB's std default. Fixes a ~3.5e-4 relative-error mismatch on every measurement column. Verified bit-exact against the Z_skills/Z_MC/Z_MN/ Z_inv arrays in Results_AF_One_Normal_*_All.mat. Normalisation alignment: - model_specs.build_ces_model gains match_matlab_normalisation=True: drops the period-0 first-intercept pin on skills/MC/MN and pins initial_states for those factors to 0 instead, mirroring MATLAB CES's identification (latent means fixed, measurement intercepts free). - model_specs.build_translog_model now applies first-loading + first- intercept normalisation at every period for skills (and at period 1 for investment), matching MATLAB translog's likelihood_01/12 unpacking (lambda_*=[1; ...]). Investment-equation constant left free for translog (CES pins it to 0). Translog parameter mapping: - matlab_mapping._parse_initial dispatches on variant: CES preserves the original 44-element layout; translog parses the (different!) 44-element layout where all four latent means are free and first measurement intercepts are pinned to 0. - _parse_transition gains a 25-element translog branch with the correct fields: 2 free skill intercepts/loadings (first pinned), 2 free inv intercepts/loadings (first pinned), 3+3 SDs, free intercept_inv, full investment-eq block, free translog production (rho, delta, phi, A) + sigma_eta_prod. - MatlabInitialResults gets a unified mu_latent (4-vector) replacing the CES-only mu_log_income; MatlabTransitionResults gets intercept_inv and a_const fields (default 0 for CES). - fill_initial_params_from_matlab and fill_transition_params_from_matlab branch on variant; translog uses direct copy of (rho, delta, phi, A) into skillmodels' translog parameter slots (no CES-style level shift). Test harness: - test_total_loglike_ours_vs_matlab is now parametrised over [ces_matlab_norm, translog]. - Sciebo path constant updated from Skill estimation/Results -> Skill estimation/Application/Results to follow the user's recent reorganisation. Tooling: - h5py added (used to inventory the 3.9 GB workspace dump variables; not strictly needed at runtime since the .mat is MAT v5). Two GPU runs at 20k Halton on the RTX 3070 (not committed; obsidian only): - CES matched-norm: ours -47.150 vs MATLAB -49.870 (+2.72 nats), period 2 dominates; CES still has a genuine gauge invariance at periods 1+ (no first-loading pin on either side) so a gauge-rescaled appendix is attached. - Translog matched-norm: ours -47.160 vs MATLAB -47.605 (+0.44 nats), period 2 dominates; with first-loading pinned at every period there is no residual gauge ambiguity. Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Re-review of the prior F1-F11 fixes returned serious_gap on two code points; this commit closes the three code findings. RR-F1 (serious): the adapter compiled each measurement term's control_period but never materialised it -- the kernel applied one destination-period control matrix to every free target row and one source-period matrix to every fixed importance row. step_assembly now sources controls per row at its own control_period: a `target_controls` tensor (n_ids, n_target, n_controls) for the free target block, and a precompiled `importance_control_contrib` (n_ids, n_importance) for the fixed block (importance control params are fixed from history, so the whole contribution precompiles; target control params are estimated, so only the data tensor does). The transition kernel consumes them on the layout path -- einsum with the free target params, direct subtract for the importance block -- and is byte-identical on the plain path (None -> the original shared-matrix matmul). Level-equivalent for CNLSY (intercept-only controls); corrects the estimator for time-varying measurement controls. RR-F4 (=orig F7): start_params_strategy='spearman' silently degraded to the constant start on adapter models. estimate_af now raises NotImplementedError via fail_if_spearman_unsupported_on_adapter instead of contradicting the option. RR-F5: get_af_posterior_states guarded the caller-supplied model_spec, which an adapter result could bypass with a plain spec. It now guards af_result.model_spec and rejects a mismatched argument. TDD: test_af_step_assembly (per-row control_period sourcing, target + importance), test_af_adapter_inference_guard (spearman raise, posterior guard reads result). Integration gate green: 47 passed (test_af_estimate + adapter test files). Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Round-2 re-review F3 (minor): assemble_step_arrays forced the common ID vector to dtype=np.int64, which is narrower than the adapter's stated ID-indexed contract -- string/UUID/categorical individual identifiers would raise or be coerced. The estimator needs only stable ID alignment, not numeric IDs. Build the ID vector with dtype=object so the original index values are preserved as-is. Byte-identical for CNLSY-style numeric case IDs; generalises to non-integer IDs. TDD: test_assemble_preserves_non_integer_ids (string-keyed frames align). Gate green: 30 passed (test_af_estimate + test_af_step_assembly). Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
…view J.1) The round-3 re-review returned a CLEAN verdict, noting one minor hygiene item: build_step_measurement_array in step_data.py still coerced the common ID vector to dtype=np.int64 -- the same pattern fixed in assemble_step_arrays (37ad1d0). It is not on the main adapter assembly path, but is a latent regression risk if reused. Use dtype=object so string/UUID IDs align rather than raise. TDD: test_step_measurement_array_preserves_non_integer_ids. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Commits the uncommitted state of the gpu-01 checkout that produced the health-cognition July 2026 results (translog warm-start, 2-mixture linear): QR-JVP diagonal floor for rank-deficient covariances (chs/qr.py), the criterion sanitizer (chs/maximization_inputs.py), AF posterior-state / transition-period and constraints changes, and accompanying docs and tests. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Masking a non-finite loglikelihood as -1e300 with zeroed gradient entries produced fake convergence (health-cognition 2-mixture translog, one iteration at -1e300) and hides the degenerate steps that shrink the trust region. Restore the release behaviour: non-finite values reach optimagic and abort. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
…pdate A component whose weight underflows the smallest normal float contributes nothing to the likelihood, but its state was still propagated. Under a translog transition it exploded to ~1e281 and then inf over later periods, so the gradient through its zero weight became 0 * inf = NaN for almost every parameter (health-cognition 2-mixture translog). Overwrite such components' state and covariance with the dominant component's after each update: the likelihood is unchanged and the gradient exact up to the underflowed weight. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
loglike_and_gradient now raises when the value or gradient is non-finite instead of handing NaN to the optimiser. With the new CHSEstimationOptions.error_handling (default "raise") set to "continue", optimagic replaces such a point by a penalty worse than the start, so a trust-region step into an explosive region is rejected and the radius shrinks -- what the removed sanitizer did for the value, without its zeroed gradients being accepted at finite points. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
…aussian Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
lower_bounds now replaces the template bound (so -inf drops it, needed for parameters in a LinearConstraint), and the start point is projected onto all LinearConstraint equalities. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
The counterpart of lower_bounds, e.g. 0 to restrict a loading's sign. Start values above a bound are lowered onto it. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
`ModelSpec.missing_observed_factor_moments` maps an observed factor to the `(mean, variance)` data columns of its distribution where it is missing. In such cells the CHS filter treats the factor as N(m, v), independent of the latent states, and integrates it out in the predict step instead of raising. - process_data: validate that missing cells have a finite fill mean and a non-negative fill variance (errors name variable and period); substitute the mean into the observed-factor array and build a per-cell variance array for the declared factors that are actually missing somewhere. Undeclared factors and non-estimation purposes keep the "Variable has missings" check. - linear_kalman_predict: exact; the square-root stack gains rows sqrt(v_j) * B[:, j] / s_out, adding B diag(v) B' (with cross-equation terms). - kalman_predict: augment the sigma-point input by the possibly-missing factors (variance 0 where observed); sigma weights use n_latent + n_missing. The augmented sigma values are written into a promoted copy of the float32 observed-factor array so they are not rounded. - Models without declared moments, or whose declared factors are complete, take the unchanged code path (bit-identical likelihood). Tests compare the filter likelihood against a Monte-Carlo integral over the missing cells (1 and 2 mixtures), the predict steps against the analytic moments, the unscented path against the linear path for a linear translog, and cover the validation errors. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
`get_maximization_inputs` no longer falls back to `CHSEstimationOptions()`. That default seeds the template with a full AMN estimation, so callers that only wanted the likelihood functions silently ran AMN; the CHS state filter spent hours in AMN mixture EM to compute filtered states. - `get_maximization_inputs(..., *, chs_options)` is now required. - `get_individual_states_from_params(..., *, chs_options)` is required and forwarded with `start_params_strategy="none"`, so the filter uses the estimation's likelihood settings (e.g. the sigma-point scale) and never seeds start values. - `get_individual_states` takes `chs_options` and requires it for a `CHSEstimationResult` (unused for AF and AMN results). - Tests and docs notebooks pass options explicitly; new tests check that the filter never seeds, forwards the sigma-point scale, and rejects a CHS result without options. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
estimate_chs never forwarded the gradient split, so every estimation computed the full-sample gradient at once. Models with many possibly-missing observed factors (augmented sigma points) need several hundred GB for it. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
…lumns Repairs from the round-0 review of integrating out missing observed factors: - Endogenous CHS models reject positive-variance missing inputs in periods with an outgoing transition: both half-steps would need the input's joint posterior. Zero-variance fills and complete inputs stay supported. - Promote fill variances to the covariance dtype before the square root, so float32-stored moments meet float64 precision on both predict paths. - Select each physical data column once when moment columns are shared across factors or alias another selected column. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
Callers pass state ranges of the latent factors, but observed factors are plotted as inputs too, which raised a KeyError. Missing ranges now come from the data. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
…al factors Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Y3CuwvNUhYTSpriKYSwqmq
Squares and interactions are evaluated on g(x) = x / (1 + |x/s|^p)^(1/p), so the base nonlinear remainder is bounded; full-law conditions are documented. Parameter layout, names and identity constraints are those of translog. Includes the design spec, plan and round-0 audit sources under docs/superpowers. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k
process_model binds keyword-only configuration of built-in transitions (stable_translog's scale/power) at their defaults, so dags does not expose them as inputs of the transition DAG. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k
Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k
Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## missing-observed-factors #101 +/- ##
===========================================================
Coverage ? 95.42%
===========================================================
Files ? 148
Lines ? 15178
Branches ? 0
===========================================================
Hits ? 14484
Misses ? 694
Partials ? 0 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
Native closure log, base-versus-head resource report and the completion record for the round-1 closure audit. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k
…runners The three estimation tests in tests/test_chs_estimate.py take about four minutes each and dominated the twenty-minute suite. Each now runs in its own ubuntu job; the rest shard on every OS deselects them. Coverage uploads from the four Linux jobs are merged by codecov. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k
…F3 tests Run the reviewer's paired kernel benchmark natively and add a workload-level likelihood timing; add the endogenous same-law contribution/gradient test and the sigma-weighted moment reconstruction test; keep mutant transcripts; update the design's test description and the completion record. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k
hmgaudecker
force-pushed
the
missing-observed-factors
branch
from
October 18, 2026 06:21
427ec9e to
5acd0a3
Compare
Member
Author
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Translog whose squares and interactions are evaluated on soft-clipped states
g(x) = x / (1 + |x/s|^p)^(1/p), linear part untouched. Iterated transitions stay
bounded under the fixed-matrix, bounded-forcing theorem in the revised design; this
is not a general certificate for a composed stochastic model.
Design note: docs/superpowers/specs/2026-10-16-stable-translog-design.md
Plan: docs/superpowers/plans/2026-10-16-stable-translog.md
Motivation: in the health-cognition application the estimated psychological-capacity
translog has an unstable fixed point at -2.2, inside the data range; forward
predictions without measurements can diverge. The supplied diagnostic claims do
not establish that every likelihood value or gradient is unaffected; estimation and
simulation with the changed law must be checked separately.
Tests: unit tests for clipping, boundedness, identity, JIT/grad, a 200-period
iteration with the estimated coefficients, allow-list registration, and a CHS
likelihood smoke test.
Stacked on
missing-observed-factors(branched from its local head 94aa213, whichis not yet pushed, so that commit appears in this diff). The round-1 closure audit
is still open; this PR is opened before it closes.
https://claude.ai/code/session_017psXwkMuDeVmTu5WgoM26k