FRET imaging and statistical signal processing reveal positive and negative feedback loops regulating the morphology of randomly migrating HT-1080 cells

Cell migration is an essential process in a wide variety of biological phenomena, such as embryonic development, wound repair, immune surveillance, tumor cell invasion and metastasis in mammalian cells (Lauffenburger and Horwitz, 1996). Extracellular cues elicit the various intracellular responses in the organization of both the actin and the microtubule cytoskeletons. It is widely accepted that the major driving force of migration is the extension of a leading edge protrusion or lamellipodium, cell body contraction, and detachment of adhesions at the cell rear (Rajakulendran et al., 2009). All these steps involve the assembly, disassembly and reorganization of the actin cytoskeleton, which are coordinated both in space and time (Mitchison and Cramer, 1996; Rajakulendran et al., 2009).

The molecular mechanisms of cell migration have been extensively studied in tissue culture cells. Rho GTPases and lipid kinases have been identified as key regulators of cell migration. The best-characterized function of Rho GTPases is in the regulation of actin dynamics: RhoA increases the actomyosin contractility, whereas Rac1 and Cdc42 induce the polymerization of actin to form lamellipodial and filopodial protrusions, respectively (Ridley, 2001). Phosphatidyl inositols, which are composed of a D-myo-inositol-1-phosphate linked through its phosphate group to diacylglycerol, are phosphorylated and dephosphorylated by phosphatidyl inositol kinases and phosphatases, generating distinct phosphoinositide species (Saarikangas et al., 2010). Among them, phosphatidylinositol 3,4,5-trisphosphate [PtdIns(3,4,5)P₃] and phosphatidylinositol bisphosphates [PtdIns(3,4)P₂ and PtdIns(4,5)P₂] regulate reorganization of the actin cytoskeleton by providing a membrane-binding platform for a variety of proteins, including guanine nucleotide exchange factors (GEFs) for Rho GTPases and actin filament nucleating proteins (Saarikangas et al., 2010). With the advent of a myriad of fluorescent proteins, genetically encoded biosensors have been increasingly used to visualize the activities of Rho GTPases and phosphatidylinositols (Aoki et al., 2008; Miyawaki, 2003). In particular, live-cell imaging with biosensors based on the principle of Förster (or fluorescence) resonance energy transfer (FRET) have uncovered the spatiotemporal dynamics of these molecules in migrating cells (Itoh et al., 2002; Kurokawa and Matsuda, 2005; Nishioka et al., 2008).

Morphological polarization (i.e. a clear distinction between cell front and rear) is required for cell migration to enable cells to turn intracellularly generated forces into net cell body translocations (Lauffenburger and Horwitz, 1996). It has been suggested that intrinsic mechanisms of intracellular signaling systems with their kinetic fluctuations or noise suffice to activate this front–rear polarization during the locomotion of eukaryotic motile cells. For example, spatially homogenous stimulation with a chemoattractant, at least for neutrophils, induces a front–rear polarity with a change in filamentous actin (F-actin) distribution from azimuthal symmetry around the cell rim to concentration within a particular region (Coates et al., 1992; Lauffenburger and Horwitz, 1996). In addition, even in the absence of concentration gradients of external stimuli, ordered patterns of cell shape are observed in Dictyostelium, flies and mammalian cells, such as directional migration, oscillation behavior of membranes and laterally traveling waves along the cell boundary (Döbereiner et al., 2006; Machacek and Danuser, 2006; Maeda et al., 2008; Weiner et al., 2007). The Turing's reaction–diffusion model is one of the most supported models for the formation of ordered patterns in migrating cells (Gierer and Meinhardt, 1972; Kondo and Miura, 2010; Otsuji et al., 2010). On the basis of a framework in the reaction–diffusion model, self-enhancing positive and negative feedback regulations play a pivotal role in the autonomous emergence of an ordered pattern. It is suggested that PtdIns(3,4,5)P₃, Rac1 and the actin cytoskeleton are part of a positive feedback loop that can initiate polarization (Peyrollier et al., 2000; Srinivasan et al., 2003; Wang et al., 2002; Weiner et al., 2002). However, contrary to these models, rapid activation of endogenous Rac1 proteins triggered effective actin polymerization but failed to activate phosphoinositide 3-kinase (PI3K) to generate PtdIns(3,4,5)P₃ or induce cell polarization (Inoue and Meyer, 2008). Furthermore, the molecular mechanisms of such negative feedback loops remain largely unknown. We have attempted to identify molecular mechanisms of positive and negative feedback loops underlying the random migration of mammalian cells using FRET imaging with statistical signal processing.

