Applications · Alanine dipeptide kinetics

Agent-Guided Variance Reduction for Protein Kinetics

An agentic method search for importance-sampled alanine dipeptide kinetics

Michael M. Kim · Department of Mechanical Engineering, Stanford University

Read the full write-up (PDF) Watch the 2-minute explainer
Explore
Illustration: walkers rattling in two basins of a stylized landscape. Most excursions fall back; a few cross.

Overview

An agent chooses the estimator

Estimating protein conformational kinetics from importance-sampled trajectories is limited by weight degeneracy: when biased transition paths are reweighted back to the true dynamics, the variance of the path weights grows exponentially with path length and a handful of paths come to dominate the estimate. The remedy is resampling, but “resample” names a family of methods, not a single method, and which member works is system dependent.

This project asks which weight-control scheme tames that variance for a real peptide, and answers it as an agentic method search. Claude Code proposes an estimator variant, implements it in the sampler, launches a rate calculation, reads back the weight statistics and effective sample size against a brute-force reference, and iterates. On the β/C5 ↔ C7eq isomerization of alanine dipeptide, the loop compared sequential importance sampling (SMC), sequential importance resampling (SIR), and a branching resampled random walk (RRS). RRS is the only scheme that converges in both directions, recovering the brute-force rate to within 3–6% where the naive committor estimate is up to 6× wrong.

Contributions

  1. An exact yardstick

    A committor-based, reweighted rate estimator for alanine dipeptide with an exact brute-force reference in both directions.

  2. A controlled comparison

    SMC, SIR, and RRS weight control, compared against that reference.

  3. An agentic method search

    A demonstration that an agentic method search can select a variance-reduction scheme that a human would otherwise tune by trial and error.

Selected participant

Built with Claude: Life Sciences

  • July 7–13, 2026
  • Global · Virtual
  • Researcher Track

I was selected to take part in Built with Claude: Life Sciences, a global hackathon hosted by Anthropic in partnership with the Gladstone Institutes and produced by Cerebral Valley. Roughly 500 participants were chosen from nearly 6,000 applicants worldwide and given access to Claude Science and Claude Code, Anthropic API credits, and real research datasets, competing for a $100,000 prize pool. Over the week the cohort submitted more than 290 projects, pushing over 12,000 commits and 3 million lines of code. The work below is my submission to the Researcher Track.

Video 1 · Explainer

The project in two minutes

A two-minute walkthrough of the project: the weight-degeneracy problem in importance-sampled protein kinetics, the agentic search over weight-control schemes driven by Claude Code, and the branching resampled random walk (RRS) that recovers the brute-force β/C5 ↔ C7eq isomerization rate.

Explainer · 2 min

The problem

Rare events and the path weight

  • Hill relation
  • Unbiased for any I
  • Path weight W(Γ)
  • Effective sample size

The rare event problem

The rates of protein conformational changes control biological function, yet they are among the hardest quantities to compute. The transitions that matter (a loop reorganizing, a dihedral flipping, a domain rotating) are rare on the timescale of the fastest molecular motions, so a direct molecular dynamics trajectory spends essentially all of its steps inside a metastable region and almost none crossing between regions. Waiting for enough spontaneous crossings to estimate a rate can require billions of steps, and the cost grows exponentially with the barrier height.

The standard way out is to stop waiting. We decompose a long trajectory into short paths $\Gamma = \{\mathbf{x}_1, \ldots, \mathbf{x}_N\}$, each initiated when the trajectory leaves region A and terminated when it enters region A or B. Paths ending in A are failure paths; those ending in B are success paths. With $\mathcal{K}_\tau(\mathbf{y}\,|\,\mathbf{x})$ the stopped transition kernel at lag $\tau$, the density of a success path factorizes,

Path density $$ p(\Gamma_{\rm AB}) = \rho_{\rm A}(\mathbf{x}_1) \prod_{i=1}^{N-1} \mathcal{K}_\tau(\mathbf{x}_{i+1}\,|\,\mathbf{x}_i), $$

