A mathematical model of intercellular signaling during epithelial wound healing

Recent experiments in epithelial wound healing have demonstrated the necessity of Mitogen-activated protein kinase (MAPK) for coordinated cell movement after damage. This MAPK activity is characterized by two wave-like phenomena. One MAPK “wave” that originates immediately after injury, propagates deep into the cell layer, and then rebounds back to the wound interface. After this initial MAPK activity has largely disappeared, a second MAPK front propagates slowly from the wound interface and continues into the tissue, maintaining a sustained level of MAPK activity throughout the cell layer. It has been suggested that the first wave is initiated by reactive oxygen species (ROS) generated at the time of injury. In this paper, we develop a minimal mechanistic diffusion-convection model that reproduces the observed behavior. The main ingredients of our model are a competition between ligand (e.g., Epithelial Growth Factor) and ROS for the activation of Epithelial Growth Factor Receptor (EGFR) and a second MAPK wave that is sustained by stresses induced by the slow cell movement that closes the wound. We explore the mathematical properties of the model in connection with the bistability of the MAPK cascade and look for traveling wave solutions consistent with the experimentally observed MAPK activity patterns.

💡 Research Summary

The paper addresses a striking phenomenon observed in epithelial wound healing: two distinct waves of MAPK (Mitogen‑activated protein kinase) activity. The first wave appears immediately after injury, spreads rapidly through the cell layer, and then rebounds toward the wound edge, disappearing within a short time. After this acute response, a second, slower‑moving MAPK front emerges from the wound margin and propagates steadily across the tissue, maintaining a sustained level of MAPK activation throughout the epithelium.

To capture these dynamics, the authors construct a minimal mechanistic model that combines diffusion, convection, and nonlinear receptor activation. The model rests on three biological premises. First, the epidermal growth factor receptor (EGFR) can be activated either by reactive oxygen species (ROS) generated at the moment of injury or by a soluble ligand such as epidermal growth factor (EGF) that is secreted continuously by cells surrounding the wound. ROS and ligand compete for the same receptor binding sites; ROS appears as a short‑lived, high‑amplitude pulse that diffuses quickly, whereas the ligand forms a slowly varying background. The binding probabilities are represented by Hill‑type functions, providing a nonlinear activation term for downstream MAPK signaling.

Second, the MAPK cascade itself is modeled as a bistable switch. Low MAPK activity is self‑inhibitory, while once a threshold is crossed the system flips to a high‑activity state due to positive feedback loops (e.g., Raf‑MEK‑ERK amplification). This bistability is introduced through a cubic reaction term that yields two stable fixed points separated by an unstable one. The presence of two stable states is essential for the formation of traveling fronts: a high‑MAPK region can invade a low‑MAPK region and vice‑versa, depending on the local stimulus.

Third, the authors incorporate mechanical stress generated by collective cell migration. As cells crawl to close the wound, they generate traction forces that feed back onto EGFR activation, effectively providing a second, slowly varying stimulus that sustains MAPK activity after the ROS pulse has vanished. This stress‑induced activation is modeled as a convection term proportional to the cell‑velocity field, which itself is linked to the spatial gradient of MAPK activity (cells move faster where MAPK is high).

Mathematically, the system consists of coupled partial differential equations for ROS concentration (R), ligand concentration (L), EGFR activation (E), and MAPK activity (M). The ROS equation includes diffusion (coefficient D_R) and a rapid decay term; the ligand equation features slower diffusion (D_L) and a source term at the wound edge. EGFR activation obeys a Michaelis‑Menten‑type kinetics with competitive inhibition by ROS and ligand. MAPK dynamics follow a reaction‑diffusion‑convection equation: ∂M/∂t = D_M∇²M + v·∇M + f(M, E), where f encodes the bistable cubic nonlinearity.

The authors perform a linear stability analysis to identify parameter regimes where the homogeneous low‑MAPK state becomes unstable in the presence of a ROS pulse, giving rise to a fast traveling wave. By imposing a reflecting boundary at the wound edge, they derive an analytical expression for the wave speed of the first pulse, showing it scales with √(D_R·k_R), where k_R is the ROS‑induced activation rate. Numerical simulations confirm that this wave propagates deep into the tissue, reflects, and then decays as ROS is depleted.

For the second, slower wave, the analysis focuses on the convection term generated by cell movement. When the stress‑induced activation exceeds a critical value, the high‑MAPK state becomes self‑sustaining and a traveling front moves outward with speed proportional to the cell‑migration velocity v_c. Because diffusion of MAPK is relatively weak (small D_M), the front retains a sharp profile, matching experimental observations of a persistent MAPK “front” that advances at a few micrometers per hour.

Parameter sweeps reveal that the coexistence of both waves requires (i) a sufficiently high initial ROS burst to trigger the first front, (ii) a ligand concentration that maintains a basal EGFR activation, and (iii) a cell‑migration speed that supplies enough mechanical feedback to keep the high‑MAPK state alive. If any of these conditions are weakened—e.g., ROS scavenging, EGFR inhibition, or impaired cell motility—the model predicts loss of one or both waves, consistent with pharmacological experiments that block ROS or EGFR signaling.

In summary, the paper delivers a concise yet powerful mathematical framework that unifies biochemical (ROS, growth factor) and biomechanical (cell‑generated stress) cues into a single set of diffusion‑convection equations. By exploiting the intrinsic bistability of the MAPK cascade, the model reproduces the experimentally observed biphasic MAPK activity: an early, ROS‑driven wave that quickly traverses and rebounds, followed by a later, stress‑driven wave that persists and drives coordinated cell migration. The analysis not only clarifies the mechanistic underpinnings of wound‑induced MAPK dynamics but also offers quantitative predictions (wave speeds, threshold ROS levels, required migration velocities) that can guide future experimental interventions and therapeutic strategies aimed at modulating wound healing.