We aimed to elucidate the molecular networks of Rho GTPases and phosphatidylinositols in randomly migrating mammalian cells by time-lapse imaging. We used human HT-1080 fibrosarcoma cells, which exhibit a high motility and invasiveness, as a model system of spontaneous random migration (Hoffman, 2005; Itoh et al., 2002). All experiments presented below employed time-lapse imaging of HT-1080 cells 1 hour after seeding on collagen-coated glass-bottomed dishes to induce spontaneous migration. The HT-1080 cells expressing biosensors were imaged every 1 minute for 2 hours in order to minimize phototoxicity and photobleaching.

We first examined the activity of the Rho GTPases Rac1, Cdc42 and RhoA in migrating HT-1080 cells using the FRET biosensors Raichu–Rac1, Raichu–Cdc42 and Raichu–RhoA, respectively. As reported previously (Itoh et al., 2002; Kraynov et al., 2000; Nalbant et al., 2004), Rac1 and Cdc42 were locally activated at lamellipodia and membrane ruffles at the leading edge, and inactivated at the retracting tail of migrating cells (Fig. 1A,B; supplementary material Movies 1, 2). By contrast, activation of RhoA was observed at both the leading edge and retracting tail (Fig. 1C; supplementary material Movie 3). This result was also consistent with the previous findings (Kurokawa and Matsuda, 2005; Pertz et al., 2006).

The spatiotemporal localization of phosphatidylinositols was also visualized with the FRET biosensors (Aoki et al., 2005; Nishioka et al., 2008) Pippi–PtdIns(3,4,5)P₃, Pippi–PtdIns(3,4)P₂ and Pippi–PtdIns(4,5)P₂. PtdIns(3,4,5)P₃ accumulated at the leading edge of migrating HT-1080 cells (Fig. 1D; supplementary material Movie 4), supporting the previous report that PtdIns(3,4,5)P₃ recruits and activates guanine nucleotide exchange factor(s) for Rac1 and Cdc42 (Aoki et al., 2005; Schmidt and Hall, 2002). PtdIns(3,4)P₂ and PtdIns(4,5)P₂ were also found to accumulate at the front of migrating HT-1080 cells (Fig. 1E,F; supplementary material Movies 5, 6). The increase in PtdIns(3,4)P₂ and PtdIns(4,5)P₂ at the lamellipodia of migrating cells could be attributed to the rapid turnover of these lipids (Nishioka et al., 2008).

We also examined the accumulation of F-actin in HT-1080 cells with Lifeact–EGFP, a marker of F-actin (Riedl et al., 2008). To correctly evaluate the accumulation of F-actin at membrane ruffles and pseudopodia (Dewitt et al., 2009), DsRed-express fused with a nuclear export sequence (ERed–NES) was coexpressed with Lifeact–EGFP, and used as a reference of the volume of cytoplasm. There was a higher ratio of Lifeact–EGFP signal to ERed–NES signal at the leading edge of migrating cells, clearly indicating the actin polymerization (Fig. 1G; supplementary material Movie 7). Notably, a negative control of the FRET biosensor, which included a Cdc42- and Rac1-interacting binding (CRIB) domain of PAK1 and RhoA, did not show any spatial or temporal change of FRET signal in HT-1080 cells (Fig. 1H).

Of note, lateral propagation of membrane protrusions and retractions were frequently observed around the cell periphery of migrating HT-1080 cells (Fig. 1, arrows), suggesting that common mechanisms underlie the ordered pattern of membrane dynamics, as observed in the migration of other eukaryotic cells (Döbereiner et al., 2006; Machacek and Danuser, 2006; Maeda et al., 2008).

To quantitatively evaluate the morphological changes, the edge velocity of membrane displacement was quantified at the whole cell periphery by a modified mechanical model (Machacek and Danuser, 2006). Briefly, the distance between a window on the cell periphery at time T and that at time T+1 was measured using time-lapse images (Fig. 2A,B; supplementary material Fig. S1). Then, the edge velocity was calculated and represented in a heat map as a function of the time and the sampling window around the cell periphery (see Materials and Methods for details; Fig. 2C). Although there was apparent heterogeneity of the pattern in the edge velocity map between different migrating HT-1080 cells, we often recognized a periodic pattern, as shown in Fig. 2C. In almost all cases, migrating cells exhibited one protrusion and one retracting tail, indicating the acquisition of front–back polarity even in the absence of a concentration gradient of external stimuli.

We next classified the migration pattern of HT-1080 cells by applying auto-correlation analysis to the edge velocity map (Döbereiner et al., 2006; Maeda et al., 2008). The auto-correlation function (ACF) in Fig. 2D has a wave-like pattern, demonstrating the occurrence of the lateral membrane wave along the cell boundary. With the temporal ACF (ΔWindow no.=0, dashed line in Fig. 2D), three parameters, Δt₀, Δt_min and Δt_max, were defined as shown in Fig. 2E. The pattern of the ACF map and threshold criteria of Δt₀ value (Δt₀>30 min) were utilized to classify the pattern of morphological dynamics into four categories: wave-like, oscillatory, directional, and non-classifiable mode (Maeda et al., 2008) (see Materials and Methods for more details; Fig. 2F; supplementary material Movie 8). Almost half of HT-1080 cells demonstrated the wave-like mode of migration (Table 1). Furthermore, the velocity of the lateral membrane wave, V_L, was calculated as 4.2±2.1 μm/minute by dividing the length of the cell periphery by Δt_max. These results suggested the existence of an intrinsic signaling network that established an ordered pattern of morphological dynamics without external stimuli. Notably, the expression of FRET biosensors per se did not substantially perturb the mode of cell migration in HT-1080 cells (supplementary material Fig. S2). Therefore, in the following study, we examined a molecular mechanism for generating wave-like motion in randomly migrating HT-1080 cells.

