4. Physical model

This document specifies the physics: what each layer models at its most accurate (reference) level. How each layer is computed, including cheaper numerical levels, is in 06-compute.md; the fidelity lattice (§3 there) lists the levels per layer. Every layer is a component with a defined interface (07-platform.md).

L1: Electrical (cable)

• Multicompartment neurons built from EM meshes or skeletons. Compartment length is set by the λ-rule (≤ 0.1 λ at the relevant frequency).
• Correct EM shrinkage and diameter errors using priors fitted where light-microscopy and EM morphologies overlap. Diameter uncertainty is propagated as a parameter, not ignored.
• Graded, non-spiking dynamics are first class. Most C. elegans neurons and many fly interneurons are non-spiking or plateau-generating.
• Production fly runs use an impedance-preserving reduction (~10–50 compartments) that keeps the synapse-to-output transfer impedances H2 relies on (06-compute.md §3).

L2: Ion channels (genome-anchored)

• A channel library maps each channel gene to a kinetic model: HH or Markov, temperature-dependent (Q₁₀). Sources: invertebrate channel literature and ICGenealogy-style curated kinetics. Genes without measured kinetics use homology-based priors with wide uncertainty.
• Channel presence is probabilistic, not a hard filter. Single-cell transcriptomes have false negatives, and CeNGEN sets expression below its detection thresholds to zero. For each gene g and type c, the prior on density is spike-and-slab: present with probability π_gc, with a lognormal slab whose location comes from expression level.
  • π_gc is computed from the evidence: expression across CeNGEN's four threshold levels (or Fly Cell Atlas detection rates), detection sensitivity for transcripts of that abundance, and reporter, antibody or electrophysiology data where they exist.
  • Hard exclusion (π_gc = 0) needs stronger evidence than non-detection: a validated negative reporter or knock-in, or absence of the current in that type in electrophysiology with a pharmacological or genetic control.
  • For computation, genes with π_gc < π_min (default 0.05, set at M0) are dropped from the type's mechanism signature. The dropped prior mass is recorded per type, and the truncation is checked on the worm by re-running with π_min = 0.
  • Expression level remains a prior on density, not a value.
• Density varies along the neurite as a per-type parametric function of path distance and compartment class (axon, dendrite, soma, initial segment).
• Stochastic channel gating (Markov / Langevin) in small compartments where channel counts are low. Worm neurons are small, so channel noise matters there.

L3: Chemical synapses (per synapse, per site)

• Every EM synapse is an individual site at its measured location. There is no aggregation into edge weights in the reference model.
• Presynaptic side: Ca²⁺ nanodomain → release probability (Hill), vesicle pools with replenishment, short-term facilitation and depression, stochastic binomial release. Graded release for non-spiking neurons.
• Postsynaptic side: receptor composition from the transcriptome, using Markov receptor models (e.g. nicotinic vs. muscarinic ACh; GluCl glutamate-gated Cl⁻, which is inhibitory in insects; GABA_A/B).
• Synapses can be merged for computation only under stated lumpability conditions; otherwise merging is an approximation with measured error (06-compute.md §3). Sharing presynaptic neuron, receptor type and postsynaptic compartment is necessary but not sufficient.
• Neurotransmitter identity comes from EM-based predictions (Eckstein et al. 2024) and is cross-checked against expression of transmitter synthesis and transport genes.

L4: Electrical synapses (latent)

• Gap-junction edges are latent variables:
  • Worm: EM-annotated, plus uncertainty.
  • Fly: absent from EM.
• Prior = f(innexin co-expression between the two types, membrane contact area in EM, light-microscopy innexin maps). Rectification is set by innexin pairing. Conductance is inferred.

L5: Volume transmission (neuropeptides and monoamines)

• Release from dense-core vesicles depends on activity and Ca²⁺ with slow kinetics.
• Reaction–diffusion in extracellular space. Peptide diffusion lengths over 0.1–10 s are tens of µm, much larger than the EM extracellular microgeometry, so the geometry enters through effective diffusivity (tortuosity, volume fraction) per neuropil region, computed from EM. Full EM-voxel geometry is the reference level used to validate this (06-compute.md §3).
• Uptake and degradation by peptidases and transporters.
• Targets are GPCRs, using peptide–receptor maps (C. elegans: Beets et al. 2023, Ripoll-Sánchez et al. 2023). Downstream second-messenger cascades (Gs/Gi → cAMP/PKA, Gq → PLC/IP₃/DAG) phosphorylate channels and change gating or density.

L6: Ion homeostasis, glia and energy

