Predicting Individual Cell Division Events from Single-Cell ERK and Akt Dynamics

Predictive determinants of stochastic single-cell fates have been elusive, even for the well-studied mammalian cell cycle. What drives proliferation decisions of single cells at any given time? We monitored single-cell dynamics of the ERK and Akt pathways, critical cell cycle progression hubs and anti-cancer drug targets, and paired them to division events in the same single cells using the non-transformed MCF10A epithelial line. Following growth factor treatment, in cells that divide both ERK and Akt activities are significantly higher within the S-G2 time window (∼8.5-40 hours). Such differences were much smaller in the pre-S-phase, restriction point window which is traditionally associated with ERK and Akt activity dependence, suggesting unappreciated roles for ERK and Akt in S through G2. Machine learning algorithms show that simple metrics of central tendency in this time window are most predictive for subsequent cell division; median ERK and Akt activities classify individual division events with an AUC=0.76. Surprisingly, ERK dynamics alone predict division in individual cells with an AUC=0.74, suggesting Akt activity dynamics contribute little to the decision driving cell division in this context. We also find that ERK and Akt activities are less correlated with each other in cells that divide. Network reconstruction experiments demonstrated that this correlation behavior was likely not due to crosstalk, as ERK and Akt do not interact in this context, in contrast to other transformed cell types. Overall, our findings support roles for ERK and Akt activity throughout the cell cycle as opposed to just before the restriction point, and suggest ERK activity dynamics are substantially more important than Akt activity dynamics for driving cell division in this non-transformed context. Single cell imaging along with machine learning algorithms provide a better basis to understand cell cycle progression on the single cell level.