To identify a hierarchical connection among Rho GTPases, phosphatidyl inositols and morphological change, we used a cross-correlation analysis and characterized their temporal correlation (Machacek et al., 2009; Tsukada et al., 2008). The cross-correlation function (CCF) between edge velocity (Fig. 3A) and Rac1 activity (Fig. 3B) was calculated with Rac1 activity being fixed as described in the Materials and Methods (Fig. 3C). Temporal CCF was also extracted (Fig. 3C, dotted line) to analyze the lag time of signaling and morphological change (Fig. 3D; supplementary material Fig. S3; Table S1).

The positive peak value of the cross correlation coefficient was obtained with a time delay of −2.2 minutes [95% confidence interval (CI), −3.0 to −1.4 minutes] from Rac1 activity, indicating that the edge velocity was increased before the increase in Rac1 activity and vice versa (Fig. 3E). However, taking the interval of image acquisition (every 1 minute) into account, the Cdc42 activity was slightly delayed or nearly synchronized with morphological change (time delay of −1.3 minutes, 95% CI, −2.2 to −0.5 minutes; supplementary material Fig. S4A). RhoA activity did not demonstrate any correlation with edge velocity (peak value of 0.082; supplementary material Fig. S4B). This was probably because the high RhoA activity at both protrusion and retraction during cell migration marginally contributed to the overall correlation, leading to an offset of the net correlation coefficient value (Kurokawa and Matsuda, 2005). Furthermore, F-actin accumulation was coupled to the edge velocity with a time delay of −0.25 minutes (95% CI, −1.5 to 1.0 minute; Fig. 3F). We time-lapse imaged Rac1 activity and F-actin accumulation every 30 seconds in migrating HT-1080 cells, demonstrating the time delay values of −2.6 minutes (95% CI, −3.3 to −1.5 minutes) and −0.9 minutes (95% CI, −1.4 to −0.4 minutes) for Rac1 activity and F-actin, respectively (supplementary material Fig. S5, Table S1). These values are consistent with those obtained by imaging with the time interval of 1 minute. Of note, we previously demonstrated that EGF-induced Rac1 activation preceded the lamellipodial extension in Cos1 cells (Kurokawa et al., 2004). Thus, the role of Rac1 appears to be different in growth factor-induced and stochastic membrane protrusions.

Cross-correlation analysis revealed that PtdIns(3,4,5)P₃ accumulation was almost synchronized with membrane protrusion (time delay of −1.1 minutes, 95% CI, −1.8 to −0.5 minutes; Fig. 3G). The downstream metabolite, PtdIns(3,4)P₂, demonstrated a time delay of −1.6 minutes in comparison with the edge velocity (95% CI, −2.3 to −0.8 minutes; supplementary material Fig. S4C). PtdIns(4,5)P₂ accumulation was delayed 1.4 minutes (95% CI, −2.0 to −0.7 minutes; supplementary material Fig. S4D). In migrating MDCK cells, phosphatidyl inositol 5-kinase and PLC-γ, but not PI3K, mainly contribute to the distribution of PtdIns(4,5)P₂ at lamellipodia of the leading edge, suggesting a possible role of phosphatidyl inositol 5-kinase and PLC-γ in random migration in HT-1080 cells (Nishioka et al., 2008). As expected, the negative control FRET biosensor did not show any direct correlation with edge velocity (Fig. 3H).

In brief, the cross-correlation analysis showed a sequential relationship between signaling molecules and morphological change in the following temporal order: (1) PtdIns(3,4,5)P₃ accumulation and morphological change (polymerized actin localization); (2) PtdIns(3,4)P₂ accumulation and Cdc42 activation; and (3) Rac1 activation.