• Dynamic intra- and extracellular ion concentrations: Na⁺/K⁺-ATPase, NCX, PMCA, SERCA/ER stores, Ca²⁺ buffers.
• The calcium indicator (GCaMP) is part of the model. It buffers Ca²⁺ and is the readout that experiments observe. The model includes indicator binding kinetics and a fluorescence forward model, so we compare against raw fluorescence and never against deconvolved spikes.
• Glia: K⁺ spatial buffering and neurotransmitter uptake, using glia geometry where EM has it.
• Optional ATP budget to constrain activity levels.

L7: Slow regulation and plasticity (Regulatory Closure, H*)

This layer defines the developmental outcome that H* predicts: the complete regulatory dynamics, the initial state, the developmental schedule, and the rule that selects an equilibrium. The fixed-point solver in 06-compute.md §5 is a shortcut for this definition. It is valid only under the conditions stated there.

Regulatory program

For each cell type c, a regulatory program R_c includes:

• sensors: a set of k low-pass filtered local signals, mainly Ca²⁺ in different bands, after Liu et al. 1998 and O'Leary et al. 2014;
• target set points s*_c ∈ ℝᵏ, which may depend on neuromodulatory state;
• gains G_c, which map sensor errors to rates of change of the regulated quantities;
• leak D_c: decay of regulated quantities that the error signal does not drive;
• the expression cascade (mRNA → protein → membrane insertion): a stable filter between the regulatory state and the parameters;
• synaptic scaling rules, regulated in the same way as channel densities;
• an initial-state distribution p₀,c (below).

Regulatory dynamics (canonical form)

For neuron i of type c, let z_i ∈ ℝⁿ be the regulatory state. Its coordinates are the integrator variables of the n regulated quantities: channel and receptor densities in log coordinates, synaptic scaling factors, and fast-path shifts of voltage-dependence. Biophysical parameters follow from it through the steady state of the expression cascade, θ_i = h_c(z_i). Every rule, in both families below, is written as

dz_i/dt = G_c(z_i) · e_i  −  D_c (z_i − z_c^ref)  +  O(‖e_i‖²)
e_i     = E_closed-loop[ s*_c(m_t) − σ_i(t) ]          (k-vector)

where:

• σ_i(t) ∈ ℝᵏ are the cell's own sensor signals;
• m_t is its local neuromodulatory state;
• the expectation is over closed-loop activity (brain + body + environment) at fixed z.

Rules may include terms of higher order in e. Without leak, these terms vanish wherever e = 0, so they do not change the set of equilibria; they can change which equilibrium is reached ("When history matters", below).

Averaging is valid when regulation is slow relative to the mixing of closed-loop activity, including behavioral-state switching: ε = τ_mix / τ_reg ≪ 1, where τ_reg is the fastest regulatory timescale (the inverse of the largest eigenvalue magnitude of the loop gain, below). The averaged dynamics then follow the true regulated trajectory with error O(ε) over times O(1/ε). Coordinates with ε not small, typically parts of the fast path, are simulated in time together with activity, not averaged.

Rule families

Two families are fitted and compared, both in the canonical form:

1. Mechanistic: calcium-sensor rules after Liu et al. (1998) and O'Leary et al. (2014). Both are constant-gain integral control: G_c is constant (in log coordinates for Liu's multiplicative rule) and D_c = 0.
2. Learned local: small neural networks per family of cell types parameterize G_c(z), D_c, the sensor filters and the higher-order terms. Their inputs are restricted to the cell's own signals (voltage, Ca²⁺ bands, second messengers) and its own regulatory state. The number of sensors k is a hyperparameter.

The locality restriction is the hypothesis; the functional form is not. If only the learned family fits, the mechanistic rule was wrong but H* survives.

Two regulation timescales. Both families include:

• a fast path: shifts in channel voltage-dependence (e.g. by phosphorylation) over seconds to minutes, with gains G^f and a nonzero leak;
• a slow path: channel density through expression, over hours.

The fast path explains immediate partial compensation after a perturbation and relaxes as the slow path takes over; the slow path explains persistent change. F3 time courses are fitted with both.

Set points can depend on neuromodulatory state. Neuromodulators from L5 can shift a type's set points (e.g. between fed and starved states), so homeostasis maintains a state-appropriate target instead of fighting modulation. This is a per-type function of local GPCR signaling, so it stays within the locality claim.

Initial state and developmental schedule

• Initial state. Regulation starts at the onset of each neuron's terminal differentiation, from z_i(t₀) ~ p₀,c. By default p₀,c is lognormal (Gaussian in z), with its mean from the earliest available expression profile of the type and its covariance Σ₀,c inferred as part of Φ_c. Neurons and animals draw independently. Under H*, within-type and between-animal variability of parameters therefore come partly from p₀,c (05-inference.md).
• Developmental schedule. The schedule specifies, as functions of developmental time, the wiring W(t), body size and mechanics, and the rearing environment. It is an explicit model input; it includes the rearing environment and is uncertain (below).
• Developmental outcome (definition). The adult parameters are θ_i = h_c(z_i(t_adult)), where z solves the regulatory dynamics from z(t₀) ~ p₀,c along the schedule, in the closed loop. This replaces "fitting θ". Everything else in this section characterizes this outcome.