Introduction
The mammalian cell cycle is in large part driven by growth factor activation of the Ras-ERK 1-5 and the PI3K-Akt 2,6-11 pathways. Growth factors cause auto-phosphorylation of receptor tyrosine kinases (RTKs). For the ERK pathway, RTK phosphorylation recruits the guanine exchange factor SOS to the membrane, catalyzing the exchange of GDP for GTP bound to Ras, initiating Raf activation 6,[12][13][14][15][16] . This in turn activates the MEK-ERK phosphorylation cascade. When activated, the effector kinases ERK1/2 translocate from the cytoplasm to the nucleus and activate transcriptional regulators such as Elk-1 and CREB 17,18 19 . These transcriptional regulators induce immediate early genes (IEGs) like c-fos 17,18 that then contribute to the expression of cyclin D1 4,6,8,[19][20][21] , a key step in S-phase entry 22 .
Transient ERK activity maintains G2 arrest, whereas sustained ERK activity promotes escape by reducing p53 levels, and inducing the expression of pro-mitotic factors such as Plk1 and cyclin B 41 . Akt activity also contributes to G2-M progression as its inhibition is associated with reduced cyclin B levels, promoting Chk1 activity and G2 arrest 42 .
These observations motivate a closer look at determining how ERK and Akt dynamics are informative of cell cycle completion after the canonical restriction point.
On a single cell level, both ERK and Akt activity dynamics have substantial cellto-cell and dynamic variation, exhibiting complex pulses and more simple steady activity 1,27,[45][46][47][48][49][50][51][52][53] . Such variation, when coupled with the observations that cell cycle progression is also heterogeneous 1,54,55 , have prompted investigations into the correlation between dynamics and cell cycle fate in single cells. What determines proliferation on a single cell level? What relative contributions do ERK and Akt activity have to the decision of individual cells to divide? Much prior work has focused on activity dynamics. Both Ras-ERK [56][57][58] and PI3K-Akt 59 exhibit biphasic growth factor-induced activation dynamics, with a transient peak followed by sustained activity hours later. The dynamics of each phase contributes differently towards driving progression to S-phase and is cell type dependent [59][60][61] . Live-cell imaging and analysis of recently divided sister cells reveal that time-integrated ERK activity has some predictive power of the timing to S-phase entry 1 . 5 Time-integrated ERK dynamics were found to influence proliferation decisions in daughter cells 50 . Predicting PC12 cell differentiation/proliferation decisions required both ERK and Akt activity dynamics to best define the decision boundary between these two cell fate outcomes 62 . Yet, the extent to which both ERK and Akt activities throughout the cell cycle are predictive of division in single cells remains unclear.
Here, we use live-cell imaging to pair measurements of growth factor-induced ERK and Akt activity to cell division outcomes in the same single cells. We aim to assess the extent to which these activities are associated with cell cycle progression beyond S-phase entry, and to evaluate their ability to predict single cell division responses jointly in the well-established non-transformed breast epithelial MCF10A cell line, a model system that is commonly used to study epithelial signaling biology and cell division control 26,[63][64][65][66][67][68] . We found that following treatment of synchronized cells with growth factors EGF and insulin, both ERK and Akt activity are significantly higher within the S-G2 interval in dividing cells. Such differences were much smaller in the pre-Sphase window, which is traditionally associated with ERK and Akt activity dependence [59][60][61] , suggesting unappreciated roles for ERK and Akt in S through G2.
These higher activities could classify division events with AUC=0. 76. Surprisingly, ERK activity dynamics alone enable AUC=0.74, suggesting Akt activity dynamics contribute little to the decision governing cell division in this context. Interestingly, we found that ERK and Akt activities are less correlated in cells that divide. Network reconstruction experiments demonstrated that this correlation behavior was not due to crosstalk, as ERK and Akt do not interact in this context, in contrast to other cell types 69 . Overall, our findings support roles for ERK and Akt activity throughout the cell cycle as opposed to 6 just before the restriction point, and suggest ERK activity dynamics are substantially more important than Akt activity dynamics for driving cell division in this nontransformed context. 7 46,70 or Akt 48 kinase translocation reporters (KTRs), but not yet both at once, and paired those single cell dynamics to cell division events (Fig. 1A). We first verified that cell cycle progression and division are related to ERK and Akt activity dynamics in MCF10A cells using small molecule inhibitor experiments ( Fig S1). KTR-expressing cells were G0-synchronized by serum and growth factor starvation for 24 hours. After acquiring 1 hour of baseline ERK or Akt activity, cells were treated with EGF and insulin, growth factors that promote cell division in MCF10A cells 71 . Images were acquired every 15 minutes for 48 hours, and then single-cell data for kinase activity and division outcome were extracted using custom-built image processing pipelines (see Methods).

Predictability of Cell Division Events from
Dynamic regimes of KTR specificity were determined using two independent (four total) MEK and Akt inhibitors (Fig. S2). ERK KTR was found to be specific in all regimes explored here, whereas the Akt KTR was found to be specific >~ 1 hour after EGF and insulin co-stimulation.
Single cell traces of ERK or Akt activity (thin lines) along with the population median (bold line) show rapid activation following growth factor treatment, which largely persists for the duration of the experiment, without recognizable pulsing (Fig. 1B,C). In (blue) dividing cells, population median ERK and Akt activity dynamics are higher throughout the cell cycle compared to non-dividing cells, with larger differences evident for ERK. In the pre-S-phase entry window (~< 8 hours after growth factor treatment), 8 there are slight differences between dividing and non-dividing cells in terms of population median ERK and Akt dynamics. These differences grow larger in the subsequent 8.5-40 hour interval post growth factor addition, which largely corresponds to S and G2 phases. These trends were also evident with 10-fold less concentration of growth factors (Fig. 1C). These results suggest that ERK and Akt activity may have importance after S through G2 phase.
To assess the statistical significance of this finding, we calculated the median ERK or Akt activity for individual single cells within the 8.5 -40 hour window postgrowth factor treatment, and then compared median activity between dividing and nondividing cells with the rank-sum test (Fig. 1D-E). Individual dots in the boxplot represent the median ERK or Akt activity calculated within the 8.5 -40 hour interval in a single cell. These median single cell activities were significantly different in dividing vs. nondividing populations (Fig. 1E). Yet, there is substantial overlap in the two populations.
We evaluated whether single cell median ERK or Akt dynamics are predictive of cell division using a logistic regression classification model, and ROC analysis to quantify the outcome. Both ERK and Akt dynamics have some predictive power for cell division under high and low growth factor conditions, with high growth factor conditions having slightly elevated predictive power, as quantified by the area under the ROC curve (AUC). ERK dynamics have more predictive power than Akt dynamics. Yet, the best achieved AUC is 0.74, indicating there are other factors driving differences in cell division fate.