and the reactive probability is the integral over all success paths, $P_{\rm AB} = \int p(\Gamma)\, d\Gamma$. The rate then follows from the Hill relation (Hill, 2012; Baudel et al., 2023):

Hill relation $$ r_{\rm AB} = \Phi_{\rm A}\, P_{\rm AB}, \qquad \mathrm{MFPT}_{\rm AB} = 1 / r_{\rm AB}, $$

where $\Phi_{\rm A}$ is the flux of trajectories escaping region A (the inverse mean time between successive first exits). The flux is cheap, because escapes from a region are not rare. All of the difficulty is concentrated in $P_{\rm AB}$: the probability that an escaped trajectory commits to the far region rather than falling back.

Importance sampling and the path weight

For a genuine rare event the success probability is tiny compared with that of a failure path, so we modify the transition kernel with an importance function $I(\mathbf{x})$ (Cai et al., 2002; de Koning et al., 2005) that biases the walker toward B:

Biased kernel $$ \mathcal{K}'_\tau(\mathbf{y}\,|\,\mathbf{x}) = Z_\tau^{-1}(\mathbf{x})\, \frac{I(\mathbf{y})}{I(\mathbf{x})}\, \mathcal{K}_\tau(\mathbf{y}\,|\,\mathbf{x}), \qquad Z_\tau(\mathbf{x}) = \int_\Omega d\mathbf{y}\, \frac{I(\mathbf{y})}{I(\mathbf{x})}\, \mathcal{K}_\tau(\mathbf{y}\,|\,\mathbf{x}), $$

where the normalization factor $Z_\tau(\mathbf{x})$, computed on the fly by Monte Carlo integration, keeps $\mathcal{K}'_\tau$ a probability density. Propagating a walker under $\mathcal{K}'_\tau$ and reweighting back to the true dynamics attaches to each path a correction factor built from those same normalizations, and this defines the path weight:

Path weight $$ p'(\Gamma_{\rm AB}) = \frac{I(\mathbf{x}_N)}{I(\mathbf{x}_1)} \left[\prod_{i=1}^{N-1} Z_\tau^{-1}(\mathbf{x}_i)\right] p(\Gamma_{\rm AB}), \qquad W(\Gamma) := \prod_{i=1}^{N-1} Z_\tau(\mathbf{x}_i). $$

Setting $I(\mathbf{x}) = 0$ on A and $I(\mathbf{x}) = 1$ on B, the optimality condition $\int I(\mathbf{y})\, \mathcal{K}_\tau(\mathbf{y}\,|\,\mathbf{x})\, d\mathbf{y} = I(\mathbf{x})$ is exactly the backward equation whose solution is the committor. When $I$ equals this optimal importance function every $Z_\tau = 1$, and the bias enhances success paths by the known factor $I(\mathbf{x}_N)/I(\mathbf{x}_1)$. When $I$ is imperfect, the weights $W(\Gamma)$ restore the unbiased statistics exactly. The reweighted estimator of the reactive probability is

Reweighted estimator $$ \widehat{P}_{\rm AB} = \frac{1}{N_p} \sum_{i=1}^{N_p} I(\mathbf{x}_1^{(i)})\, W(\Gamma_i)\, \mathbf{1}\!\left[\Gamma_i \in \{\Gamma_{\rm AB}\}\right], $$

which is unbiased for any importance function, including an imperfect learned committor. This is the crucial licence the method operates under: errors in the guide degrade efficiency, not correctness. A neural network that is merely a decent approximation of the committor still yields a statistically exact rate, provided we can afford the variance it induces.

The weight-variance problem

Unbiasedness is not enough. Because $W(\Gamma)$ is a product of $N-1$ fluctuating factors, its variance grows exponentially in path length,

Variance growth $$ \mathrm{Var}[W] \sim e^{cN}, \qquad c > 0, $$

the classic weight-degeneracy pathology of sequential importance sampling. The efficiency of the estimator is governed by the effective sample size,

Effective sample size $$ \mathrm{ESS} = \frac{\left(\sum_i W_i\right)^2}{\sum_i W_i^2}, $$