Rac1 induces actin reorganization and, consequently, lamellipodial formation (Hall, 1998). According to this scheme, the delay in Rac1 activation shown in Fig. 3 implies a positive feedback regulation from actin polymerization to Rac1 in migrating HT-1080 cells. To explore this assumption, we examined the effects of the F-actin depolymerizers latrunculin B and cytochalasin D on PtdIns(3,4,5)P₃ accumulation and Rac1 activity. Latrunculin B (Fig. 4A,B) and cytochalasin D (supplementary material Fig. S6A,B) treatment did not induce any change in the level of PtdIns(3,4,5)P₃, but were associated with a decrease in Rac1 activity (Fig. 4C,D). These results indicated that actin polymerization activated Rac1 by a positive feedback mechanism. Cdc42 activity was not suppressed by latrunculin B or cytochalasin D treatment (data not shown), although the reason for this remains unknown. Therefore, we focused on only Rac1 activity in the following study. Next, we investigated another positive feedback mechanism from Rac1 to PI3K (Aoki et al., 2007; Weiner et al., 2002) by the inducible translocation technique (Inoue et al., 2005). Rapamycin triggers heterodimerization formation of a Lyn N-terminal sequence-tagged FRB (LDR) with an FKBP-fused protein. Upon rapamycin treatment, FKBP–Tiam1 was translocated to the plasma membrane, and consequently activated Rac1 within 2 minutes, followed by the formation of lamellipodia (Fig. 4E,G). However, rapid activation of Rac1 did not increase the level of PtdIns(3,4,5)P₃ (Fig. 4F,H). These findings demonstrated that the positive feedback pathway consisting of Rac1 and PI3K did not contribute to the signaling network in migrating HT-1080 cells as reported in HL-60 neutrophils (Inoue and Meyer, 2008).

The results of auto-correlation analysis, which showed a wave-like pattern of membrane dynamics, prompted us to investigate a negative feedback loop involving Rac1, phosphatidyl inositols and morphological change. We found that treatment with a PI3K inhibitor, LY294002, induced an immediate and transient decrease in Rac1 activity to below the basal level, and suppressed lamellipodial protrusions (Fig. 4J,L). Interestingly, Rac1 activity and the morphological repression returned to the basal levels within 20 minutes even though PtdIns(3,4,5)P₃ was sustained at low levels for at least 30 minutes after treatment with LY294002 (Fig. 4I,K). We confirmed these findings using another PI3K inhibitor, PIK-93 (Knight et al., 2006) (supplementary material Fig. S6C,D).

The transient suppression of Rac1 by PI3K inhibition might be caused by adaptation, i.e. the ability of the signaling network to reset itself after responding to a stimulus. A recent study revealed that an adaptation emerges from only two major topologies of a signaling network: a negative feedback loop with a buffering node and an incoherent feed forward loop with a proportioner node (Ma et al., 2009). This prompted us to explore possible kinetic models of PI3K, Rac1 and actin polymerization modules by considering the following two constraints: (1) the inhibition of actin polymerization induces a sustained decrease in Rac1 activity, and (2) Rac1 activity and actin polymerization adapt to PI3K inhibition. Because Rac1 activation did not affect PI3K, we set PI3K as an input. Actin polymerization was set as an output for the membrane protrusion. In preliminary simulations, Rac1 was located between PI3K and actin polymerization (supplementary material Fig. S7). None of the models derived from this topology could fulfill the aforementioned constraints because actin polymerization cannot positively regulate Rac1.

Therefore, four models of possible networks with an additional ‘node’ were constructed to reproduce the experimental results obtained by cross-correlation and perturbation analyses (Fig. 5A). In these models, we assumed that actin polymerization corresponded to the edge velocity. The numerical simulation in every proposed model reproduced the effect of actin polymerization inhibitors and PI3K inhibitors on the change in Rac1 activity – namely, sustained decrease and adaptation, respectively (Fig. 5B, Model A; supplementary material Figs S8–S10 Models B–D, Tables S2–S5).

To distinguish these four models, we searched for candidate molecules that might act as nodes in the signaling scheme. In particular, we focused on negative regulators for actin reorganization and Rac1 signaling, and found that myosin light chain kinase (MLCK) and/or its target signaling played a role in a buffering node of the network. To visualize MLCK accumulation, the high molecular mass isoform of human MLCK1 tagged with GFP at the N-terminus, hereafter referred to as EGFP–Long-MLCK, was expressed in HT-1080 cells (Poperechnaya et al., 2000). EGFP–Long-MLCK demonstrated clear accumulation at not only the retracting tail, but also the leading edge of migrating HT-1080 cells (Fig. 6A; supplementary material Movie 9). The temporal CCF demonstrated a positive peak of correlation efficiency at the time delay of approximately −4.0 minutes (95% CI, −4.9 to −3.1 minutes; Fig. 6B, arrow). Furthermore, the negative correlation coefficient peak was obtained at the positive time delay of approximately 0.5 minutes (95% CI, −0.8 to 1.8 minutes; Fig. 6B, arrowhead). These delay timing values did not correlate with the expression level of EGFP–Long-MLCK, negating the artificial effect of over-expression of Long-MLCK on the migration in HT-1080 cells (supplementary material Fig. S11). To validate the relationship between MLCK and Rac1, we examined whether rapid activation of Rac1 induced translocation of MLCK to the periphery of the membrane. As expected, rapamycin-induced heterodimerization between LDR and FKBP–Tiam1-activated Rac1, which in turn increased the level of localization of EGFP–Long-MLCK at the peripheral membrane (Fig. 6C, arrows) with a delay of 4.6 minutes (Fig. 6D). In addition, treatment with an MLCK inhibitor, ML-7, induced a rapid increase in Rac1 activity (Fig. 6E,F). Thus, a reaction between Rac1 and MLCK could be the buffering node, which was postulated in Model A (Fig. 5). Finally, we examined the involvement of Rho-associated protein kinase (ROCK), which is well known to regulate myosin contractility (Totsukawa et al., 2004). ROCK inhibitor, Y27632, was added to the medium of HT-1080 cells expressing Raichu–Rac1 to determine whether Y27632 increased Rac1 activity as ML-7 does. Unexpectedly, however, Y27632 had little effect on Rac1 activity, although Y27632 induced typical morphological changes in ROCK inhibition: long protrusions and a sickle-cell shape (supplementary material Fig. S12). These results strongly supported the negative feedback loop in Model A, in which the MLCK pathway acted as a buffering node for the negative feedback.

