Spatial Turing-type Pattern Formation in a Model of Signal Transduction Involving Membrane-based Receptors Coupled by G Proteins

Abstract

In this paper, a model of signaling pathways involving G proteins is investigated. The model incorporates reaction-diffusion mechanisms in which various reactants participate inside and on the extra-cellular surface membrane. The messenger molecules may diffuse over the surface of the cell membrane and signal transduction across the cell membrane is mediated by membrane receptor bound proteins which connect the genetically controlled biochemical intra-cellular reactions to the production of the second messenger, leading to desired functional responses. Dynamic and steady-state properties of the model are then investigated through weakly nonlinear stability analysis. Turing-type patterns are shown to form robustly under different delineating conditions on the system parameters. The theoretical predictions are then discussed in the context of some recently reported experimental evidence.

Keywords

signal transduction Turing pattern weakly nonlinear stability G proteins membrane based receptors

Introduction

The ability to respond appropriately to signals in the environment is essential for the survival of any organism. A very sophisticated mechanism has therefore evolved for detecting external signals, transducing them into internal signals and eliciting cellular responses. As evidenced by several earlier investigations (Levchenko and Iglesias, 2002; Rapple et al. 2002; Iglesias, 2003; Krishnan and Iglesias, 2003), the cells can sense spatial gradients as temporal changes in receptor occupancy and change certain behavioral patterns in response. Several biological systems have such ability to sense the direction of external chemical sources and respond by transient activation of specific intracellular signaling pathways. This ability to adapt to different levels of external stimuli, so that it is the gradient of signaling molecule rather than the average signal value that determines the response, is a common feature of most chemotactic signaling systems (Levchenko and Iglesias, 2002).

The amoeboid organism Dictyostelium discoideum has been widely recognized as a useful model system for the study of chemotaxis (Rapple et al. 2002; Iglesias, 2003). Several recent research reports have proposed and studied mathematical models of signal transduction applied to Dictyostelium discoideum and other eukaryotic cells (Levchenko and Iglesias, 2002; Rapple et al. 2002; Krishnan and Iglesias, 2003). In his study of feedback control in intracellular signaling pathways in Dictyostelium discoideum, Iglesias (Iglesias, 2003) stated that, for cells to sense and respond to change in their environment, they must first have external sensors on receptors for each of the different stimuli to which it needs to respond. These external messages are then relayed to a series of internal reactants, which in turn trigger key cellular functions. A healthy functioning cell signaling mechanism, termed signal transduction, is essential for the well-being of the life form. Abnormalities of signal transduction pathways have been linked to the development of many serious disorders, such as cancer for example.

Mostly, a cancer cell is a cell that has escaped the controls that maintain its normal differentiated function within the regulatory mechanism of the body (Norman and Litwack, 1997). This therefore underlines the role of signal transduction in the loss of organismic control. In other words, cancer may be viewed as an aberration of the signal transduction process since it derives from a cell that has lost the ability to respond normally to controls from outside, or inside, the cell (Norman and Litwack, 1997). It is thus not surprising that hormones and their receptors are intimately related to carcinogenesis. Many tumors produce ectopical amounts of biologically active hormones that create dysfunctions of the signal transduction process leading to abnormal effects. Various tumors have been shown to secrete ACTH (Adrenocorticotrophic hormone) and to cause hypercortisolism, even when the tumor is undetectably small for many years (Norman and Litwack, 1997).

Additionally, hormones and antihormones are used to treat certain types of cancer. Many cancers are related to the status of hormones in the body. Hormone dependencies of a cell usually are a feature of the presence of the cognate receptor in the cell, while hormone independence becomes a feature of tumor cells that no longer express the appropriate receptor. An avenue for cancer treatment is to remove the grand responsible for the offensive secretion of a hormone. Another avenue would be to utilize appropriate hormones as chemotherapeutic agents. For an example, tamoxifen can interfere with offensive effects of estrogen and result in the inhibition of cellular growth of the tumor. For another example, Vasopressin has been proposed for its potential effect of slowing down the flow of blood that tumors depend on for growth.

Since hormones and their receptors are so closely related to carcinogenesis, better understanding of signal transduction mechanisms have been subject of recent intense research (Norman and Litwack, 1997; Iglesias, 2003). However, signaling pathways are highly nonlinear processes, involving feedbacks and cross-talk among interconnected components, subject to control by many independent events, making clear-cut description highly difficult. New techniques and approaches are needed to understand these complicated networks. According to J. Michael Bishop, “In order to fully understand these pathways, we need a convenient and powerful model to compliment the experimental research” (Iglesias, 2003).

Our goal is thus to extend current understanding of the signal transduction processes by creating a single model explaining adaptation and gradient sensing. Not unlike chemotaxis, subcellular response calls for the cells to detect often exceedingly shallow and changing gradients of extracellular substances and regulating a complex response in accordance with the direction and the value of these gradients. Such highly complex and integrated response needs an explanation and modeling in order to translate biochemical observations and clinical evidence into a set of predictions of dynamical and steady-state properties of the system.

Based on observation and principles of the signal transduction pathways proposed in earlier investigations (Norman and Litwack, 1997; Spiegel, 2000; Levchenko and Iglesias, 2002; Rapple et al. 2002; Iglesias, 2003), we may arrive at a model which describes gradient sensing and adaptation. The resulting model incorporates reaction-diffusion mechanisms in which various reactants participate inside and on the surface of the cell membrane. The messenger molecules diffuse over the extra-cellular membrane surface in two dimensions, while some transport of molecules across cell membrane may take place to a certain extent. Signal transduction across the cell membrane is mediated by membrane receptor bound proteins which connect the genetically controlled biochemical reactions in the cytosol to the production of the second messenger, eliciting desired intracellular responses. Dynamic and steady-state properties of the model are then investigated through the application of weakly nonlinear stability analysis. We show that Turing-type patterns will be formed robustly under different physiological conditions reflected by the system parameters.

The Turing mechanism has often been put forward as a model for certain aspects of morphogenesis (Glansdorff and Prigogine, 1971; Meinhardt, 1982; Koch, 1994; Kondo and Asai, 1995; Sawai et al. 2000; Raunch and Millonas, 2004; Pansuwan et al. 2005) such as pre-patterning in the embryo (Raunch and Millonas, 2004). It has also been a basis for models of self-organization in several other physical systems (Glansdorff and Prigogine, 1971), supported by verifiable observations in a real chemical system (Castets et al. 1990). The investigations involved the possibility of an instability occurring in purely dissipative systems of chemical reactions far from equilibrium and the transport process of diffusion but no hydrodynamic motion. If the system involves two chemical species, commonly termed an activator and an inhibitor, the existence of diffusive instabilities of this sort requires an autocatalytic reaction for the activator and a diffusive advantage for the inhibitor (Turing, 1952). That is, it is necessary that the activator species diffuses significantly more slowly than the inhibitor.

A more recent work by Rauch and Millonas, (2004) showed the formation of Turing-type patterns in trans-membrane signal transduction, by numerical simulations of the models based on their spectrum analysis of the linear stability theory. However, their work involves signalling in the context of a population of cells which are considered as point sources rather than at the subcellular level. Also, a great deal of useful information may be lost through the assumptions of linear stability analysis. In this paper, we are able to predict rhombic and hexagonal planforms in messenger hormone patterns on the cell membrane through nonlinear stability analysis. The theoretical results are then discussed in the context of experimental observations reported in recent literatures.