which collapses from $N_p$ toward $1$ as a single weight comes to dominate. Proteins are close to the worst case for this pathology. The fast, stiff vibrational modes (bond stretches, angle bends) make the per-step factors $Z_\tau(\mathbf{x}_i)$ noisy at every step without advancing the slow $(\phi,\psi)$ reorganization we actually care about. Weight noise accumulates on the fast clock while progress accumulates on the slow one. Degeneracy therefore sets in quickly, and a raw importance-sampled estimator is unusable in practice.

Weight control

Three weight-control schemes

  • SMC · SIR · RRS
  • One estimator
  • One learned importance function

Resampling interrupts the exponential growth by periodically replacing the weighted population with a nearly equally weighted one of the same expectation. We compare three members of that family. All three share the same estimator and the same learned importance function, and differ only in how weights are controlled, so the comparison isolates exactly one design choice.

  • SMC Schematic

    Sequential importance samplingNo resampling

    Weights accumulate unchecked. This is the baseline, and the scheme most exposed to the exponential variance growth $\mathrm{Var}[W] \sim e^{cN}$.

    Trials completed in the search, of 1000

    β/C5 → C7eq1000
    C7eq → β/C531†

    Degenerates. Stalls after 31 of 1000 trials in the reverse direction.

  • SIR Schematic

    Sequential importance resamplingGlobal resampling

    The whole population is resampled by low-variance systematic resampling whenever $\mathrm{ESS} < \theta N_p$ (here $\theta = 0.5$), and every weight is reset to the mean. Global resampling is aggressive: it bounds the weights, but it can over-prune walker diversity.

    Trials completed in the search, of 1000

    β/C5 → C7eq118†
    C7eq → β/C5533†

    Unstable. Per-trial spread 2–3× the mean; both directions terminated early.

  • Selected
    RRS

    Resampled random walkBranching

    Each walker is treated individually. Whenever its weight leaves a target band $[w_{\min}, w_{\max}]$, it is branched into an integer number of unit-weight offspring,

    Offspring count $$ n_{\rm off} = \lfloor W \rfloor + \mathbf{1}\!\left[u < W - \lfloor W \rfloor\right], \qquad u \sim \mathcal{U}(0,1), $$

    so overweight walkers split and underweight walkers undergo Russian-roulette pruning, with survivors reset to unit weight (Shi, 2015; Vorba and Křivánek, 2016; Talbot, 2005). Branching keeps the weight distribution concentrated near unity while preserving the estimator’s expectation.

    Trials completed in the search, of 1000

    β/C5 → C7eq1000
    C7eq → β/C51000

    Converges. All 1000 trials in both directions.

† Run terminated before 1000 trials (degeneracy / non-convergence). The meters show trials completed; the full numbers are in Table 1.

The agentic loop

Claude Code in the loop

  • Claude Code
  • One routine: process_walkers
  • Logged end to end

Why an agent

Exploring this design space by hand is slow: each candidate is a small code change followed by a long run and a careful reading of variance diagnostics. So rather than tuning the estimator by hand, we treated the choice of estimator as a search problem and drove it agentically with Claude Code.

Why it fits

Conveniently, the SMC/SIR/RRS logic all lives in a single process_walkers routine, which means a new variant is a localized diff: exactly the kind of change an agent can make reliably and then evaluate honestly against ground truth.

Claude CodeClosed loop

Propose

Suggest a weight-control variant, informed by why the previous one failed: degeneracy, over-pruning, or non-convergence.

Implement

Write it as a branch of the sampler’s process_walkers routine: a localized, reviewable diff.

Run

Launch the biased rate calculation together with the matched brute-force reference.

Analyze

Read back the rate-versus-trial log, the path-weight histogram, and the ESS, and compare against the exact rate.

Figure 1. The agentic method search. Claude Code drives a closed loop over estimator designs: propose a weight-control scheme, implement it as a localized change to the sampler, run the biased rate calculation against a matched brute-force reference, and read back weight-variance and effective-sample-size diagnostics to decide the next variant. Claude also generated the diagnostics (weight distributions, convergence traces, MFPT bar charts) that appear as the figures below. The full propose–implement–run–analyze loop was logged end to end, so the search is auditable and reproducible. The loop selected RRS for this protein.

