Discrete models of autocrine cell communication in epithelial layers

Abstract

Pattern formation in epithelial layers heavily relies on cell communication by secreted ligands. Whereas the experimentally observed signaling patterns can be visualized at single-cell resolution, a biophysical framework for their interpretation is currently lacking. To this end, we develop a family of discrete models of cell communication in epithelial layers. The models are based on the introduction of cell-to-cell coupling coefficients that characterize the spatial range of intercellular signaling by diffusing ligands. We derive the coupling coefficients as functions of geometric, cellular, and molecular parameters of the ligand transport problem. Using these coupling coefficients, we analyze a nonlinear model of positive feedback between ligand release and binding. In particular, we study criteria of existence of the patterns consisting of clusters of a few signaling cells, as well as the onset of signal propagation. We use our model to interpret recent experimental studies of the EGFR/Rhomboid/Spitz module in Drosophila development.

FIGURE 1

(A) The geometry of the problem is motivated by cell communication in Drosophila oogenesis. Epithelial cells cover the large oocyte. EGF receptors are uniformly distributed across the epithelial layer and are absent on the surface of the oocyte. (B) Ligands diffuse in a thin gap between the epithelial layer and a reflective surface. (C) The two model layers considered in this article—periodic arrays of squares and hexagons. (D) Main processes in ligand binding and transport. (E) Ligand binding stimulates ligand release. Receptor activation leads to the degradation of a factor inhibiting the transcription of the ligand-releasing protease. In the absence of inhibition, the protease is synthesized and generates the secreted ligand.

FIGURE 2

The steady-state ligand field due to a single ligand-releasing cell. (A) Square cells, L_x = 5 × 10⁻⁴ cm. (B) Hexagonal cells, 2L_v = 5 × 10⁻⁴ cm. Other parameters: h = 5 × 10⁻⁵ cm, k_e = 0.1 min⁻¹, k_off = 0.1 min⁻¹, R₀ = 1 × 10⁴ molecules/cell surface, D = 1 × 10⁻⁷ cm² s⁻¹, k_on = 0.1 nM⁻¹ min⁻¹, and Q_s = 100 molecules/cell/min. The fields were computed by solving Eq. 12 using an adaptive mesh finite element package (FEMLAB).

FIGURE 3

Normalized coupling coefficients (, see text for details) plotted as a function of cell-to-cell distance, k, defined as the ratio of the Euclidian distance divided by the cell size. (A) Effect of ligand-receptor affinity k_on (square cells, L_x = 5 × 10⁻⁴ cm, D = 1 × 10⁻⁷ cm² s⁻¹). (B) Effect of ligand diffusivity D (hexagonal cells, 2L_v = 5 × 10⁻⁴ cm; k_on = 0.1 nM⁻¹ min⁻¹). All other parameters are as in Fig. 2.

FIGURE 4

Stable localized patterns generated by the positive feedback. Existence of the localized three-cell patterns for squares (A) and hexagons (B) as a function of ligand release rate and ligand-receptor affinity. The shaded area corresponds to the stable localized pattern. Crossing the upper boundary ignites the neighboring cells. Crossing the lower boundary leads to the extinction of the pattern. (C) A one-dimensional cut (Q_s = 300 molecules/cell/min) through the diagram in A showing the disconnected region of existence of the localized pattern. Localized patterns are realized both for low and high ligand-receptor affinities. (D) Critical rate of ligand release necessary to destabilize the localized pattern as a function of the number of cells in it, k_on = 0.1 nM⁻¹ min⁻¹. A–D: C_T = 500 molecules/cell. Other parameters are as in Fig. 2.

FIGURE 5

Response of a cellular layer to a localized perturbation in protease release. (A) The structure of the perturbation. B–D, The ligand fields induced by the perturbation computed for different ligand-release rates. (B) Q_s = 100 molecules/cell/min. Ligand field 40 min after the perturbation. The perturbation decays. (C) Q_s = 200 molecules/cell/min. The perturbation generates a stable pattern of the same structure. (D) Q_s = 300 molecules/cell/min. The perturbation generates an “ignition” wave (see text for details). All parameters are as in Fig. 4 A.

FIGURE 6

(A) Solution of Eq. 12 for the square cells was derived on the two-dimensional semi-infinite domain with the origin in the center of the square. The dimensionless area of the square is . (B) Solution of Eq. 12 for the hexagons was derived on the two-dimensional, semi-infinite domain with origin in the center of the hexagon. The dimensionless area of the hexagon is a = 3 × r × v. (C) Indexing used for the computation of coupling coefficients for square cells. (D) Indexing used for the computation of coupling coefficients for hexagonal cells.