Predictability of Cell Division Events Using Measurements of Both ERK and Akt
Dynamics in the Same Single Cells. As shown above, ERK and Akt activity dynamics 9 alone contain information about subsequent cell division. Would simultaneous measurements of both ERK and Akt activity dynamics in the same single cell improve cell division predictions? To answer this question, we performed a similar experiment as  (Table 1). Thus, we conclude that the above 10 analysis focused on median activity is likely to be sufficient for assessing how predictive ERK or Akt dynamics are for cell division events in the studied context.
Inferring the Topology of the ERK-Akt Network. Information content is related to correlation, so we investigated the extent to which ERK and Akt dynamics in the same single cells were correlated, looking across every cell and every time point (Fig. 3A, replicate Fig. S3D). Interestingly, in dividing cells, single cell ERK and Akt dynamics within the 8.5-40 hour window are significantly less correlated than in non-dividing cells, at both high and low growth factor doses. Network topology can strongly influence correlated behaviors. In different studies, ERK and Akt have been reported to exhibit very different network behavior, such as cross-pathway activation, inhibition 6,20,21,69,[74][75][76][77][78][79] and non-interaction 23,69,80 . Factors such as cell type and growth factor context can influence these discrepant network topologies 3,69 . Previous work conducted in panel of growth factors and cell lines show varying probabilities of forming an interaction network edge between ERK and Akt 69 . The differences in network edge formation can affect downstream signaling and cell fate decisions 3 . Could ERK and Akt network topology be dynamic, and give insight into the division-related correlated behaviors observed above?
To reconstruct the ERK and Akt network in MCF10A cells, we implemented recent theory from our lab that specifies a sufficient experimental design for such tasks, based on perturbation time course data 81 . Specifically, for this 2-node network, three time course experiments should be done: response of ERK and Akt activity to EGF and Insulin co-treatment with (i) no inhibitor; (ii) an ERK pathway inhibitor; and (iii) an Akt pathway inhibitor. Additionally, we wanted to understand whether the network would be 11 different in the acute phase of growth factor treatment from a serum starved state vs.
the "chronic" condition where ERK and Akt activities are steady over time, particularly because these time regimes seem to have different biological information encoded for cell division decisions.
In the acute regime (Fig. 3B), MCF10A cells expressing either ERK or Akt KTR were seeded, serum and growth factor starved for 24 hours, and then pretreated with Akt activation (Fig. 3C). Akt inhibition ablates Akt activation and has a negligible effect on ERK activation (Fig. 3C). Although the Akt KTR may reflect kinase activity other than Akt in this acute stimulus regime, the fact that complete inhibition of the ERK pathway has negligible impact on the Akt KTR readout means that ERK does not impact Akt or the others. These results show that in the acute stimulus regime, ERK and Akt exhibit negligible crosstalk after treatment with EGF and insulin.
In the chronic regime, cells were pretreated with either EGF and insulin for 30 minutes followed by 30 minutes of baseline acquisition, leading to robust ERK and Akt activation ( Fig. 3D-E, Replicates in S5). Akt inhibition reduces Akt activity, as expected, but negligibly affects ERK activity. MEK inhibition reduces ERK activity, as expected, 12 but does not appreciably affect Akt activity. These conclusions are also consistent when a second set of MEK and Akt inhibitors are used (Fig. S4). Therefore, in the chronic regime ERK and Akt also do not exhibit appreciable cross pathway interactions after EGF and insulin co-treatment. We conclude it is unlikely that crosstalk interactions account for correlations that change in dividing vs. non-dividing cells. 13