Setup

Alanine dipeptide, with an exact answer

  • 22 atoms
  • Ac-Ala-NHMe in vacuum
  • AMBER14 · OpenMM
  • 300 K
  • Seed 2026

The learned importance function

A SchNet-style graph neural network parameterizes the committor that guides the walkers.

The importance function $I(\mathbf{x})$ is a SchNet-style graph neural network (Schütt et al., 2017) acting on ten chemically distinct heavy atoms selected from the 22-atom topology, with 64 features, three interaction layers, a 20-Gaussian radial basis (cutoff 0.9 nm), and translational/rotational invariance built in. A single scalar head yields a bidirectional committor via $I_+ = \sigma(f)$ and $I_- = \sigma(-f)$, so the same network guides both the forward and reverse transitions. Each interaction layer updates atom features through a continuous-filter convolution,

Continuous-filter convolution $$ \mathbf{h}_i^{(\ell+1)} = \mathbf{h}_i^{(\ell)} + \sum_{j \neq i} \mathbf{W}^{(\ell)} \mathbf{h}_j^{(\ell)} \odot \boldsymbol{\phi}^{(\ell)}\!\left(\mathbf{e}(r_{ij})\right). $$

Training minimizes the log-space residual of the optimality condition, $\ln\!\left[\int I(\mathbf{y}) \mathcal{K}_\tau(\mathbf{y}\,|\,\mathbf{x})\, d\mathbf{y}\right] - \ln I(\mathbf{x})$, over swarm-propagated lags, using a bootstrapped target network updated by exponential moving average ($\tau_{\rm ema} = 5 \times 10^{-3}$), a replay buffer, and gradient clipping. The production importance function is an ensemble average of checkpoints from epochs 1000–1100, which reduces the estimator noise of any single network.

InputAtomic positions $\{\mathbf{x}_i\}_{i=1}^{22}$
Atom selection$\{\mathbf{x}_i\}_{i \in S}$,   $|S| = 10$
Pairwise distances$r_{ij} = \|\mathbf{x}_i - \mathbf{x}_j\|$
Atom embedding$\mathbf{h}_i^{(0)} = \mathbf{E}[s_i]$
RBF expansion$e_k(r_{ij}) = e^{-\gamma (r_{ij} - \nu_k)^2}$
Interaction layer × L
Linear map$\mathbf{W}^{(\ell)}\mathbf{h}_j^{(\ell)}$
⊙
Filter network$\boldsymbol{\phi}^{(\ell)}(\mathbf{e}(r_{ij}))$
Continuous-filter convolution, added as a residual (⊕)
Per-atom readout$f_{\rm RO}(\mathbf{h}) = \mathbf{w}_2^{\mathsf{T}}\, \mathrm{SiLU}(\mathbf{W}_1 \mathbf{h} + \mathbf{b}_1) + b_2$
Importance function$I(\{\mathbf{x}_i\}; \theta) = \sigma\!\left[\sum_i f_{\rm RO}\!\left(\mathbf{h}_i^{(L)}\right)\right]$
Figure 2. SchNet-style graph neural network parameterizing the committor $I(\{\mathbf{x}_i\}; \theta)$ for alanine dipeptide. Ten chemically distinct atoms are selected from the 22-atom topology; their pairwise distances feed a per-atom embedding and a radial-basis expansion. Each of the $L$ interaction layers updates atom features through a continuous-filter convolution. A per-atom readout, summed and passed through a sigmoid, yields the importance function.
Potential of mean force of alanine dipeptide in the backbone dihedrals phi and psi, with the beta/C5 and C7eq metastable regions marked by rings of sampled points Free energy in (φ, ψ)
Figure 3. The β/C5 ↔ C7eq benchmark. (a) Potential of mean force in the backbone dihedrals $(\phi, \psi)$; the two metastable regions β/C5 and C7eq (rings of white points) are separated by a free-energy barrier. (b,c) The two metastable conformations of alanine dipeptide (carbon black, nitrogen blue, oxygen red, hydrogen white).

