Extended Subsystem Couplings
The subsystem models developed in Energy Metabolism Models and Neuroendocrine and Autonomic Models have been expanded with mast cell dynamics (Mast Cell Activation Dynamics), coagulation (Coagulation and Microvascular Dynamics), mitochondrial quality control (Mitochondrial Dynamics and Quality Control), metabolic flexibility (Lactate Kinetics and Metabolic Flexibility), BH₄ competition (Tetrahydrobiopterin Competition Model), cerebrovascular autoregulation (Cerebral Blood Flow Autoregulation), and central sensitization (Central Sensitization and Pain Amplification). This section formalizes the additional couplings these subsystems introduce into the integrated model.
1 Mast Cell–Energy–Autonomic Triangle
Mast cell degranulation couples to both the autonomic and energy subsystems simultaneously. Histamine-mediated vasodilation reduces MAP, triggering compensatory sympathetic activation (baroreflex). The sympathetic surge increases cardiac energy demand while the mast cell inflammatory mediators impose an immune energy cost (immune energy). The model reveals a three-way interaction invisible to verbal reasoning: a mast cell flare simultaneously depletes the energy budget (mediator clearance + immune activation + sympathetic compensation), worsens orthostatic intolerance (vasodilation + tachycardia), and drives neuroinflammation (TNF-\(\alpha\), IL-6 from mast cells crossing the BBB). The integrated model quantifies the total energy cost of a single mast cell flare as \(E_\text{MC} = \int_0^{\tau_\text{flare}} (J_{\text{immune,MC}} + J_{\text{SNS,comp}} + J_\text{clearance}) , d t\), which for moderate flares consumes 15–25% of the daily energy budget, explaining why ME/CFS patients with MCAS comorbidity have substantially reduced activity tolerance.
2 Coagulation–Oxygenation–Mitochondrial Cascade
Microclot burden reduces effective oxygen delivery (o2 microclot), which reduces maximal ETC flux (\(J_{\text{CIV,max}} \propto \text{DO}_2^\text{eff}\)). Reduced ETC flux increases ROS production (ros production), which damages endothelial cells, increasing tissue factor expression and promoting further coagulation. This cascade creates a slow positive feedback loop operating on timescales of days to weeks—a feedback loop that only becomes apparent when the coagulation and energy models are coupled. The model predicts that the cascade has a threshold: below a critical microclot burden \(M_c^*\), fibrinolysis keeps pace with formation; above \(M_c^*\), the feedback loop becomes self-sustaining. This threshold is patient-specific, depending on fibrinolytic capacity (PAI-1 levels), inflammatory state, and baseline mitochondrial reserve.
3 Connective Tissue and Vascular Compliance
Hypermobile Ehlers–Danlos syndrome (hEDS) is overrepresented among ME/CFS patients. The biomechanical consequences of connective tissue laxity enter the model through modified vascular compliance parameters. Increased venous compliance (\(C_v\)) amplifies gravitational blood pooling during orthostatic stress:
\[ \Delta V_\text{pool} = C_v \cdot \rho g \cdot h_\text{column}, \quad C_v^\text{hEDS} = \kappa \cdot C_v^\text{normal}, \quad \kappa \in [1.3, 2.0] \tag{1}\]
where \(C_v\) is venous compliance (volume/pressure), \(\rho g h_\text{column}\) is the gravitational hydrostatic pressure over column height \(h_\text{column}\), and \(\kappa\) is the compliance amplification factor in hEDS. The increased pooling volume enters the orthostatic model (orthostatic) as a larger effective \(\Delta P_\text{grav}\), requiring greater compensatory sympathetic activation to maintain MAP. The model predicts—uniquely through quantitative analysis—that hEDS patients have a narrower energy envelope because the energy cost of maintaining upright posture is amplified: greater sympathetic drive means higher cardiac work, and the compensatory tachycardia consumes additional ATP. For \(\kappa = 1.5\), the model predicts a 10–20% reduction in energy available for activity during upright hours, independent of any mitochondrial or immune dysfunction. This provides a quantitative explanation for why hEDS patients with ME/CFS tend toward greater severity: the connective tissue deficit imposes a permanent energy tax on postural maintenance.
Additionally, increased joint laxity raises proprioceptive error, which the autonomic system compensates for through increased muscle co-contraction. This adds a term \(J_\text{proprioceptive} = e_\text{co-contract} \cdot (1 - \text{stability} \\/ \text{stability}_0)\) to the energy demand, further narrowing the energy envelope.
4 Epigenetic State as a Slow Variable
DNA methylation and histone modification changes have been documented in ME/CFS Vega, Vernon, and McGowan (2021) and likely contribute to disease persistence. The preceding models treat disease parameters (e.g., \(\alpha_\text{CI}\), \(k_\text{exh}\), \(n_F\)) as fixed constants. In reality, these parameters are partially determined by the epigenetic state, which itself evolves on timescales of months to years in response to the disease state.
The epigenetic state vector \(\mathbf{\mathcal{E}} = (\mathbf{\mathcal{M}}, \mathcal{A})\) comprises a per-locus methylation vector \(\mathbf{\mathcal{M}} = (m_1, m_2, ..., m_n)\) and a histone acetylation index \(\mathcal{A}\). For full specification of \(\mathbf{\mathcal{M}}\) — including locus classes (ProB RepSeqs, ProA RepSeqs, gene-region), the DNMT3B zero-sum constraint, the irreversibility threshold \(m_i^\text{crit}\), the derived \(B_\text{strength}\) variable, tissue/cell-type indexing, and histone mark coupling — see the formal model in Chapter Formal Causal Hierarchy Analysis (Per-Locus Dynamics: Vector Model for Bidirectional Methylation). The scalar \(overline(\mathcal{M})\) used in this chapter is retained as a first approximation for the integrated ODE system; it represents the net deviation \(||\mathbf{\mathcal{M}} - \mathbf{\mathcal{M}}^\text{baseline}||\) and should be interpreted as a compressed summary, not a complete epigenetic description. \(overline(\mathcal{M})\) and \(\mathcal{A}\) modulate key disease parameters:
\[ \begin{aligned} \alpha_\text{CI}(\mathcal{E}) &= \alpha_{\text{CI,0}} \cdot (1 - \epsilon_\text{CI} \cdot overline(\mathcal{M})) \\ k_\text{exh}(\mathcal{E}) &= k_{\text{exh,0}} \cdot (1 + \epsilon_\text{exh} \cdot overline(\mathcal{M})) \\ J_\text{biogenesis}(\mathcal{E}) &= J_{\text{biogenesis,0}} \cdot (1 + \epsilon_\text{bio} \cdot \mathcal{A}) \end{aligned} \tag{2}\]
where the \(\epsilon\) coefficients quantify the sensitivity of each parameter to epigenetic state. The epigenetic state evolves in response to the disease state itself:
\[ \begin{aligned} \frac{d overline(\mathcal{M})}{d t} &= k_\text{DNMT} \cdot \mathbf{C}_\text{pro} \cdot (1 - overline(\mathcal{M})) - k_\text{demeth} \cdot overline(\mathcal{M}) \\ \frac{d \mathcal{A}}{d t} &= k_\text{HAT} \cdot \frac{[\text{AcCoA}]}{K_\text{AcCoA} + [\text{AcCoA}]} \cdot (1 - \mathcal{A}) - k_\text{HDAC} \cdot \mathcal{A} - k_\text{SIRT} \cdot \frac{[\text{NAD}^+]}{K_\text{SIRT} + [\text{NAD}^+]} \cdot \mathcal{A} \end{aligned} \tag{3}\]
Scalar approximation caveats. The \(d overline(\mathcal{M})/d t\) equation treats \(overline(\mathcal{M})\) as driven upward by DNMT3A/B activity (gain-model), with passive demethylation as reversal. This approximation omits three structural features captured in the full per-locus vector model (Chapter Formal Causal Hierarchy Analysis, Per-Locus Dynamics: Vector Model for Bidirectional Methylation): (1) DNMT3B redistribution — in cancer, loss at ProB repeats and gain at ProA/gene-region are coupled through a shared enzyme pool (\(\sum_i \Delta m_i = 0\) at the allocation level), not independent; for ME/CFS, the mechanism is proposed as more targeted (HSAT2-specific loss-meCpG, Geneviève Fourel pers. comm. May 2026); (2) irreversibility threshold — loci falling below \(m_i^\text{crit}\) lose the MeCP2/MBD-DNMT1 maintenance loop and require active, SAM-dependent de novo remethylation, which passive demethylation kinetics cannot represent. Low-density CpG regions are empirically more vulnerable to maintenance failure (Tiedemann et al. 2024); DNMT1 operates via density-dependent allosteric switch (Kimura and Sasaki 2012); (3) direction-dependent therapeutic response — the scalar model implies demethylation is always therapeutic and methyl donors are always neutral, while the vector model shows that ProB loss-dominant patients require remethylation (methyl donors) and are harmed by demethylating agents (5-azacitidine, HDAC inhibitors, TET activators). The scalar model’s hysteresis and intervention-window predictions remain qualitatively valid for gain-dominant patterns but may invert for loss-dominant patterns.
DNMT3B substrate specificity. The current model treats methylation as CpG-only for simplicity. However, DNMT3B also methylates non-CpG sites (CAG trinucleotides) in early development and in neurons. A complete model would distinguish CpG from CAG methylation at DNMT3B target sites; the vector model supports this by tracking per-site entries with distinct dynamics.
Provisionality. The ProA/ProB repeat framework (Bonnet, Hulo, Fourel et al. 2026 preprint (Bonnet et al. 2026)) on which the loss-dominant scenario rests is a computational genomics hypothesis developed in cancer, not peer-reviewed. Its extension to ME/CFS is our extrapolation, not endorsed by those authors.
Basis of the gain-model formulation. In the gain-model formulation above, pro-inflammatory cytokines promote DNA methylation (DNMT upregulation), acetyl-CoA availability limits histone acetylation (HAT activity), and SIRT1/SIRT3 (NAD⁺-dependent deacetylases) modulate histone state. This model produces hysteresis—a phenomenon that can only be identified through dynamical systems analysis: the disease state promotes epigenetic changes that stabilize the disease parameters, which in turn maintain the disease state. Even if the original trigger (e.g., viral infection) is fully resolved, the epigenetic state \(\mathbf{\mathcal{E}}\) may have shifted to values that keep the system in the disease attractor. Recovery then requires not just removing the trigger but reversing the epigenetic changes—a process that operates on the slow timescale \(\tau_\text{epi} ~ 1 \\/ k_\text{demeth} ~\) months to years.
The hysteresis width—the difference between the perturbation required to trigger disease onset and the intervention required for recovery—is determined by the epigenetic coupling strengths \(\epsilon\). The model predicts that early intervention (before epigenetic consolidation) requires substantially less therapeutic intensity than late intervention, formalizing the clinical intuition of an “intervention window” in the first months after disease onset.
5 Updated Whole-Body State Vector
The extended model incorporates the new subsystem variables:
- Energy metabolism (12 variables): original 8 \(+\) \([L]\) (lactate), \(\Phi\) (fuel selection), \(M_h\), \(M_d\) (mitochondrial mass)
- Immune system (17 variables): original 12 \(+\) \(\text{MC}_r\), \(\text{MC}_p\), \(\text{MC}_d\) (mast cells), \([\text{His}]\) (histamine), PGD₂
- Cytokines (6 variables): unchanged
- Coagulation (3 variables): \([\text{Fbn}]\), \(M_c\) (microclot burden), \([\text{Pls}]\) (plasmin)
- Neuroendocrine (10 variables): original 8 \(+\) \([\text{BH}_4]\), \([\text{NO}]\)
- Autonomic/cardiovascular (5 variables): original 4 \(+\) CBF
- Pain/sensitization (2 variables): \(\mathcal{S}\) (sensitization state), \(\mathcal{F}\) (fiber density)
- Sleep (2 variables): unchanged
- Gut (5 variables): unchanged
- Epigenetic (n+1 variables): \(\mathbf{\mathcal{M}}\) (n per-locus methylation entries) and \(\mathcal{A}\) (histone acetylation). For the integrated ODE, the scalar summary \(overline(\mathcal{M})\) is used as a compressed approximation; see Per-Locus Dynamics: Vector Model for Bidirectional Methylation for the full vector specification.
The total system now comprises approximately 64 state variables. The additional 19 variables (relative to the initial 45-variable model) capture mast cell/histamine dynamics, coagulation, mitochondrial quality control, metabolic flexibility, neurovascular regulation, pain processing, epigenetic memory, and GI motility/SIBO dynamics. While the system is larger, its modular structure (each subsystem has \(\leq 5\) coupling variables with others) maintains computational tractability: the Jacobian is sparse, enabling efficient numerical integration and bifurcation analysis.