In this study, we have revealed a signaling network comprising positive and negative feedback regulations between PI3K, actin polymerization, Rac1 and MLCK in randomly migrating HT-1080 cells (Fig. 7A). On the basis of the auto- and cross-correlation analyses, we can summarize the relative delay timing in comparison to edge velocity change (Fig. 7B) and the process of signaling and membrane protrusion and retraction as follows (Fig. 7C). In the first phase of protrusion, the edge velocity takes a positive value, indicating membrane extension and actin polymerization (0 minutes). Simultaneously or slightly after membrane protrusion, PtdIns(3,4,5)P₃ accumulates at the tip of the protrusion. In the second phase, after the initiation of membrane extension (2 minutes), Rac1 activity increases. In the third phase (4 minutes), MLCK accumulates at the tip of the leading edge in a Rac1-dependent manner, and consequently decelerates the edge velocity through actin reorganization. According to the results of the auto-correlation analysis, the lateral membrane wave revolves around the cell periphery with a periodicity of 40 minutes. Thus, the edge velocity, Rac1 activity and MLCK localization reach a peak at 10, 12 and 14 minutes, respectively. Meanwhile, the edge position of the cell periphery is described by integrating the edge velocity value, and peaks at 20 minutes after the start of the increase in edge velocity (Fig. 7C, gray line). The latter half of the process, i.e. retraction, represents the same sequential process in which the edge velocity and molecular activities show negative and low values, respectively.