The system

As a molecular proof of principle we take alanine dipeptide (Ac-Ala-NHMe, 22 atoms) in vacuum, the standard benchmark for conformational sampling and rare-event dynamics. Molecular-mechanics parameters come from the AMBER14 force field and the dynamics are propagated with OpenMM tooling. The system is evolved by constrained overdamped Langevin dynamics, integrated with the Euler–Maruyama scheme,

Euler–Maruyama step $$ \mathbf{x}(t + \Delta t) = \mathbf{x}(t) - \gamma^{-1}\mathbf{M}^{-1}\nabla U(\mathbf{x}(t))\, \Delta t + \sqrt{2 \gamma^{-1}\beta^{-1}\Delta t}\; \mathbf{M}^{-1/2}\boldsymbol{\xi}, \qquad \boldsymbol{\xi} \sim \mathcal{N}(\mathbf{0}, \mathbb{I}), $$

at $T = 300$ K, damping $\gamma = 1.0$ ps$^{-1}$, and timestep $\Delta t = 0.01$ fs, with hydrogen bond-length constraints imposed by SHAKE-like projection.

The regions

The conformation is characterized by the backbone dihedrals $(\phi, \psi)$, and the metastable regions are disks of radius $10^\circ$ in Ramachandran space. The $\alpha_{\rm R}$ region acts as an absorbing guard so that indirect crossings are not miscounted. All runs use random seed 2026.

  • Region Aβ/C5centered at (−150°, 160°)
  • Region BC7eqcentered at (−75°, 60°)
  • Absorbing guardαRat (−70°, −40°)

A benchmark with an exact answer

Alanine dipeptide places the two slow states in a two-dimensional dihedral plane whose free-energy surface we can draw explicitly. The β/C5 and C7eq regions sit in separate free-energy wells; because the peptide is small, unbiased MD still yields a converged first-passage rate in both directions. Every biased estimate can therefore be checked against an exact number rather than against itself, which is what makes the system usable as a testbed for a method search at all.

The estimator and the yardstick

How each rate is computed, and the exact numbers every estimate is checked against.

Rate estimator

The escape flux $\Phi$ is measured by counting successive first exits from a region during an equilibration run ($10^4$ escapes). The reactive probability $P_{\rm AB}$ is estimated by launching $N_p = 1000$ importance-guided walkers, each carrying a swarm of $n_{\rm swarm} = 100$ micro-steps used to estimate the normalization factor $Z_\tau$, accumulating weights, applying the chosen resampling scheme, and evaluating the estimator; the calculation is repeated for 1000 independent trials. RRS uses the weight band $[1, \infty)$ (Russian-roulette pruning of sub-unit weights) and SIR the threshold $\theta = 0.5$.

Brute-force reference

For an exact reference we launch $10^4$ unbiased walkers from each region and record first-passage times to the other, giving

  • β/C5 → C7eq$r_{\rm AB}^{\rm BF} = 1.08 \times 10^{-3}$ fs$^{-1}$MFPT $=$ 926 fs
  • C7eq → β/C5$r_{\rm BA}^{\rm BF} = 4.48 \times 10^{-4}$ fs$^{-1}$MFPT $=$ 2234 fs

These are the yardstick for every estimator below.

Results

SMC degenerates, SIR is unstable, RRS converges

Weight degeneracy is severe

The path-weight distributions make the central difficulty concrete. Although the vast majority of weights sit near the mean $\langle W(\Gamma) \rangle$, the distribution is extraordinarily heavy tailed, with rare paths reaching weights orders of magnitude larger. A single such path can outweigh the other $10^5$ combined, precisely the effective-sample-size collapse anticipated by the theory. Any scheme that does not actively control this tail produces a rate dominated by a few paths.

