Skip to content

Commit ae4bee8

Browse files
hmgaudeckerclaude
andcommitted
Mahler-Yum example: split the PRNG key once; drop an unsupported causal claim
Round-3 Pro audit findings, both verified at source. F2. `create_inputs` drew `type_indices` from `key` and then split that same already-consumed key for the health, adjustment-cost and productivity-shock draws. JAX's contract is split-before-consume, and its own key-reuse checker rejects this exact pattern: the type draw and the other three are not guaranteed independent. The root key is now split once into four children, each consumed exactly once. No bias was measured at this calibration (a million-draw correlation test found none), but independence is what the initial distribution assumes, so it should hold by construction rather than by luck. This changes every initial draw, so `tests/test_mahler_yum_2024.py`'s expected labor-supply distribution and mean-wealth profile at seed=32/n=10000 are now stale and need regenerating on a GPU. Not done here. F5. The module docstring claimed the Fortran's old-habit continuation was a shortcut FORCED by the continuous effort search landing off the habit grid, since the code could not interpolate along that axis. That is refuted by the authors' own simulator: `mod_simul.f90` carries `effhabit` as `real(8)`, sets the next habit to the continuous golden-section effort (line 634), and calls `bilini(agrid, fgrid, ...)` at the continuous incoming habit throughout. The two-dimensional interpolation the solver would have needed already existed a few files away. The mechanism (the effort objective interpolates over assets only) stands; the motive does not follow from the source and is no longer asserted. Also softened the numerical-scope note: the moments' numerical derivative not converging is a property of this implementation, but the measurements do not isolate WHICH numerical choice causes it. The grid is a consistent explanation, not an identified one. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01PVSLmeMo6iWUTcNER7bvM6
1 parent fe483fd commit ae4bee8

1 file changed

Lines changed: 29 additions & 10 deletions

File tree

src/lcm_examples/mahler_yum_2024/__init__.py

Lines changed: 29 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -19,6 +19,16 @@
1919
*simulation* meanwhile does advance the habit (`mod_simul.f90:634`), so its
2020
solve and simulate are mutually inconsistent.
2121
22+
What their solver omits is an interpolation along the habit axis: with
23+
`puregrid = 2` the chosen effort lands off the 30-point habit grid (`nf = 30`),
24+
so evaluating `EV` there needs one, and inside the effort objective they
25+
interpolate over assets only (`lini(agrid, EV(:,h), ap)`). That is the mechanism;
26+
the source does not reveal the motive, and it does not support calling the
27+
omission *forced*. Their own simulator carries `effhabit` as `real(8)`, sets the
28+
next habit to the continuous golden-section effort (`mod_simul.f90:634`), and
29+
calls `bilini(agrid, fgrid, ...)` at that continuous habit throughout — the
30+
two-dimensional interpolation the solver would have needed already existed.
31+
2232
`next_lagged_effort` below implements eq. (6). Consequently this model must NOT
2333
be expected to reproduce the paper's published model column at the published
2434
Table-II estimates: those estimates were fit against the other recursion.
@@ -44,10 +54,12 @@
4454
5-node equal-probability rule integrates it exactly only when a single branch
4555
dominates over the whole support; the node count needs a refinement check. Note
4656
that simulated adjustment shares are NOT restricted to multiples of 1/5 —
47-
simulation draws the cost continuously. (ii) The discrete grids make the
48-
simulated moments locally flat in the parameters, which is why MSM inference is
49-
withheld in the replication. That is a property of this solver, not of the
50-
paper's model: the authors' continuous searches do not share it.
57+
simulation draws the cost continuously. (ii) The simulated moments' numerical
58+
derivative does not converge across finite-difference steps, which is why MSM
59+
inference is withheld in the replication. The instability is a property of this
60+
implementation, not of the paper's model — the authors searched continuously and
61+
estimated successfully — but the measurements do not isolate *which* numerical
62+
choice produces it; the grid is a consistent explanation, not an identified one.
5163
"""
5264

5365
import dataclasses
@@ -933,18 +945,25 @@ def create_inputs(
933945
}
934946

935947
td = _build_type_distribution()
936-
key = jax.random.key(seed)
948+
# Split the root key ONCE and consume each child exactly once. Deriving
949+
# `type_indices` from `key` and then splitting that same `key` reuses a
950+
# consumed key, which JAX's key-reuse checker rejects: the type draw and the
951+
# three draws below are then not guaranteed independent. No bias was measured
952+
# at this calibration, but the guarantee is what the initial distribution
953+
# relies on, so it has to hold by construction rather than by luck.
954+
type_key, health_key, adjustment_key, shock_key = jax.random.split(
955+
jax.random.key(seed), num=4
956+
)
937957
type_indices = np.asarray(
938958
jax.random.choice(
939-
key,
959+
type_key,
940960
jnp.arange(len(td)),
941961
(n_simulation_subjects,),
942962
p=jnp.array(td["probability"].to_numpy()),
943963
)
944964
)
945965

946-
keys = jax.random.split(key=key, num=3)
947-
health_draw = jax.random.uniform(keys[0], (n_simulation_subjects,))
966+
health_draw = jax.random.uniform(health_key, (n_simulation_subjects,))
948967
health_thresholds = np.asarray(td["health_threshold"])[type_indices]
949968

950969
# Round the continuous initial habit to the NEAREST grid point. `searchsorted`
@@ -978,15 +997,15 @@ def create_inputs(
978997
"productivity_shock": np.asarray(
979998
shock_gridpoints[
980999
jax.random.choice(
981-
keys[2],
1000+
shock_key,
9821001
jnp.arange(5),
9831002
(n_simulation_subjects,),
9841003
p=stationary_shock_dist,
9851004
)
9861005
]
9871006
),
9881007
"adjustment_cost": np.asarray(
989-
jax.random.uniform(keys[1], (n_simulation_subjects,))
1008+
jax.random.uniform(adjustment_key, (n_simulation_subjects,))
9901009
),
9911010
}
9921011
)

0 commit comments

Comments
 (0)