G Protein-Coupled Signal Transduction

Since G proteins couple receptors for many hormones to effectors that regulate second messenger metabolism, we shall use heterotrimeric G protein pathways as a model of signal transduction mechanisms. G proteins couple hundreds of receptors for hormones, growth factors, neurotransmitters, odorants and other extracellular ‘first messengers’ to effectors such as adenylyl cyclase (AC), phospholipase Cβ, and various ion channels (Spiegel, 2000). In recent years, an increasing number of human disorders, particularly endocrine diseases, have been shown to be caused by mutations in either G protein or G protein-coupled receptors (GPCRs).

G protein-coupled signal transduction involves the following components: GPCRs, the G proteins, G protein-regulated effectors such as AC, and a family of proteins that activates G protein guanosine triphosphatase (GTPase) activity (Spiegel, 2000).

G proteins are heterotrimers composed of three subunits, α, β and γ. The latter two subunits form a tightly, but noncovalently associated functional unit, the β/γ dimer. In the resting state, the α-subunit tightly binds guanosine diphosphate (GDP) and is associated with the β/γ dimer. When a hormone or other first messenger binds to a receptor, the receptor causes the G protein to exchange GDP for the nucleotide guanosine triphosphate (GTP) which activates the G protein. The GTP-bound α-subunit dissociates from the β/γ dimer. In Figure 1, a schematic description of the pathways, based on a scheme shown in a report by A.M. Spiegel, (2000), is shown, where the arrows between the GTP-bound α-subunit and effector and between the β/γ dimer and effector indicate the interaction between the effector(AC) and the respective subunits. Thus, G proteins are the true signal transducers since they respond to the occupancy of the receptor and modulate the activity of the catalytic subunit of the membrane-bound adenylate cyclase(AC) enzyme (Norman and Litwack, 1997). Under physiologic conditions, effector (AC) activation by G protein subunits is transient and is terminated by the GTPase activity of the α-subunit. The latter converts bound GTP to GDP, thus returning the α-subunit to its inactivated state with high affinity for the β/γ dimer which reassociates to again form the heterotrimer (Spiegel, 2000).

Figure 1

The G protein GTPase cycle (Figure is based on a scheme in the report by A.M. Spiegel (2000)). (1) Synthesis and targeting of components, (2) Receptor activation by signaling molecule, (3) Receptor activation of G protein, (4) G protein subunits-effector interaction and (5) GTPase activity terminates the interaction and returns a subunit to inactivated state.

The physical interaction of the G protein subunits, with the membrane-bound adenylate cyclase (AC) activates the catalytic activity of AC so that the substrate ATP is converted to produce cyclic adenosine mono-phosphate (cAMP) which functions as a second messenger. Intracellular responses are subsequently elicited. It activates protein kinase A, freeing its active C subunits. These liberated C subunits in turn function as an active kinase to phosphorylate proteins, which thereby amplifies the initial hormonal signal. This amplification effect is an important factor in the signaling pathways and will be incorporated into our model construction later (Equation (2.13)).

To arrive at a model for the signal transduction process, we are guided initially the work of Iglesias, (2003) and let R stand for the concentration of regulators in their inactive forms. We may think of these as the inactivated units of adenylate cyclase. The number of activated regulators is then represented by R^*(t). We assume that the activation and inactivation are modulated by a pair of agents A, and I, respectively. In our particular pathway of interest, the G proteins are the true signal transducers that switch effectors on and off. According to Norman and Litwack, (1997), the α-subunits of the G protein take on two functions: G_αs of the amount A, which activates the AC, and the active G_αi of the amount 7, which inhibits the activity of AC. Both G_αs and G_αi are activated by guanosine triphosphate (GTP) and both also functions as a GTPase. The GTPase activity thus endows the G protein with a turn on-off mechanism (Norman and Litwack, 1997). Therefore, R transforms to R^* according to the following reaction equation.

R ⇄_{A}^{I} R^{*}

(2.1)

Using mass action dynamics, an equation for the above reaction is

\frac{d R^{*}}{d t} = - k_{- r} I R^{*} + k_{r} A (t) R (t)

(2.2)

Assuming that the total number R_T = R + R^* of the regulators remains constant, (2.2) becomes

\frac{d R^{*}}{d t} = - [k_{- r} I + k_{r} A (t)] R^{*} (t) + k_{r} A (t) R_{T}

(2.3)

On the other hand, the activating agents and inhibiting agent I are regulated by the external signal which is proportional to the membrane surface concentration S of the signaling hormone, such as the external cAMP, giving rise to the following equations:

\frac{d A}{d t} = - k_{- a} A + k_{a} S

(2.4)

\frac{d I}{d t} = - k_{- i} I + k_{i} S

(2.5)

where the first terms in (2.4) and (2.5) are the rate of removal by natural means and the last terms are rates of their synthesis.

It is commonly assumed (Iglesias, 2003) that in such a process, the activated regulator R^* and A equilibrate relatively quickly while I has the slower dynamics. Using equation (2.3), one finds that R^* will equilibrate to the value

R^{*} = \frac{R_{T} A}{A + K_{R} I}

(2.6)

where

K_{R} = \frac{k_{- r}}{k_{r}}

(2.7)

From setting

\overset{\cdot}{A} = 0

in (2.4), one finds

A= \frac{k_{a}}{k_{- a}} S

(2.8)

which transforms (2.6) into

R^{*} = \frac{\tilde{k}}{{\tilde{b}}_{1} S + I} S

(2.9)

where

\tilde{k} = \frac{R_{T} k_{a}}{K_{R} k_{- a}}

and

{\tilde{b}}_{1} = \frac{k_{a}}{k_{- a} K_{R}}

Now, as the subunits of G protein interact with AC aggregates, they activate AC to synthesize intra-cellular cAMP, the level of which is denoted by C(t), described by the following equation:

\frac{d C}{d t} = - k_{- c} C + k_{c} R^{*} R^{*} + k'_{c}

(2.10)

in which the first term corresponds to the removal rate. The last two terms correspond to its synthesis, k′_C being the apparent zero order synthesis rate. The second term is the synthesis rate arising from the interaction between the activating subunits G_αs of the G protein of the amount A given by (2.8) among the activated catalytic subunits of AC whose amount is R^* as given in (2.9).

This enzyme is also found to equilibrate quickly to