Equilibrium selection

Let M = {z : e(z) = 0} be the set where every neuron's averaged sensors meet its set points. With k sensors and n regulated quantities per neuron, M has dimension n − k per neuron. Meeting the set points therefore does not determine the parameters when k < n. What the dynamics select from M depends on the structure of the rule.

Case 1: integral control with constant gain (G_c constant, D_c = 0). This includes the mechanistic family.

• Conserved quantities. For every t, z_i(t) − z_i(t₀) lies in range(G_c), for any sensor nonlinearity, network coupling or schedule. Equivalently, the n − k quantities N_cᵀ z_i are conserved, where N_c spans null(G_cᵀ).
• Slice coordinates. Write z_i(t) = z_i(t₀) + G_c a_i(t), where a_i = ∫ e_i dt ∈ ℝᵏ is the integrated error. Regulation restricted to this slice is the k-dimensional system da_i/dt = e_i(z(t₀) + G a).
• Selection rule. The outcome is the point of M on the slice through the initial state, z(t₀) + range(G). Generically this intersection is a set of isolated points.
• Stability. Let J = ∂E[σ]/∂z be the sensitivity of the averaged sensors to the regulatory state of all neurons, including network coupling. An equilibrium on the slice is stable if every eigenvalue of the loop gain J G has positive real part. A root that fails this test is not a developmental outcome, even though it satisfies the sensor equations.
• Uniqueness. If the symmetric part of J G is uniformly positive definite on a convex region of the slice that contains the trajectory, the equilibrium in that region is unique and attracts every trajectory in it (the slice dynamics are then strongly monotone).
• Path independence. When the slice equilibrium is unique, the adult outcome depends on the initial state and on the adult wiring and environment that define M. It does not depend on the developmental schedule or on transient manipulations of activity. If the slice holds several stable equilibria, the schedule selects among them.

Case 2: leak (D_c ≠ 0).

• Linear quantities conserved by the dynamics are those in the left null space of [G_c D_c]. If D_c has full rank there are none, and equilibria are generically isolated.
• Set points are not met exactly. The steady-state sensor error satisfies G_c e = D_c (z − z_c^ref).
• Memory of the initial state and of earlier activity decays at rates set by D_c. If the slowest of these is comparable to the animal's age, the adult state is not an equilibrium and must be computed by transient simulation.

Case 3: state-dependent gain G_c(z).

• Involutive case. If the Lie bracket of every pair of columns of G_c lies in their span, trajectories stay on k-dimensional integral manifolds (leaves), by the Frobenius theorem. Case 1 holds with the leaf through z(t₀) in place of the slice. This is automatic for k = 1.
• Non-involutive case. Otherwise, the states reachable from z(t₀) form the orbit generated by the columns and their brackets, which can have dimension up to n. The intersection of this orbit with M then has positive dimension, and where the trajectory ends depends on the time course of e(t). This is holonomy: identical initial states and identical adult conditions can give different adult parameters after different histories.

When history matters

Activity history during development can leave a lasting effect on adult parameters, after conditions are restored, only through:

1. several stable equilibria on the same slice or leaf;
2. a non-involutive gain G_c(z) (Case 3);
3. terms of higher order in e; a term in e² changes the state along directions outside range(G) whatever the sign of the error;
4. active bounds on regulated quantities (saturation of expression), which break the conservation law while active;
5. state outside the regulatory program, such as associative plasticity;
6. a leak whose decay is slow compared with the animal's age (Case 2): a memory that fades rather than a lasting change.

In particular, constant-gain mechanistic rules predict no lasting effect of a transient developmental manipulation unless mechanism 1 or 4 applies. F7 tests the equilibrium shift during a sustained manipulation (F7a) separately from the effect that remains after it ends (F7b) (02-hypotheses.md).

Noise. Under Case 1 the slice conservation holds for every realization of activity, so fluctuations in e(t) only jitter the state around the slice equilibrium. Under mechanisms 2 and 3 the state is not confined in this way. Whether fluctuations then produce within-type variability that grows with age is an open question for the demonstrator. It is not assumed as a prediction: a first synthetic test with a non-involutive rule (n = 3, k = 2) showed no growth.

Identifiability from equilibria and from transients

• Gauge. Write G_c = U_c C_c, where U_c is an orthonormal basis of range(G_c) and C_c ∈ GL(k). In Case 1 every equilibrium depends on G_c only through range(G_c), so replacing G_c by G_c C for any invertible C changes no endpoint. This includes multiplying all regulation rates by a common factor.
• Equilibrium data. Steady-state data, including within-type covariance, carry zero Fisher information about C_c: overall regulation rates, their ratios within the subspace, and their mixing. They constrain s*_c, range(G_c) (a point on the Grassmannian Gr(k, n)) and p₀,c.
• Transient data. C_c and the leak D_c are identified only by transients: F3 recovery time courses and F7 time courses. M3 reports these two groups of parameters separately (09-roadmap.md).