(a) β/C5 → C7eq Log-frequency histogram of path weights for the beta/C5 to C7eq direction, showing a heavy tail far above the mean weight (b) C7eq → β/C5 Log-frequency histogram of path weights for the C7eq to beta/C5 direction, showing a heavy tail far above the mean weight
Figure 4. The weight-variance problem. Histograms (log frequency) of path weights $W(\Gamma)$ under RRS, for (a) β/C5 → C7eq and (b) C7eq → β/C5. Even after branching control, the distribution is heavy tailed: rare paths carry weights far above the mean $\langle W(\Gamma) \rangle$ (dashed), illustrating why naive reweighting is unusable and why weight control is essential.

The method search

Final estimates of each weight-control scheme against the brute-force reference.

The pattern the agent uncovered is stark. SMC (no resampling) completes the A→B direction but stalls after only 31 of 1000 trials in the reverse direction, where weight degeneracy inflates the population and the compute budget is exhausted before convergence. SIR (global resampling) never settles: its per-trial relative standard deviation is 2–3× the mean, and both directions were terminated early (118 and 533 trials) as the estimate wandered rather than converged; aggressive global resampling destroys walker diversity here. RRS is the only scheme that runs all 1000 trials in both directions and converges cleanly, with per-trial spread well below SIR’s and no catastrophic tail behavior.

Table 1Weight-control schemes versus brute force
Scheme Direction Trials Rate (fs−1) σ/μ MFPT (fs)
SMC no resampling β/C5 → C7eq 1000 1.02 × 10−3 0.70 978
C7eq → β/C5 31 † 3.91 × 10−4 0.38 2555
SIR global resample β/C5 → C7eq 118 † 1.22 × 10−3 2.18 822
C7eq → β/C5 533 † 5.25 × 10−4 3.15 1906
RRS branchingSelected β/C5 → C7eq 1000 1.05 × 10−3 0.92 950
C7eq → β/C5 1000 4.22 × 10−4 0.66 2372
Brute force exact β/C5 → C7eq n/a 1.08 × 10−3 n/a 926
C7eq → β/C5 n/a 4.48 × 10−4 n/a 2234

Table 1. Final estimates over up to 1000 trials ($N_p = 1000$ walkers, ensemble importance function, seed 2026). “Trials” is the number completed before convergence or termination; σ/μ is the per-trial relative standard deviation of the rate, a proxy for weight variance. Only RRS (highlighted) completes both directions and matches the reference. † Run terminated before 1000 trials (degeneracy / non-convergence).

The headline validation

RRS recovers the brute-force rate; the naive committor does not

The running-rate traces for the selected RRS estimator show the reweighted estimate rIS settling onto the brute-force line rBF within a few hundred trials, far above the naive committor estimate rNS.

This is the headline validation, and it is worth stating precisely: using the same learned committor two different ways gives very different answers. Applied naively (flux times the mean committor over escape configurations, “N.S.”), the estimate is badly biased: MFPTAB = 926 fs against a naive 5412 fs, and MFPTBA = 2234 fs against a naive 4028 fs, i.e. up to ~6× too slow, because an imperfect committor does not satisfy the exact flux relation. Using that same committor as an importance-sampling guide and reweighting with RRS (“I.S.”) corrects this to 950 fs (+2.6%) and 2372 fs (+6.2%), both within error of the brute-force bars (“B.F.”). The reweighting turns an imperfect learned object into a quantitatively correct rate (the property promised by the unbiased estimator), provided the weight variance is controlled, which is exactly what RRS delivers.

(a) β/C5 → C7eq · convergence
Running rate estimate versus trial for the beta/C5 to C7eq direction, converging onto the brute-force reference line
Blue: running estimate rIS with its confidence band. Red dashed: brute-force reference rBF. Light dashed: naive committor estimate rNS.
(c) β/C5 → C7eq · MFPT
Bar chart of mean first-passage time for the beta/C5 to C7eq direction comparing naive committor, RRS importance sampling, and brute-force MD
Mean first-passage time from the naive committor (N.S., 5412 fs), the RRS importance-sampling estimate (I.S., 950 fs), and brute-force MD (B.F., 926 fs).
(b) C7eq → β/C5 · convergence
Running rate estimate versus trial for the C7eq to beta/C5 direction, converging onto the brute-force reference line
Blue: running estimate rIS with its confidence band. Red dashed: brute-force reference rBF. Light dashed: naive committor estimate rNS.
(d) C7eq → β/C5 · MFPT
Bar chart of mean first-passage time for the C7eq to beta/C5 direction comparing naive committor, RRS importance sampling, and brute-force MD
Mean first-passage time from the naive committor (N.S., 4028 fs), the RRS importance-sampling estimate (I.S., 2372 fs), and brute-force MD (B.F., 2234 fs).