C = \frac{k_{C}}{k_{- C}} R^{*} R^{*} + \frac{k'_{C}}{k_{- C}}

(2.11)

Using (2.8) in (2.11), one obtains

C = \frac{{\tilde{b}}_{2} S^{2}}{{({\tilde{b}}_{1} S + I)}^{2}} + K_{C}

(2.12)

at equilibrium, where

{\tilde{b}}_{2} = \frac{k_{C} k_{a} \tilde{k}}{k_{- C} k_{- a}}

and

K_{C} = \frac{k'_{C}}{k_{- C}}

This cAMP (C) in turn acts as a second messenger and amplifies the initial hormone signal S(t,x,y). Thus, the rate equation for the signaling hormone level at the point (x,y) on the extracellular membrane-surfaces at time t should read as follows.

\frac{\partial S}{\partial t} = - {\tilde{a}}_{3} S - \frac{{\tilde{b}}_{3} S}{{\tilde{b}}_{4} + S} + k_{S} C (t) + \tilde{μ} \nabla^{2} S

(2.13)

where the first term is the rate of removal by natural means, while the second accounts of its transport through the cell membrane causing its disappearance from the extra-cellular membrane surface. Since it is reasonable to expect this to saturate as the level of hormone increases, the absorption term assumes the form of a Hill-type saturation function, b₃ being the maximum absorption rate and b₄ the half-saturation constant. The third term accounts for the signal amplification arising from the synthesis of cAMP mentioned earlier. We have incorporated the diffusion effect of S over the 2–D cell-membrane surface by the last term in (2.13), where $Δ^{2} = \frac{\partial^{2}}{\partial x^{2}} + \frac{\partial^{2}}{\partial y^{2}}$ , the rate of disappearance through the cell membrane having been accounted for by the second term.

Substituting (2.12) into (2.13), and introducing dimensionless variables $\hat{I} = \frac{I}{[G_{i}] T}, \hat{S} = \frac{S}{[G_{s}] T},$ , where [G_i]_T and [G_s]_T are the total concentrations per cell of G_αi and G_αs, respectively, in μmol m^–2, we arrive at our model consisting of the following 2 equations

\frac{\partial \hat{I}}{\partial τ} = - a_{1} \hat{I} + a_{2} \hat{S}

(2.14)

\begin{matrix} \frac{\partial \hat{S}}{\partial τ} = - a_{3} \hat{S} - \frac{b_{1} \hat{S}}{b_{2} + \hat{S}} + \frac{a_{4} \hat{S}}{{(a_{5} \hat{S} + \hat{I})}^{2}} \\ + a_{6} + μ \nabla^{2} \hat{S} \end{matrix}

(2.15)

where

\begin{array}{l} τ = \frac{t}{t_{b}}, a_{1} = t_{b} {[G_{i}]}_{T} k_{- i}, a_{2} = t_{b} {[G_{S}]}_{T} k_{i}, \\ a_{3} = t_{b} {[G_{S}]}_{T} {\tilde{a}}_{3}, a_{4} = \frac{t_{b} {[G_{S}]}_{T}^{2}}{{[G_{i}]}_{T}^{2}} {\tilde{b}}_{3}, a_{5} = \frac{{\tilde{b}}_{1} {[G_{S}]}_{T}}{{[G_{i}]}_{T}}, a_{6} = t_{h} k_{S} K_{c}, \\ b_{1} = t_{b} {\tilde{b}}_{3}, b_{2} = \frac{{\tilde{b}}_{4}}{{[G_{i}]}_{T}}, μ = \frac{t_{b} \tilde{μ} {[G_{S}]}_{T}}{d^{2}}, \end{array}

, with t_b being the characteristic ligand binding time, and d the inner plus outer membrane thickness. Typical values of t_b and d are 0.01 s (Erban and Othmer, 2005) and 20 nm (Institute for Biomolecular Design, http://redpoll.pharmacyualberta.ca).

In the next section, we shall carry out a weakly non-linear stability analysis on (2.14)-(2.15) in order to show the existence of rhombic and hexagonal planform solutions to our model following the technique discussed by Wollkind et al.(1994) and reviewed by Stephenson and Wollkind (1995).

Nonlinear Stability Analysis

In order to apply the technique of weakly nonlinear stability theory, we let

F (\hat{I}, \hat{S}) = - a_{1} \hat{I} + a_{2} \hat{S}

(3.1)

\begin{matrix} G (\hat{I}, \hat{S}) = - a_{3} \hat{S} - \frac{b_{1} \hat{S}}{b_{2} + \hat{S}} \\ + \frac{a_{4} {\hat{S}}^{2}}{{(a_{5} \hat{S} + \hat{I})}^{2}} + a_{6} \end{matrix}

(3.2)

which transforms (2.14)-(2.15) into

{\hat{I}}_{t} = F (\hat{I}, \hat{S})

(3.3)

{\hat{S}}_{t} = G (\hat{I}, \hat{S}) + μ \nabla^{2} \hat{S}

(3.4)

By expanding F(Î,Ŝ) and G(Î, Ŝ) into Taylor's series about the steady state (I₀,S₀) of (3.3) and (3.4) and letting the perturbations

i = \hat{I} - I_{0},

, and

s = \hat{S} - S_{0},

, we obtain the following system:

\begin{matrix} {(\begin{matrix} \begin{matrix} i \\ s \end{matrix} \end{matrix})}_{t} = (\begin{matrix} f_{1} & f_{2} \\ g_{1} & g_{2} \end{matrix}) (\begin{matrix} i \\ s \end{matrix}) + (\begin{matrix} o & o \\ g_{3} & g_{4} \end{matrix}) (\begin{matrix} i^{2} \\ s^{2} \end{matrix}) \\ + (\begin{matrix} o \\ g_{5} \end{matrix}) (i s) + (\begin{matrix} o \\ μ \nabla^{2} \hat{S} \end{matrix}) \end{matrix}

(3.5)

where

\begin{matrix} f_{1} \equiv - a_{1}, f_{2} \equiv a_{2}, g_{1} \equiv - \frac{2 a_{4} S_{0}^{2}}{{(a_{5} S_{0} + I_{0})}^{3}}, \\ g_{2} \equiv - a_{3} - \frac{b_{1} b_{2}}{{(b_{2} + S_{0})}^{2}} + \frac{2 a_{4} I_{0} S_{0}}{{(a_{5} S_{0} + I_{0})}^{3}}, \\ g_{3} \equiv \frac{3 a_{4} S_{0}^{2}}{{(a_{5} S_{0} + I_{0})}^{4}}, \\ g_{4} \equiv \frac{b_{1} b_{2}}{{(b_{2} + S_{0})}^{3}} + \frac{a_{4} I_{0} (- 2 a_{5} S_{0} + I_{0})}{{(a_{5} S_{0} + I_{0})}^{4}} and \\ g_{5} \equiv - \frac{2 a_{4} (- a_{5} S_{0}^{2} + 2 S_{0} I_{0})}{{(a_{5} S_{0} + I_{0})}^{4}} \end{matrix}

(3.6)

A Rhombic Planform Analysis

In order to investigate the possibility of occurrence in our model of rhombic-type patterns, we shall consider a rhombic planform solution of (3.5) of the form (referring the readers to Stephenson and Wollkind, (1995) for the motivation for the choices of functions therein)

\begin{matrix} \tilde{v} (x, y, t) \sim A_{1} (t) \cos (q_{c} x) (\begin{matrix} 1 \\ 1 \end{matrix}) \\ + B_{1} (t) \cos (q_{c} z) (\begin{matrix} 1 \\ 1 \end{matrix}) \\ + A_{1}^{2} (t) [{\tilde{v}}_{2000} + {\tilde{v}}_{2020} \cos (2 q_{c} x)] \\ + A_{1} (t) B_{1} (t) [{\tilde{v}}_{1111} \cos (q_{c} (x + z)) \\ + {\tilde{v}}_{111 (- 1)} \cos (q_{c} (x - z))] \\ + B_{1}^{2} (t) [{\tilde{v}}_{0200} + {\tilde{v}}_{0202} \cos (2 q_{c} z)] \\ + A_{1}^{3} (t) [{\tilde{v}}_{3010} \cos (q_{c} x) \\ + {\tilde{v}}_{3030} \cos (3 q_{c} x)] \\ + A_{1}^{2} (t) B (t) [{\tilde{v}}_{2101} \cos (q_{c} z) \\ + {\tilde{v}}_{2121} \cos (q_{c} (2 x + z)) \\ + {\tilde{v}}_{212 (- 1)} \cos (q_{c} (2 x - z))] \\ + A_{1} (t) B^{2} (t) [{\tilde{v}}_{1210} \cos (q_{c} x) \\ + {\tilde{v}}_{1212} \cos (q_{c} (x + 2 z)) \\ {\tilde{v}}_{121 (- 2)} \cos (q_{c} (x - 2 z))] \\ + B_{1}^{3} (t) [{\tilde{v}}_{0301} \cos (q_{c} z) \\ + {\tilde{v}}_{0303} \cos (3 q_{c} z)] \end{matrix}

(3.7)

where

\begin{array}{l} \tilde{v} (x, y, t) = (\begin{matrix} i (x, y, t) \\ s (x, y, t) \end{matrix}), {\tilde{v}}_{j l m n} = (\begin{matrix} i_{j l m n} \\ s_{j l m n} \end{matrix}) and z = x \\ \cos (ϕ) + y \sin (ϕ) \end{array}

with the amplitude equations:

\begin{matrix} \frac{d A_{1} (t)}{d t} \sim σ A_{1} (t) - A_{1} (t) [α_{1} A_{1}^{2} (t) \\ + β_{1} B_{1}^{2} (t)] \end{matrix}

(3.8)

\begin{matrix} \frac{d B_{1} (t)}{d t} \sim σ B_{1} (t) - B_{1} (t) [β_{1} A_{1}^{2} (t) \\ + α_{1} B_{1}^{2} (t)] \end{matrix}

(3.9)

as the zero^th order system of the most dangerous mode (Stephenson and Wollkind, 1995), α₁ and β₁ being the Landau constants to be determined later, as well as σ.

In (3.7), we are using the notation $\tilde{ν}$ _jlmn for the coefficient of each term in (3.7) of the form $A_{1}^{j} (t) B_{1}^{l} (t) \cos (q_{c} (m x + n z))$ . Please see the work of Stephenson and Wollkind, (1995) for more detail of the technique. On substituting this solution (3.7) into (3.5), we obtain a sequence of vector systems, each of which corresponds to one of these terms. We now catalogue the solutions for the first-order system. In particular, the first order system which corresponds to j = m = 1, l = n = 0 is

(\begin{matrix} σ \\ σ \end{matrix}) = (\begin{matrix} f_{1} & f_{2} \\ g_{1} & g_{2} \end{matrix}) (\begin{matrix} 1 \\ 1 \end{matrix}) + (\begin{matrix} 0 \\ - μ q_{c}^{2} \end{matrix})

(3.10)

Hence,

\begin{array}{l} σ = g_{1} + g_{2} - μ q_{c}^{2} \equiv σ_{0}, and μ = \frac{1}{q_{c}^{2}} (g_{1} + \\ g_{2} - f_{1} - f_{2}) \end{array}

The solutions for the four second-order systems, which can be solved in a straight forward manner, are as follows:

\begin{matrix} i_{2000} = \frac{g_{3} + g_{4}}{2 (- g_{1} + \frac{g_{2} f_{1}}{f_{2}})}, \\ s_{2000} = - \frac{f_{1}}{f_{2}} i_{2000}, \\ i_{2020} = \frac{g_{3} + g_{4}}{2 (- g_{1} - \frac{f_{1}}{f_{2}} (- g_{2} + 4_{μ}))}, \\ s_{2020} = \frac{(2 σ - f_{1})}{f_{2}} i_{2020}, \\ i_{1111} = \frac{(g_{3} + g_{4})}{- g_{1} - \frac{f_{1}}{f_{2}} [- g_{2} + μ (1 + 2 \cos (ϕ))]}, \\ s_{1111} = - \frac{f_{1}}{f_{2}} i_{1111}, \\ i_{111 (- 1)} = \frac{(g_{3} + g_{4})}{- g_{1} - \frac{f_{1}}{f_{2}} [- g_{2} + μ (1 - 2 \cos (ϕ))]}, a n d \\ s_{111 (- 1)} = - \frac{f_{1}}{f_{2}} i_{111 (- 1)} \end{matrix}

From the third-order systems, we can find the Landau constants α₁ and β₁ as

\begin{matrix} α_{1} = - \frac{f_{2}}{g_{1} + f_{2}} [(2 g_{3} + g_{5}) i_{2000} \\ + (2 g_{4} + g_{5}) s_{2000} + g_{3} i_{2020}] \end{matrix}

(3.11)

and

\begin{matrix} β_{1} = - \frac{f_{2}}{g_{1} + f_{2}} [(2 g_{3} + g_{5}) i_{2000} \\ + (2 g_{4} + g_{5}) s_{20000} \\ + g_{3} (i_{1111} + i_{111 (- 1)}) \\ + g_{4} (s_{1111} + s_{111 (- 1)})] \end{matrix}

(3.12)

leaving the derivation to the Appendix.

Having developed these formulae for Landau constants, we now turn our attention to the rhombic planform amplitude equations (3.8)–(3.9) which possess the following equivalence classes of critical points $(A_{1}^{0}, B_{1}^{0})$ when $q_{c}^{2} = 1$ :

I : A_{1}^{0} = B_{1}^{0} = 0

(3.13)

I I : {(A_{1}^{0})}^{2} = \frac{σ}{α_{1}}, B_{1}^{0} = 0

(3.14)

I I I : A_{1}^{0} = B_{1}^{0} w i t h {(A_{1}^{0})}^{2} = \frac{σ}{α_{1} + β_{1}}

(3.15)

By assuming that α₁ > 0 and α₁ + β₁ > 0, we investigate the stability of the critical points in (3.13)-(3.15) by seeking a solution of our amplitude equations (3.8)–(3.9) of the form

A_{1} (t) = A_{1}^{0} + \in c_{0} \exp (p t) + O (\in^{2})

(3.16)

B_{1} (t) = B_{1}^{0} + \in c_{0} \exp (p t) + O (\in^{2})

(3.17)

with

| \in | < 1

. One finds the following associated roots for p.

I : p_{1, 2} = σ,

(3.18)

I I : p_{1} = - 2 σ, p_{2} = (1 - \frac{β_{1}}{α_{1}}) σ,

(3.19)

I I I : p_{1} = - 2 σ, p_{2} = \frac{2 (β_{1} - α_{1}) σ}{α_{1} + β_{1}},

(3.20)

which yield the stability criteria for each critical point

(A_{1}^{0}, B_{1}^{0})

that I is stable for σ < 0, II is stable for σ > 0, α₁ >0 and

(1 - \frac{β_{1}}{α_{1}}) < 0

, and III is stable for σ > 0, α₁ + β₁ > 0 and

(\frac{β_{1} - α_{1}}{α_{1} + β_{1}}) < 0

. According to Stephenson and Wollkind (1995), I and II represent the homogeneous and striped states, respectively, while III can be identified with a rhombic pattern possessing characteristic angle φ.

We now investigate the critical points II and III, when stable, of our amplitude equations in relation to Turing patterns of interest. To the lowest order, the solution of the model associated with these critical points is given by the deviation

s \sim A_{1}^{0} \cos (\frac{2 π x}{λ_{c}}) + B_{1}^{0} \cos (\frac{2 π z}{λ_{c}})

(3.21)

where z = xcos(φ) + ysin(φ) and

λ_{c} = \frac{2 π}{q_{c}}

. The contour plot of this deviation function with

A_{1}^{0} > 0

and

B_{1}^{0} = 0

relevant to critical point II is shown in the (x, y) plane in Figure 2. Clearly, such alternating light and dark parallel bands produced by the critical point II should be identified with a striped Turing pattern as anticipated above. In order to make an analogous interpretation of critical point III, we consider the deviation function (3.21) with

A_{1}^{0} = B_{1}^{0} > 0

. We generate the contour plot of s in Figure 3, which exhibits the rhombic pattern.

Figure 2

Contour plot of s for the critical point II in (3.14) with a₁ = 0.6, a₂ = 0.8, a₃ = 0.5, a₄ = 0.4, a₅ = 0.005, a₆ = 0.001, b₁ = 0.0002, b₂ = 0.55, μ = .00262, and $ϕ = \frac{π}{2}$ .

Figure 3

Contour plot of s for the critical point III in (3.15) with a₁ = 0.1, a₂ = 0.2, a₃ = 0.5, a₄ = 0.4, a₅ = 0.05, a₆ = 0.001, b₁ = 0.001, b₂ = 1.5, μ = 0.0959972, and $ϕ = \frac{π}{3}$ .

A Hexagonal Planform Analysis

In order to investigate the possibility of occurrence in our model of hexagonal-type patterns, we shall consider a hexagonal planform solution of (3.5) of the form (Stephenson and Wollkind, 1995)

\begin{matrix} \tilde{v} (x, y, t) \sim A_{2} (t) \cos (q_{c} x) (\begin{matrix} 1 \\ 1 \end{matrix}) \\ + B_{2} (t) \cos (\frac{1}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y) (\begin{matrix} 1 \\ 1 \end{matrix}) \\ + A_{2}^{2} (t) [{\tilde{v}}_{2000} + {\tilde{v}}_{2040} \cos (2 q_{c} x)] \\ + A_{2} (t) B_{2} (t) \\ [{\tilde{v}}_{1111} \cos (\frac{1}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y) \\ + {\tilde{v}}_{1131} \cos (\frac{3}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y)] \\ + B_{2}^{2} (t) [{\tilde{v}}_{0200} + {\tilde{v}}_{0220} \cos (q_{c} x) \\ + {\tilde{v}}_{0202} \cos (\sqrt{3} q_{c} y) \\ + {\tilde{v}}_{0222} \cos (q_{c} x) \cos (\sqrt{3} q_{c} y)] \\ + A_{2}^{3} (t) [{\tilde{v}}_{3020} \cos (q_{c} x) + {\tilde{v}}_{3060} \cos (3 q_{c} x)] \\ + A_{2}^{2} (t) B_{2} (t) \\ [{\tilde{v}}_{2111} \cos (\frac{1}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y) \\ + {\tilde{v}}_{2131} \cos (\frac{3}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y) \\ + {\tilde{v}}_{2151} \cos (\frac{5}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y)] \\ + A_{2} (t) B_{2}^{2} (t) [{\tilde{v}}_{1200} + {\tilde{v}}_{1220} \cos (q_{c} x) \\ + {\tilde{v}}_{1240} \cos (2 q_{c} x) + {\tilde{v}}_{1202} \cos (\sqrt{3} q_{c} y) \\ + {\tilde{v}}_{1222} \cos (q_{c} x) \cos (\sqrt{3} q_{c} y) \\ {\tilde{v}}_{1242} \cos (2 q_{c} x) \cos (\sqrt{3} q_{c} y)] \\ + B_{2}^{3} (t) [{\tilde{v}}_{0311} \cos (\frac{1}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y) \\ {\tilde{v}}_{0331} \cos (\frac{3}{2} q_{c} x) \cos (\frac{\sqrt{3}}{2} q_{c} y) \\ {\tilde{v}}_{0313} \cos (\frac{1}{2} q_{c} x) \cos (\frac{3 \sqrt{3}}{2} q_{c} y) \\ {\tilde{v}}_{0333} \cos (\frac{3}{2} q_{c} x) \cos (\frac{3 \sqrt{3}}{2} q_{c} y)] \end{matrix}

(3.22)

with the amplitude equations

\begin{matrix} \frac{d A_{2} (t)}{d t} \sim σ A_{2} (t) - α_{0} B_{2}^{2} (t) \\ - A_{2} (t) [α_{1} A_{2}^{2} (t) \\ + α_{2} B_{2}^{2} (t)] \end{matrix}

(3.23)

\begin{matrix} \frac{d B_{2} (t)}{d t} \sim σ B_{2} (t) - 4 α_{0} A_{2} (t) B_{2} (t) \\ - B_{2} (t) [2 α_{2} A_{2}^{2} (t) \\ \frac{1}{4} (α_{1} + 2 α_{2}) B_{2}^{2} (t)] \end{matrix}

(3.24)

as the zero^th order system of the most dangerous mode (Stephenson and Wollkind, 1995).

In (3.22), we are employing the notation ${\tilde{v}}_{j l m n}$ for the coefficient of each term in (3.22) of the form $A_{2}^{j} (t) B_{2}^{1} (t) \cos (\frac{m q_{c} x}{2}) \cos (\frac{n \sqrt{3} q_{c} y}{2})$ . Substituting this solution into (3.5) and proceeding in exactly the same manner as in the rhombic planform analysis, we again obtain a sequence of vector systems, each of which corresponds to one of these terms. In particular, the first order system which corresponds to j = m = 1, l = n = 0 is

(\begin{matrix} σ \\ σ \end{matrix}) = (\begin{matrix} f_{1} & f_{2} \\ g_{1} & g_{2} \end{matrix}) (\begin{matrix} 1 \\ 1 \end{matrix}) + (\begin{matrix} 0 \\ - μ q_{c}^{2} \end{matrix})

(3.25)

Hence,

σ = g_{1} + g_{2} - μ q_{c}^{2} = σ_{0},

and

\begin{array}{l} μ = \frac{1}{q_{c}^{2}} (g_{1} + \\ g_{2} - f_{1} - f_{2}) \end{array}

There are eight second-order systems, two of which contain the Landau constant α₀. We can write one of them as

M_{1} {\tilde{v}}_{0220} = α_{0} (\begin{matrix} 1 \\ 1 \end{matrix}) + (\begin{matrix} 0 \\ \frac{1}{4} (g_{3} + g_{4} + g_{5}) \end{matrix})

(3.26)

where

\begin{matrix} M_{1} \equiv (\begin{matrix} 2 σ_{0} - f_{1} & - f_{2} \\ - g_{1} & 2 σ_{0} - g_{2} + μ q_{c}^{2} \end{matrix}) \\ = (\begin{matrix} σ_{0} + f_{2} & - f_{2} \\ - g_{1} & σ_{0} + g_{1} \end{matrix}) \end{matrix}

(3.27)

Considering the adjoint linear eigenvalue problem of (3.27):

M_{1}^{T} {\tilde{v}}_{1}^{*} = σ_{1}^{*} {\tilde{v}}_{1}^{*}

(3.28)

where

σ_{1}^{*} = σ_{0}

is an eigenvalue of M₁ and

M_{1}^{T}

, we obtain

{\tilde{v}}_{1}^{*} = (\begin{matrix} g_{1} \\ f_{2} \end{matrix})

. By taking inner products of (3.26) with

{\tilde{v}}_{1}^{*}

, we find, upon making use of the adjoint property, that

\begin{matrix} σ_{0} {\tilde{v}}_{0220} . {\tilde{v}}_{1}^{*} = α_{0} (\begin{matrix} 1 \\ 1 \end{matrix}) . {\tilde{v}}_{1}^{*} \\ + (\begin{matrix} 0 \\ \frac{1}{4} (g_{3} + g_{4} + g_{5}) \end{matrix}) . {\tilde{v}}_{1}^{*} \end{matrix}

(3.29)

Then, taking the limit as σ₀ → 0 we obtain

α_{0} = - \frac{f_{2}}{4 (g_{1} + f_{2})} (g_{3} + g_{4} + g_{5})

(3.30)

The other six second-order systems can be solved in a straight forward manner, yielding solutions as follows:

\begin{matrix} i_{2000} = - \frac{\frac{1}{2} g_{3} + \frac{1}{2} g_{4} + g_{5}}{g_{1} - \frac{g_{2} f_{1}}{f_{2}}}, s_{2000} = - \frac{f_{1}}{f_{2}} i_{2000}, \\ i_{2040} = - \frac{g_{5} + \frac{1}{2} g_{3} + \frac{1}{2} g_{4}}{g_{1} - \frac{f_{1}}{f_{2}} (g_{2} - 4 μ)}, s_{2040} = - \frac{f_{1}}{f_{2}} i_{2040}, \\ i_{0200} = - \frac{g_{3} + g_{4} + g_{5}}{4 (g_{1} - \frac{g_{2} f_{1}}{f_{2}})}, s_{0200} = - \frac{f_{1}}{f_{2}} i_{0200}, \\ i_{0220} = \frac{α_{0} (g_{2} - μ) - f_{2} (α_{0} + \frac{1}{4} g_{3} + \frac{1}{4} g_{4} + \frac{1}{4} g_{5})}{f_{2} g_{1} - f_{1} (g_{2} - μ)}, \\ s_{0220} = - \frac{f_{1} i_{0220} + α_{o}}{f_{2}}, \\ i_{1111} = \frac{4 α_{0} (g_{2} - μ) - f_{2} (4 α_{0} + g_{3} + g_{4} + g_{5})}{f_{2} g_{1} - f_{1} (g_{2} - μ)}, \\ s_{1111} = - \frac{f_{1} i_{1111} + 4 α_{0}}{f_{2}}, \\ i_{1131} = - \frac{f_{2} (g_{3} + g_{4} + g_{5})}{f_{2} g_{1} - f_{1} (g_{2} - 3 μ)}, s_{1131} = - \frac{f_{1}}{f_{2}} i_{1131} \end{matrix}

Using the third-order systems, we can find the other two Landau constants α₁ and α₂ as

\begin{matrix} α_{1} = - \frac{1}{g_{1} + f_{2}} [(2 f_{2} g_{3} + f_{2} g_{5}) i_{2000} \\ + f_{2} (2 g_{4} + g_{5}) s_{2000} \\ + (f_{2} g_{3} + \frac{f_{2} g_{5}}{2}) i_{2040}] \\ - \frac{1}{g_{1} + f_{2}} (f_{2} g_{4} + \frac{f_{2} g_{5}}{2}) s_{2040} \end{matrix}

(3.31)

and

\begin{matrix} α_{2} = - \frac{1}{g_{1} + f_{2}} [8 α_{0} g_{1} i_{0220} + 8 α_{0} f_{2} s_{0220} \\ + (2 f_{2} g_{3} + f_{2} g_{5}) i_{0200} \\ + f_{2} (2 g_{4} + g_{5}) s_{0200}] \\ - \frac{1}{g_{1} + f_{2}} [(\frac{1}{2} f_{2} g_{3} + \frac{1}{4} f_{2} g_{5}) i_{1111} \\ + (\frac{1}{2} f_{2} g_{4} + \frac{1}{4} f_{2} g_{5}) s_{1111} \\ + \frac{3}{8} f_{2} g_{3} i_{1131}] \\ - \frac{1}{g_{1} + f_{2}} (\frac{3}{8} f_{2} g_{4} s_{1131}) \end{matrix}

(3.32)

also leaving the detail of the derivation to the Appendix.

Having developed these formulae for Landau constants, we now turn our attention to the hexagonal planform amplitude equations (3.23)–(3.24) which possess the following equivalence classes of critical points $(A_{2}^{0}, B_{2}^{0})$ when $q_{c}^{2} = 1$ :

I : A_{2}^{0} = B_{2}^{0} = 0

(3.33)

I I : {(A_{2}^{0})}^{2} = \frac{σ}{α_{1}}, B_{2}^{0} = 0

(3.34)

\begin{matrix} I I I^{\pm} : A_{2}^{0} = B_{2}^{0} = {(A_{2}^{0})}^{\pm} \\ = \frac{- 2 α_{0} \pm {[4 α_{0}^{2} + (α_{1} + 4 α_{2}) σ]}^{\frac{1}{2}}}{α_{1} + 4 α_{2}} \end{matrix}

(3.35)

I V : A_{2}^{0} = - \frac{4 α_{0}}{2 α_{2} - α_{1}}, B_{2}^{0} = \frac{σ - σ}{α_{1} + 2 α_{2}}

(3.36)

where we assume that

\begin{array}{l} α_{1} > 0, α_{1} + 4 α_{2} > 0, σ -_{1} \\ \equiv - \frac{4 α_{0}^{2}}{α_{1} + 4 α_{2}}, σ_{1} \equiv \frac{16 α_{1} α_{0}^{2}}{{(2 α_{2} - α_{1})}^{2}} \end{array}

and

σ_{2} \equiv \frac{32 (α_{1} + α_{2}) α_{0}^{2}}{{(2 α_{2} - α_{1})}^{2}}

The orbital stability conditions for these critical points can be posed in terms of σ. This sort of stability of pattern formation is meant in the sense of a family of solutions in the plane which may interchange with each other but do not grow or decay to a solution type from a different family. Such an interpretation depends upon the translational and rotational symmetries inherent to the amplitude-phase equations. This invariance also limits each equivalence class of critical points to a single member that must be explicitly considered. This line of reasoning was used in our cataloguing of the critical point for the rhombic planform analysis in the previous section as well. Thus, critical point I is stable in this sense for σ < 0 while the stability behavior of II and III^±, which depends upon the signs of α₀ and 2α₂- α₁ as well, has been summarized in Table 1 under the further assumption that α₁ + α₂ > 0.

Table 1

Orbital stability behavior of critical points Hand III^±

α₀	2α₂-α₁	Stable Structures
+	–, 0	III^– for σ > σ_–1
+	+	III^– σ_–1 < σ < σ₂, II for σ > σ₁
0	–	III^± for σ > 0
0	+	II for σ > 0
–	+	III⁺ for σ_–1 < σ < σ₂, II for σ > σ₁
–	–, 0	III⁺ for σ > σ_–1

We next offer a morphological interpretation of the potentially stable critical points described above relative to the Turing patterns under investigation. To the lowest order, the solution of the model associated with these critical points is given by the deviation function

s \sim A_{2}^{0} \cos (x) + 2 B_{2}^{0} \cos (\frac{1}{2} x) \cos (\frac{\sqrt{3}}{2} y)

(3.37)

Since ${(A_{2}^{0})}^{+} > 0$ for α₀ ≤ 0 and ${(A_{2}^{0})}^{-} > 0$ for α₀ > 0, we can conclude that the contour plots of(3.37) would have circular elevations at the centers of the hexagonal for critical point III⁺ when stable, and circular depression for III^–. These contour plots are shown in Figure 4 and Figure 5 for III⁺ and III^–, respectively. As in analogous plots relevant to our rhombic planform analysis, elevations appear light and depressions dark in these figures. Hence, recalling that the Turing patterns under consideration are classified by their light regions, we identify hexagonal arrays of nets or honeycombs with critical point III^– and of spots or dots with critical point III⁺.

Figure 4

Contour plot of s for the critical point III⁺ in (3.35) with a₁ = 0.1, a₂ = 1.5, a₃ = 0.6, a₄ = 0.9, a₅ = 0.001, a₆ = 0.0001, b₁ = 0.5, b₂ = 2, μ = 0.841495.

Figure 5

Contour plot of s for the critical point III^– in (3.35) with a₁ = 0.5, a₂ = 1.5, a₃ = 0.6, a₄ = 0.5, a₅ = 0.00001, a₆ = 0.01, b₁ = 0.95, b₂ = 3, μ = 1.13812.

Discussion

Since the values of various kinetic constants are not yet known, we have chosen particular sets of parametric values for illustration purposes, following the example of Krishnan and Iglesias, (2003). For example, to plot the contours in Figure 5, we used the following dimensionless parameters: a₁ = 0.5, a₂ = 1.5, a₃ = 0.6, a₄ = 0.5, a₅ = 0.00001, a₆ = 0.01, b₁ = 0.95, b₂ = 3.0 and μ = 1.13812. This dimensionless diffusion coefficient μ corresponds to the physical diffusion constant of 2.726 × 10^–10 m²s^–1, while the cAMP diffusion constant was estimated to be 2.5 × 10^–10 m²s^–1 by Rapple et al.(2002). Our choice of the parametric values took into account several physical features. The activation and deactivation of the inhibitor are more or less of the same order of magnitude. The production rate of cAMP in the absence of the activated regulator in very low so that a₆ « 1. The amplification of the signal by the second messenger is comparable with the removal rate of the signaling hormone from the system.

Comparing the planforms in Figures 4 (quantum dots) and Figure 5 (honeycombs), we observe that quantum dots are expected to occur in the signaling system in which the maximum absorbtion b₁ is not so high. On the other hand, the honeycombs pattern is expected when the maximum rate of disappearance through the membrane of the signaling hormone is relatively high. The values used here are of a comparable order of magnitude as those of the corresponding parameters in the work of Krishnan and Iglesias, (2003). Stripes and rhombic patterns, however result when the rate of diffusion over the membrane surface μ and disappearance through the membrane are noticeably lower.

It has only recently been appreciated that G protein coupled receptors (GPCR) and their associated signaling components are not randomly dispersed throughout the plasmalemma. According to Ostrom et al.(2002), expression of βARs(β-Aderenergic receptors) and AC in several cultured cell models and cardiac muscle cells is enriched in distinct microdomains of the plasma membrane. In these microdomains, one can observe retention of proteins that have particular post-translational modifications (Shaul et al. 1996). Moreover, caveolin in caveolae contains a binding domain that interacts with certain signaling molecules and thereby facilitate localization of signaling molecules in caveolae (Ostrom et al. 2002). Such compartmentation of signaling molecules contradicts the belief that components of GPCR signal transduction are randomly distributed and have extensive mobility in the plasma membrane. On the contrary, these proteins appear to be restricted to plasmalemmal microdomains, probably conductive to rapid and specific signal transduction (Anderson, 1998; Okamoto et al. 1998; Ostrom et al. 2000; Ostrom et al. 2002).

In their recent study, Ostrom et al. (2002) tested the hypothesis that expression and localization of GPCRs and isoforms of AC might be critical determinants of how vascular smooth muscle cells respond to extracellular signals. They observed that localization of AC and the components that influence AC activity might help “tailor” the ability of cells to respond to extracellular and intracellular signals by defining a precise environment in which the second messenger (cAMP) will be generated. In addition, the spatial organization of signaling molecules was found to be a likely important factor in vascular smooth muscle cell regulation, in particular with respect to regulation of cAMP formation which leads to decreased contractile tone, decreased vascular resistance and decreased blood pressure (Ostrom et al. 2002).

Because primary cells in culture can rapidly lose their differentiated phenotype, Ostrom et al.(2002) examined RASMC (rat aortic smooth muscle cells) morphology using transmission electron microscopy. Ultrastructural examination of passage 5 RASMC indicated that these cells possess both morphologic caveolae (light vesicular structures) enriched with β-aderenergic receptors, facilitating localization of signaling molecules, which might be identified with stable III⁺ pattern (hexagonal arrays of dots), and morphologic features consistent with contractile phenotype (or filaments) which may facilitate formation of stable pattern of type II (stripes) identified theoretically in the text.

Another recent report by Kriebel et al.(2003) investigated the role of AC localization in regulating the streaming formation in chemotaxis. In their work, the authors presented experimental results concerning the morphological polarization of Dictyostilium cells with a cAMP gradient, and showed that local accumulation of AC at the uropod of a cell is required for stream formation. They proposed that the asymmetric distribution of AC provides a compartment from which cAMP may act to amplify the chemical gradient. Thus, this could represent a unique mechanism that cells use to amplify chemotactic responses.

To date, little has been published in terms of pictured spatial distribution of membrane receptors. Using radioactivity-labelled drugs and so-called affinity labels, Peter Gaudeng Waser (1983) has been able to produce binding with high strength or irreversible binding on specific receptors which permitted their isolation and characterization. His electronmicrograph with high resolution shows Acetylcholine receptors in membrane fragments (Waser, 1983) which exhibit hexagonal structure of localization. Other images that show striped or rhombic structures may be found in the work of Unwin et al. (1988) and on http:\\gingi.uchicago.edu/archr.html. Since signaling molecules are expected to colocalize with their receptors, such pictured receptor complex shown in Waser's report (1983) provides experimental evidence, to a certain extent, in support of our theoretical predictions that the system of interest may operate, under certain circumstances, in such weakly nonlinear regimes.

Conclusion

In this paper, we have concentrated on the signal transduction pathways that involve membrane-based receptors which couple to adenylate cyclase. G proteins play a key role in these pathways by activating adenylate cyclase to produce cAMP which functions as a second messenger. Examples of hormones that act via the utilization of a cAMP system are ACTH, GnRH, PTH, and Vasopressin, which are significantly implicated in carcinogenesis.

For some time, various stages of cancers have been classified by means of measurement of receptor levels and their hormone dependence from biopsy specimens. Unfortunately, the indications of measurements of steroid receptor levels are not always absolute and it appears that the ability of receptors to bind steroid does not ensure the complete functioning of the receptors on the cell. Clearly, more research into this area of steroid receptors dysfunction is needed to sharpen prognostic ability. The continued identification of the presence of receptors and their hormone dependence as correlated to their specific functions will greatly enhance prognosis and treatment of such disease.

The result of our mathematical modelling underlines the need for more detailed and precise microscopic study on the spatial distribution of the protein-receptor complexes and the plasmalemma morphology, in the hope of shed more light on the precise manner in which spatial organization of multiple component signal transduction cascades may provide a means to generate signals with high fidelity and efficiency. How the different spatial patterns or gradients described in this paper are related to specific cellular function or dysfunction also remains to be investigated and characterized. More specialized experimental studies are required to provide microscopic data to test hypotheses and establish linkages among the different dynamic activities in the signal transduction cascades. Our model provides a useful context in which to explore the relevance of Turing-type pattern formation to signaling pathways.

Footnotes

Acknowledgements

This work was supported by the National Center for Genetic Engineering and Biotechnology, and the Commission on Higher Education, Ministry of Education.

Appendix

References

Anderson

1998. The caveolae membrane system. Annu. Rev. Biochem., 67: 199–225.

Castets

, Dulos

, Boissonade

, and De Kepper

1990. Experimental evidence of a sustained standing Turing-type nonequilibrium chemical pattern. Phys. Rev. Lett., 64: 2953–2956.

Erban

, and Othmer

H.G.

2005. From signal transduction to spatial pattern formation in E.Coli: a paradigm for multiscale modeling in biology. Multiscale Model. Simul., 3(2): 362–394.

Glansdorff

, and Prigogine

1971. Thermodynamics theory of structure, stability and fluctuations. Wiley, London.

Iglesias

P.A.

2003. Feedback control in intracellular signaling pathways: regulating chemotaxis in Dictyostelium discoideum. Europ. J. Contr., 9: 216–225.

Koch

A.J.

, and Meinhardt

1994. Biological pattern formation. Rev. Mod. Phys., 66: 1481–1507.

Kondo

, and Asai

1995. A reaction-diffusion wave on the skin of the marine angelfish pomacanthus. Nature, 376: 765–768.

Kriebel

P.W.

, Barr

V.A.

, and Parent

C.A.

2003. Adenylyl cyclase localization regulates sheaming during chemotaxis. Cell, 12: 549–560.

Krishnan

, and Iglesias

P.A.

2003. Analysis of the signal transduction properties of a module of spatial sensing in eukaryotic chemotaxis. Bull. Math. Biol., 65: 95–128.

10.

Levchenko

, and Iglesias

2002. Models of eukaryotic gradient sensing: application to chemotaxis of amoebae and neutrophils. Biophys. J., 82: 50–63.

11.

Meinhardt

1982. Models of biological pattern formation. Academic Press, London.

12.

Norman

A.W.

, and Litwack

1997. Hormones. Academic Press, California, U.S.A.

13.

Neubig

R.R.

, and Siderovski

D.P.

2002. Regulators of G-protein signalling as new central nervous system drug targets. Nature Rev., 1: 187–197.

14.

Okamoto

, Schlegel

, Scherer

P.E.

, and Lisanti

M.P.

1998. Caveolins, a family of scaffolding proteins for organizing “pre-assembled signaling complexes” at the plasma membrane. J. Biol. Chem., 273: 5419–5422.

15.

Ostrom

R.S.

, Liu

, Head

B.P.

, Gregorian

, Seasholtz

T.M.

, and Insel

P.A.

2002. Localization of adenylyl cyclase isoforms and G protein-coupled receptors in vascular smooth muscle cells: expression in caveolin-rich and noncaveolin domains. Mol. Pharmacol., 62: 983–992.

16.

Ostrom

R.S.

, Post

S.R.

, and Insel

P.A.

2000. Stoichiometry and compartmentation in G protein-coupled receptor signaling: implications for therapeutic interventions involving Gs. J. Pharmacol. Exp. Ther., 294: 407–412.

17.

Pansuwan

, Rattanakul

, Lenbury

, Wollkind

D.J.

, Harrison

, and Rajapakse

2005. Nonlinear stability analyses of pattern formation on solid surfaces during ion-sputtered erosion. Math. Comp. Mod., 41: 939–964.

18.

Rappel

W.J.

, Thomas

P.J.

, Levine

, and Loomis

W.F.

2002. Establishing direction during chemotaxis in eukaryotic cells. Biophys. J., 83: 1361–1367.

19.

Rauch

E.M.

, and Millonas

M.M.

2004. The role of trans-membrane signal transduction in turing-type cellular pattern formation. J. Theor. Biol., 226: 401–407.

20.

Sawai

, Maeda

, and Sawada

2000. Spontaneous symmetry breaking Turing-type pattern formation in a confined dictyostelium cell mass. Phys. Rev. Lett., 85: 2212–2215.

21.

Shaul

P.W.

, Smart

E.J.

, Robinson

L.J.

, German

, Yuhanna

I.S.

, Ying

, Anderson

R.G.

, and Michel

1996. Acylation targets endothelial nitric-oxide synthase to plasmalemmal caveolae. J. Biol. Chem., 271: 6518–6522.

22.

Spiegel

A.M.

2000. G protein defects in signal transduction. Horm. Res., 53(3): 17–22.

23.

Stephenson

L.E.

, and Wollkind

D.J.

1995. Weakly nonlinear stability analyses of one-dimensional Turing pattern formation in activator-inhibitor/immobilizer model systems. J. Math. Biol., 33: 771–815.

24.

Tang

, and Othmer

H.G.

1994. A G protein-based model of adaptation in Dictyostellium discoideum. Mathematical Biosciences, 120: 25–76.

25.

Tang

, and Othmer

H.G.

1995. Excitation, oscillations and wave propagation in a G-proein-based model of signal transduction in Dictyostellium discoideum. Phil. Trans. R. Soc. Lond. Ser., B. 349: 179–195.

26.

Turing

A.M.

1952. The chemical basis of morphogenesis. Philos. Trans. R. Soc. London, Ser., B 237: 37–72.

27.

Unwin

, Toyoshima

, and Kubalek

1988. Arrangement of the acetycholine receptor subunits in the resting and desensitized states, determined by cryoelectron microscopy of crystallized Torpedo postsynaptic membranes. J. Cell Bio., 107: 1123–1138.

28.

Waser

P.G.

1983. Cholinerge pharmakon-rezeptoren. Neujahrblatt. 83S, 50Abb., 185.

29.

Wollkind

D.J.

, Manoranjan

V.S.

, and Zhang

1994. Weakly nonlinear stability analyses of prototype reaction-diffusion model equations. Siam Review., 36(2): 176–214.