Discussion
Binary single-cell responses, like division, to perturbations such as growth factor and drug treatments, are almost universally heterogeneous even in clonally derived populations. However, predictive biochemical features, present either before the perturbation, or from dynamics following the perturbation, are seldom known. The ability to predict such binary responses would not only reflect a deep and fundamental understanding of the systems governing important cellular responses, but also have significant translational applications such as antibiotic resistance, tissue engineering, and anticancer therapy, where the fates of single cells can be of great importance.
Here, we investigated growth-factor induced cell division fates in the well-studied, nontransformed mammalian epithelial cell line MCF10A, and how they may be predicted by the dynamics of two central signaling pathways, PI-3K/Akt and Ras/ERK. Answering such questions requires single-cell, non-destructive analysis of biochemical features, in this case ERK and Akt activities, that are paired to the eventual cell division outcome.
They also must be carried out in a high-throughput manner to observe enough events to make statistically-supported conclusions. After setting up this experimental system and understanding its ranges of validity, we learned that (i) ERK and Akt activities are higher in the 8.5-40 hour window after growth factor treatment in cells that divide, suggesting Nearly all the cells we observed had relatively simple dynamics for ERK and Akt activity, a rise then a somewhat constant higher than baseline steady-state. Other recent single cell studies have reported pulsatile ERK dynamics 1,49,83,84 . Some of this may be related to differences in growth factor concentration, the reporters used, being FRET-based 85 or translocation based 70 . No live-cell imaging probe is perfect and of course has its drawbacks, some of which may be related to off-target responses, which may partly explain our Akt activity data in the "acute" phase first following growth factor treatment. For example, kinases other than Akt may recognize and phosphorylate the FOXO1-based Akt KTR docking site [86][87][88][89][90] . EGF and insulin stimulation may also 15 promote activation of such kinases including PLK1 91 , SGK and PKA 92 . Another aspect may have to do with cell-cell contact and density. In our study, cells were seeded at low density and serum/growth factor-starved prior to analysis, whereas pulsatile signaling was reported in high density environments in asynchronously cycling settings 1,84 . Yet others have found pulsatile dynamics can induce different sets of genes as compared to sustained dynamics 46 . However, phenotypic consequences, at least in terms of cell proliferation still seem to be related to simple time-integrated signaling dynamics 1,50 , similar to what we found here.
ERK and Akt activity dynamics are only a subset of the potentially important variations that drive phenotypic variability in cell division responses, as shown by the AUC=0.76 that was obtained. ERK dynamics account for nearly all this predictive power. This reinforces Akt activity as perhaps more relevant for cell maintenance and health, and more as a "checkpoint" for division but not a significant driver, at least in the studied system. As noted above, cell contacts and density are important. Such phenomena may potentially be controlled through micropatterning experiments, where cell shape and placement can be carefully controlled 93,94 . Cell "state", corresponding to different epigenetic and/or metabolic states of cells prior to the experiment, has been reasonably well documented ubiquitously, and can contribute to variability, although is difficult to assess in the "track and follow" manner that can be done with live-cell kinase reporters. Metabolic or organelle abundance variability may also contribute 95,96 . Of course, there are other pathways and biochemical correlates that are likely important, such as a balance between p53 and p21 and/or CDK2 activity 54,67,97 . Given the multitude of fluorescent proteins, and improvements in cell tracking from non-labeled 16 bright field images 98,99 , one may be able to measure more important biochemical readouts simultaneously for such purposes. There are also multiple checkpoints between growth factor treatment and cell division, such as the restriction point, and DNA damage checkpoints, that may contribute. Monitoring division with probes like the Fucci system that gives readouts of each cell cycle phase may help explore such phenomena 100 .
An interesting aspect of our study was the surprising larger differences between dividing and non-dividing cells in the time period that corresponds to S/G2 phases of the cell cycle, as opposed to pre-S-phase. The roles of growth factor signaling through ERK and Akt pathways historically focused on passing the restriction point into S-phase 101 .
Thus, our results suggest potential functional roles for ERK and Akt beyond this canonical understanding. Indeed, a recent study found time-integrated ERK activity in a mother cell's G2 phase influenced the cell cycle progression in the subsequently daughter cells 50 . The mechanisms that may be driving such functional roles are a potentially interesting area of future study.
We also studied the ERK and Akt activity network, since we found that ERK and Akt activity are less correlated with each other in dividing cells compared to non-dividing cells. We found that the observed differences in correlation are likely not a result of network topology as ERK and Akt do not appreciably interact. This lack of interaction is surprising given that some prior studies describe these pathway as exhibiting cross pathway interactions 23 80 . These studies reveal that in non-interacting pathways, differences in protein expression influence the flow of erythropoietin signaling 80 . Therefore, in our model system, it is possible that the observed differences in ERK and Akt correlation may arise from differences in protein expression across dividing and non-dividing cells.
It may also be that differences in upstream signaling capacity to ERK and Akt may be related. Characterizing the differences in protein expression level in single cells, and following cell signaling and cell division can provide insight; but this becomes a quite challenging experiment given the number of probes to be measured simultaneously.
Computational Image Analysis. While many image analysis tools exist [103][104][105] , each application still requires much novel development tuned to the problem at hand. We developed an automated image analysis pipeline using both iLastik 106 and CellProfiler 105 software packages, along with MATLAB scripts (Fig. S2C). It is available at the 1. Prior to segmentation images were flatfield corrected and background subtracted using CellProfiler. Images of nuclear localized fluorescent protein H2B-mRuby2 (ERK KTR) and NLS-mCherry (Akt KTR) were input into iLastik. Nuclei were identified using a series of features-object intensity, edge detection and texture.
2. The binary mask outputs from iLastik were input into CellProfiler to create a perinuclear ring known as the 'Cytoring' 107 which extends 10 pixels from the binary nuclei mask and into the cytoplasm (Fig S2C). Calculating the cytoplasmic to nuclear KTR fluorescence ratio provides the relative activity of the pathway of interest for that particular cell at that time point. 25 3. Segmented nuclei identified with iLastik were tracked using CellProfiler's TrackObjects module 104,105,108 and Follow Neighbors 108 . Each identified nucleus was assigned a numerical ID, which corresponds to the same cell across each timepoint. We filter tracks that are shorter than the duration of the time course to prevent quantification of cells that were transiently tracked. 4. Cell division was detected using a feature of cytoplasmic to nuclear KTR fluorescence (C/N ratio) that is unique to dividing cells. As cells divide, there is a change in morphology resulting in a rapid decrease in C/N ratio (Fig. S7). MATLAB's findpeaks function was used to detect when this steep decrease occurs. We then truncated the time series 5 timepoints before the identified peak, which is attributed to actual kinase activity. This was repeated 1000 times to define the range of correlation coefficients between the 5 th and 95 th percentiles, which was reported and rounded up to the nearest 0.01.
The hctsa package 72,73 was installed according to the published documentation. We initialized our data set from data acquired under high concentrations of EGF and insulin stimulation using the custom MATLAB script hctsa_prep_dualrep.m. The number of cells that were input into hctsa were chosen based on the smallest population of either dividing /non-dividing cells. A random permutation was then used to choose an equal number of the largest population of cells to input into hctsa. The data set was initialized using the hctsa's TS_init command, followed by TS_compute command. The processed datasets were labeled for hctsa format as 1-dividing or 0-non-dividing using 27 TS_LabelGroups. The top predictive features were identified using the raw computed hctsa values and the function TS_TopFeatures for both ERK and Akt time course data.               Figure S7