We propose that a reaction–diffusion system of the signaling network generates the lateral propagation of waves around the cell periphery (Kondo and Miura, 2010). The positive feedback loop from Rac1 to actin polymerization amplifies the increase in the edge velocity and the lamellipodial protrusion. This actin polymerization signal diffuses laterally along the cell periphery. Meanwhile, the negative feedback loop comprising MLCK as a buffering node and Rac1 as the input, respectively, retracts the membrane with a time delay. Consistent with this idea, Otsuji et al. have reported a mass-conserved reaction diffusion with positive and negative feedback loops, which reproduces characteristic behaviors of membrane dynamics, including lateral wave propagation (Otsuji et al., 2010). Our data provide a molecular basis for their mathematical modeling. Furthermore, Arai et al. have found that self-organized waves of PtdIns(3,4,5)P₃ are generated spontaneously on the membrane of Dictyostelium cells in the presence of latrunculin A, and they have developed a reaction diffusion model that recapitulates the relaxation oscillator of phosphatidylinositol dynamics with positive and negative feedback loops (Arai et al., 2010). This model differs from the present models of HT-1080 cells in terms of the signaling hubs: phosphatidylinositols play a central role in the positive and negative regulations in Dictyostelium cells, whereas Rac1 acts as a signaling hub for two feedback loops, as observed in HT-1080 cells. Of note, the velocity of the lateral membrane wave in this study, ~4 μm/minutes (Table 1), was comparable with those in mouse embryonic fibroblasts (2.4 μm/minutes), fly cells (0.72 μm/minutes), T cells (0.33 μm/minutes) (Döbereiner et al., 2006) and Dictyostelium cells (2–4 μm/minutes as roughly calculated from Arai et al. (Arai et al., 2010), implying the presence of general molecular mechanisms underlying the random migration of eukaryotic cells.

One of the key findings in this study was that a negative feedback loop from Rac1 to MLCK enables the system of actin reorganization to adapt to external perturbations. We assumed that GTP–Rac1 recruits MLCK to the edge of membrane protrusions either directly or indirectly, leading to the retraction of the peripheral plasma membrane by increasing the actomyosin contractility. Consistent with this model, MLCK localizes at the leading edge of lamellipodia in mouse embryonic fibroblasts, and regulates the periodic lamellipodial contractions (Giannone et al., 2004). The long MLCK directly binds to actin filaments in vitro and in interphase cells, and these bindings are mediated by five DXRXXL motifs (Poperechnaya et al., 2000; Smith et al., 2002). The inhibition of ROCK by Y27632 did not alter Rac1 activity (supplementary material Fig. S12), suggesting two possible molecular mechanisms of the negative feedback. First, the negative feedback from MLCK to Rac1 was mediated through myosin that was regulated by MLCK, but not ROCK. There are many reports showing spatial and functional discrepancies of myosin regulations between MLCK and ROCK (Totsukawa et al., 2004). Second, myosin-independent MLCK activity mediates the negative feedback from MLCK to Rac1. Consistent with this hypothesis, myosin-independent effect of MLCK on migration of mammalian cells has been reported (Kudryashov et al., 1999; Niggli et al., 2006). In disagreement with the negative feedback model, the rapid and global activation of Rac1 by FKBP–Tiam1 and LDR upon rapamycin stimulation resulted in global lamellipodial extensions (Fig. 4E,F). This discrepancy is presumably attributable to the difference of strength and duration in Rac1 activation between randomly migrating cells and cells expressing FKBP–Tiam1. Robust and sustained activation of Rac1 by the FKBP–Tiam1 translocation system might cause inactivation of MLCK through phosphorylation by PAK, which is a downstream effector of Rac1 (Sanders et al., 1999). However, it has been reported that PAK directly phosphorylates Ser19 of MLC and induces an increase in actomyosin contraction in endothelial cells (Chew et al., 1998; Zeng et al., 2000). Therefore, future studies will be needed to evaluate the spatial and temporal regulation of MLCK activity in migrating cells.

Evidence of a pathway from the actin cytoskeleton to Rac1 was provided by cross-correlation analysis (Fig. 3) and by the latrunculin-B-induced decrease of Rac1 activity (Fig. 4). In this context, we should note that Rac1 induces actin polymerization and formation of lamellipodia (Hall, 1998) (Fig. 4E,F). These data indicate the existence of a positive feedback loop between the actin cytoskeleton and Rac1. Although it remains unknown how the actin cytoskeleton regulates Rac1 activity, we can adduce several mechanisms of the positive feedback loop. Some GEFs acting on Rac1, such as FGDs, Frabin (Nakanishi and Takai, 2008) and Asef2 (Sagara et al., 2009), interact with F-actin directly or indirectly. These GEFs are candidate molecules for a role in the positive feedback. Furthermore, focal complexes and adhesions induced by Rac1 at lamellipodia recruit the Cas–Crk–DOCK signaling complex (Côté and Vuori, 2007; Meller et al., 2005) and lead to the further activation of Rac1 in a manner consistent with the positive feedback. Such positive feedback regulation between PtdIns(3,4,5)P₃, Rac1 and actin has been well characterized in Dictyostelium cells, neutrophils and PC12 cells (Aoki et al., 2007; Iijima et al., 2002; Merlot and Firtel, 2003; Weiner et al., 2002). In contrast to these earlier works, PtdIns(3,4,5)P₃ was not involved in the positive feedback in the present model (Fig. 7). Consistent with this, Inoue et al. have shown with HL-60 neutrophils that rapid activation of endogenous Rac1 proteins triggers effective actin polymerization in the absence of PI3K activation and increase in PtdIns(3,4,5)P₃ (Inoue and Meyer, 2008).

Cross-correlation analyses have been applied to understand how Rho GTPases coordinate the morphodynamics (Machacek et al., 2009; Tsukada et al., 2008). Tsukada et al. have demonstrated, using the edge evolution tracking (EET) method, that the change in cell edge velocity precedes the activation of Rac1 and Cdc42, with a time delay of ~6 and 8 minutes, respectively (Tsukada et al., 2008). By contrast, Machacek et al. have shown with a sophisticated image processing technique that Cdc42 and Rac1 are activated ~40 seconds after the increase in the edge velocity (Machacek et al., 2009). Differences in the values of the time delay between these studies and our study, 2.2 minutes for Rac1, could arise from two reasons: the differences in image processing algorithms, and the use or not of a classification system for the modes of cell migration. With respect to the first reason, Machacek et al. analyzed the local area of the cell periphery, i.e. only the protrusive area, with image sets acquired for a shorter time interval of 10 seconds than ours of 1 minute. We, instead, focused on the global membrane dynamics of protrusion and retraction in the entire cell periphery of migrating HT-1080 cells. This could cause larger time delay values than those obtained in Machacek's study. With respect to the second reason, we classified the modes of cell migration into four categories – wave-like, oscillatory, directional and unclassifiable – and analyzed only the cells that were categorized into the wave-like migration. By contrast, Tsukada et al. did not use any categorization. Time-delay values obtained from directionally migrating cells became much greater than those of the cells migrating in the other modes. Therefore, the use of classification resulted in shorter delay timing in our study than Tsukada's studies (Tsukada et al., 2008). Notably, we assumed that the proposed feedback regulations between Rac1, actin cytoskeleton and MLCK took place in both protrusions and retractions in migrating HT-1080 cells for the following reasons. (1) Data obtained from cross-correlation analysis at the front of migrating cells contributed to the total correlation coefficient value to a similar extent as that obtained at the back of migrating cells. (2) We did not recognize any difference in the correlative relationship of Rac1 or Cdc42 activities and edge velocity between the front and back of cell movements; namely, at the front, increase in edge velocity preceded increase in Rac1 activity, whereas, at the back, decrease in edge velocity preceded decrease in Rac1 activity. It was also noted that a correlation-based multiplexing approach should be carefully applied to a system with feedback interactions (Sabouri-Ghomi et al., 2008). In our study the delay timing obtained by correlation analysis was within 2.2 minutes for Rac1 (Fig. 3); whereas, the effect of feedback regulation appeared at 4–6 minutes later (Fig. 6). Thus, it is unlikely that the feedback regulation interferes with the correlation analysis.

By employing image processing, statistical analysis and pharmacological perturbations, this study has demonstrated the molecular mechanisms underlying the ordered membrane dynamics in randomly migrating eukaryotic cells. A pivotal question to be answered is what factors initiate the increase in edge velocity prior to Rac1 activation. Our data showed that PtdIns(3,4,5)P₃ accumulated concomitantly with the increase in edge velocity (Fig. 3). Additionally, it has been recently reported that membrane curvature induces actin polymerization through EFC and/or F-BAR proteins and the N-WASP–WIP complex (Takano et al., 2008). In line with this finding, maximal F-actin assembly was observed ~20 seconds after maximal edge advancement (Ji et al., 2008). These results suggest the involvement of lipid regulation, which serves as a trigger for membrane protrusion. Further studies will be needed to uncover the molecular mechanism of the dynamics in migratory cells.

The plasmids encoding FRET biosensors have been described as follows: Raichu-Rac1/1011x, pRaichu-Cdc42/1054x (Itoh et al., 2002), pRaichu-RhoA/1294x (Nakamura et al., 2005; Yoshizaki et al., 2003), pPippi-PtdIns(3,4,5)P₃/1747x, pPippi-PtdIns(3,4)P₂/1759X (Aoki et al., 2005; Sato et al., 2003), pPippi-PtdIns(4,5)P₂/1758X (Nishioka et al., 2008) and pRaichu-Pak-Rho/1110x, a negative control (Kurokawa and Matsuda, 2005).

HT-1080 cells were purchased from the American Type Culture Collection or the Japan Cell Resource Bank, and maintained in DMEM (Sigma-Aldrich, St. Louis, MO) supplemented with 10% fetal bovine serum (FBS). The cells were plated on 35-mm glass-bottomed dishes (Asahi Techno Glass, Tokyo, Japan), which were pre-coated with collagen type I (Nitta Gelatin Inc., Osaka, Japan). Plasmids were transfected into HT-1080 cells with Lipofectamine 2000 according to the manufacturer's protocol (Invitrogen, San Diego, CA). cDNA of Lifeact, a 17-amino-acid peptide (Riedl et al., 2008), were subcloned into the pCAGGS-EGFP vector. EGFP–Long-MLCK was a kind gift from Anne R. Bresnick (Poperechnaya et al., 2000). Expression plasmids for ERed–NES, LDR (Lyn11-targeted FRB), FKBP and FKBP–Tiam1 have been described previously (Aoki et al., 2007; Inoue et al., 2005). cDNA for FKBP–iSH2 was subcloned into pCAGGS-GST vector to obtain pCAGGS-GST-FKBP-iSH2 (Inoue and Meyer, 2008). LY294002 were obtained from Sigma-Aldrich (St. Louis, MO). Latrunculin B, cytochalasin D, ML-7 and Y27632 were purchased from Calbiochem (La Jolla, CA). PIK-93, which is a PI3K inhibitor, was obtained from Symansis NZ (Timaru, New Zealand). Rapamycin was purchased from LC Laboratories (Woburn, MA).

FRET imaging was done essentially as previously described (Aoki and Matsuda, 2009). Briefly, HT-1080 cells expressing biosensors were suspended in trypsin, plated on collagen-coated 35-mm glass-bottomed dishes, and left for approximately 1 hour. At the time that cells were being plated, the culture medium was replaced with Phenol-Red-free DMEM/F12 containing 10% FBS, and overlaid with mineral oil to prevent evaporation. Cells were imaged with an inverted microscope (IX81; Olympus, Tokyo, Japan) equipped with a 60× objective lens (Olympus), a cooled CCD camera (Cool SNAP-K4; Roper Scientific), an LED illumination system (CoolLED precisExcite; Molecular Devices), an IX2-ZDC laser-based autofocusing system (Olympus) and an MD-XY30100T-Meta automatically programmable XY stage (Sigma Koki, Tokyo, Japan). The following filters were used for the dual emission imaging studies: excitation filters, 435/20 for CFP and FRET (Olympus), S492/18x for GFP (Chroma) and 580AF20 for RFP (Omega Optical, Brattleboro, VT); dichroic mirrors, XF2034 for CFP and FRET (Omega) and 86006bs for GFP and RFP (Chroma); emission filters, 480AF30 for CFP (Omega), 535AF26 for FRET and GFP (Omega) and 635DF55 for RFP (Omega). The exposure time was 400 mseconds for CFP and FRET images when the binning was set to 4×4. After background subtraction, FRET/CFP ratio images were created with MetaMorph software (Universal Imaging, West Chester, PA), and represented by the intensity modulated display mode. In the intensity-modulated display mode, eight colors, from red to blue, are used to represent the FRET/CFP ratio, with the intensity of each color indicating the average intensity of FRET and CFP.

The velocity of the leading edge in migrating HT-1080 cells was calculated using a mechanical model with some modifications (Machacek and Danuser, 2006). The process was divided into two steps: (1) image processing and (2) velocity calculating.

1. Image processing. First, FRET or GFP images were binarized by the threshold value, which was determined arbitrarily to identify the cell edge (supplementary material Fig. S1A,B). Second, an Erode operation was applied to the binarized images three times to reduce the noise of the periphery (supplementary material Fig. S1C). Third, the outline of the eroded images was captured by line scanning to correct the data of the two-dimensional coordinates and FRET/CFP ratio values of the cell periphery (supplementary material Fig. S1D). These processes were implemented using MetaMorph software (Universal Imaging). There was no use of interior positions, only the edge positions.

2. Velocity calculating. First, the data of the XY coordinates and FRET/CFP ratio values, which were obtained by the image processing, were normalized to 200 sampling windows. Second, to quantify the edge velocity, the window that was located at an angle of 0°from the centroid was set as an initial point (Fig. 2A). Third, normal vectors were obtained along the cell boundary (supplementary material Fig. S1E). Normal vector angles were affected especially in the region with strong boundary deformation (Machacek and Danuser, 2006). Therefore, fourth, the normal vector angles were corrected by averaging angles of five backward and five forward vectors (supplementary material Fig. S1F). Lastly, by using the averaged vectors, the boundary displacement at the time T was obtained by calculating the minimum distance from the point at time T to the point at time T+1, and divided by the interval time of imaging to obtain the edge velocity (supplementary material Fig. S1G). We denoted the matrices of edge velocity and molecular activity/localization as {M(x,t)|0≤x≤200, 0≤t} and {F(x,t)|0≤x≤200, 0≤t}, respectively, where M is spatio-temporal data of morphological changes and F is spatio-temporal data of molecular activity or molecular localization. These processes were implemented using Matlab software (version R2010b; The Mathworks Inc., Natick, MA).

The auto- and cross-correlation analysis was performed essentially as described previously (Maeda et al., 2008). The auto-correlation functions of M(x,t) and F(x,t), and the cross-correlation function between M(x,t) and F(x,t) are defined as: where the numbers M_ave and F_ave indicate the average value of M(x,t) and F(x,t), respectively. The operator ST_x,t denotes the averages over time and over space. The range of calculation is the entire period of measurement in time and space from 0 to 200. Correlation analyses were implemented using Matlab software.

Random migration patterns of HT-1080 cells, i.e. wave-like, oscillatory, directional, and non-classifiable migration, were categorized on the basis of the ordered pattern of auto-correlation function in edge velocity, as shown in Fig. 2F. First, non-classifiable disordered pattern of the auto-correlation function was discriminated from the other migration patterns. Then, directional migration was classified as a value of Δt₀>30 minutes defined as a threshold criteria in the temporal auto-correlation function (ΔSpace=0, dashed line in Fig. 2D). Finally, wave-like migration was distinguished from oscillatory migration by the typical wavy shapes in the auto-correlation function.

All kinetic reactions in Models A–D (Fig. 5) were described by ordinary differential equations based on Michaelis–Menten kinetics. Equations and parameters in Models A to D are listed in supplementary material Tables S2–S5, respectively. Numerical simulation was performed using Matlab software. To achieve the adaptation, we assumed one constraint: the reaction rates and Michaelis constants of the buffering node (k_node and K_{m node}) were set to smaller values than those of their substrate in order to saturate their reaction rates (Ma et al., 2009).

We thank A. R. Bresnick for the plasmids and Y. Sakumura for the helpful discussion. Y. Inaoka, K. Hirano, R. Sakai, A. Katsumata, N. Nonaka, and A. Kawagishi are also to be thanked for their technical assistance. We are grateful to the members of the Matsuda Laboratory for their helpful input.