Predicted within-type covariance

Linearize Case 1 around the slice equilibrium. Let δz₀ be the deviation of a neuron's initial state from its type mean, and δr the deviation of its averaged sensors from the type's set points that does not come from z: per-neuron set-point jitter, and differences in input and position in the circuit. The adult deviation is then

δz* = P δz₀ + G (J G)⁻¹ δr,      P = I − G (J G)⁻¹ J
Cov(z*) = P Σ₀ Pᵀ + G (J G)⁻¹ Σ_r (J G)⁻ᵀ Gᵀ

• P is the oblique projector onto null(J) along range(G). Variability inherited from development is confined to the directions in which the sensors are insensitive.
• The sensor-controlled combinations Jz vary only as much as set points and inputs do: J Cov(z*) Jᵀ = Σ_r, independent of Σ₀.

F4 tests these structures (02-hypotheses.md).

Development follows the measured connectome series, with uncertainty

Wiring changes during development, and under H* regulation tracks it. The worm has connectomes from the first larval stage to adulthood (Witvliet et al. 2021), but each comes from a different animal. Interpolating them directly would mix age-related change with individual variation. CAKE therefore uses a probabilistic developmental wiring model:

• Each edge weight is a smooth trend in developmental time plus an animal-specific deviation. The deviation's variance is estimated from datasets at the same stage, including the adult datasets of Cook et al. (2019) and Witvliet et al. (2021).
• Burn-in runs on developmental histories sampled from this model, and adult predictions are reported as distributions.
• Under a Case 1 rule with a unique slice equilibrium, the adult outcome does not depend on the sampled history. The schedule then matters only through adult wiring, branch selection and the history mechanisms above. Burn-in reports the fraction of sampled histories that end on a different equilibrium.

This gives an extra test: H* predicts how parameters change from stage to stage, and whether adult parameters depend on the developmental path.

Associative plasticity (learning)

Distinct from homeostasis, which pulls activity back to set points, learning changes synapses in response to experience and reward. It is modeled where the rules are well characterized: dopamine-gated depression of Kenyon cell → mushroom body output neuron synapses in the fly, and the identified circuits for salt and temperature associative learning in the worm. Plasticity rules are CKL mechanisms like any other.

• Learned weights are not regulatory state and are not selected by the equilibrium solve. They are produced by simulating the training protocol.
• The regulation solve is conditioned on the current learned weights when learning and regulation act on well-separated timescales. Otherwise the two are simulated jointly in time.
• Homeostatic synaptic scaling is part of z; learned weights are not.
• Tested by E4 learned behavior.

Rearing history is an explicit input

Under H*, θ can depend on the activity each neuron experienced while it settled. The burn-in therefore needs a defined rearing environment: the distribution of sensory conditions, temperature and behavior during development. Default: the standard lab conditions of the animals that produced the data (e.g. 20 °C, fed, on agar for worms).

Sensitivity analysis re-runs burn-in under alternative rearing distributions and measures how much θ and the Emulation Ladder scores change. It reports two components separately:

• the effect of the adult environment, which changes M;
• the lasting effect of history, which needs one of the mechanisms above.

Large sensitivity is a finding in its own right: it predicts that differently reared animals should differ in measurable ways (testable, e.g. in temperature-raised worms).

L8: Body and environment (closing the sensory loop)

• C. elegans:
  • Viscoelastic continuum body (finite elements or Cosserat rod), cuticle, hydrostatic pressure, 95 body-wall muscles with excitation–contraction coupling.
  • Neuromuscular junctions.
  • Proprioception through stretch-sensitive channels in motor neurons.
  • Environment: low-Reynolds-number fluids (resistive force theory → regularized Stokeslets → full Stokes for reference), agar surface contact, and chemical and thermal gradients solved as diffusion fields.
• Drosophila:
  • Articulated body from the flybody / NeuroMechFly lineage (MuJoCo / MJX) with Hill-type muscles, motor neuron → muscle mapping from the nerve cord connectome, campaniform and chordotonal mechanosensors, halteres.
  • Compound-eye optics renderer feeding photoreceptor phototransduction models.
  • Odor plumes from CFD.
  • Flight aerodynamics (quasi-steady → CFD for reference).

Timescales and conservation

• The layers span ~10 µs (spikes) to days (regulation). Time stepping, multirate splitting and partitioning are specified in 06-compute.md §4.
• Conservation of charge, ion mass and mechanical energy is a physical requirement of every layer and is monitored continuously in every run (07-platform.md, verification).