Figure 5. RRS recovers the exact kinetics. (a,b) The importance-sampling estimate converges onto the exact rate in both directions. (c,d) The naive estimate is up to 6× too slow, while reweighted RRS matches brute force to within 3–6%.

  • +2.6%RRS mean first-passage time against brute force, β/C5 → C7eq (950 against 926 fs)
  • +6.2%the same comparison for C7eq → β/C5 (2372 against 2234 fs)
  • ~6×too slow: the naive committor estimate for β/C5 → C7eq (5412 against 926 fs)

N.S.: flux times the mean committor over escape configurations. I.S.: the same committor used as an importance-sampling guide, reweighted with RRS. B.F.: brute-force MD from 104 unbiased walkers per region. Alanine dipeptide in vacuum, constrained overdamped Langevin dynamics at 300 K, seed 2026.

Takeaway

Found, not hand-tuned

On the β/C5 ↔ C7eq isomerization of alanine dipeptide, with an exact brute-force rate available in both directions, a Claude-driven method search compared sequential importance sampling, global resampling, and a branching resampled random walk, and selected RRS. It is the only scheme that completes both directions, controls the heavy-tailed weight distribution, and recovers the true rate to within 3–6%, where the naive committor estimate is up to 6× wrong. The scientific content is the estimator and its validation; the methodological content is that the estimator was found, not hand-tuned, by an agent reading its own variance diagnostics against ground truth.

Future directions

Where this goes next

Several limitations frame the outlook.

  1. Scale

    The benchmark is a single peptide in vacuum, and the exponential weight growth is milder here than it will be for a larger, solvated system, so the real test of RRS is whether it scales.

  2. Autotune

    The B→A direction retains a slightly larger residual bias and per-trial variance than A→B, which points to the weight band $[w_{\min}, w_{\max}]$ and the ESS trigger as hyperparameters worth autotuning; these are themselves natural targets for the agent to optimize.

  3. Beyond resampling

    The agentic loop is not specific to resampling. An agent that can propose an estimator, implement it, run it, and read back a variance diagnostic against a trusted reference is a general engine for numerical-scheme discovery in molecular simulation. Turning Claude loose on the functional form of the bias, not just the resampling, is the direction I find most promising.

Keep exploring

From theory to applications

The write-up

Agent-guided variance reduction for importance-sampled protein kinetics: recovering the alanine dipeptide β/C5 ↔ C7eq isomerization rate

Michael M. KimBuilt with Claude: Life Sciences Hackathon (Researcher Track), in partnership with Gladstone Institutes, 2026

Write-up (PDF) Video Code 
@misc{kim2026agentguided,
  title={Agent-guided variance reduction for importance-sampled protein kinetics:
         recovering the alanine dipeptide beta/C5 <-> C7eq isomerization rate},
  author={Kim, Michael M.},
  year={2026},
  howpublished={Built with Claude: Life Sciences Hackathon (Researcher Track),
                in partnership with Gladstone Institutes}
}

Acknowledgements. This work was carried out during the Built with Claude: Life Sciences hackathon, in partnership with Gladstone Institutes. I thank the organizers and the Anthropic team, and gratefully acknowledge the Anthropic API credits provided for the event. The committor and rate-estimation framework builds on prior work on importance sampling for rare transition paths (accelerated Langevin dynamics and accelerated MCMC, Kim & Cai 2026). Every result uses random seed 2026 and is reproducible from the accompanying pipeline: equilibration, brute-force MD, committor training, branching-random-walk rate estimation, and plotting.