Differentiation correlates with reduced motility and increased actin intensity
We previously established that ERK1/2 inhibition induced robust, faster, and less temporally variable differentiation and fusion of primary myoblasts isolated from chick and mice, leading to the rapid formation of myotubes ex vivo within 24 h post induction (Eigler et al,
2021). Immunofluorescence staining of the differentiation markers MyoG and MyHC at different time points in cultures treated with the ERK1/2 inhibitor SCH 772984 (ERKi, 1 μM) or with DMSO as control (Morris et al,
2013) show that differentiation into fusion-competent myocytes is accompanied by the upregulation of MyoG, which initiates terminal differentiation (Figs.
1A and
EV1 and
EV2). MyHC expression begins 16 h post induction and peaks in the multinucleated myotubes at 24 h (Figs.
1A and
EV2). The number of MyoG-expressing cells increases over time, stabilizing at 14 h post induction (Figs.
1B and
EV1 and
EV2).
To characterize the dynamic behavior of differentiating myoblasts, we isolated primary myoblasts from mice co-expressing the nuclear marker tdTomato fused to a nuclear localization signal (tDTomato-NLS) (Prigge et al,
2013) and the F-actin marker LifeAct-EGFP (Riedl et al,
2010; Fig.
1C). Actin governs a range of processes during differentiation, from cell motility, elongation, and alignment to fusion and fibrillogenesis, ultimately leading to sarcomere formation post-fusion and the subsequent maturation of myotubes into myofibers (Rubinstein et al,
1976; Otey et al,
1988; Nowak et al,
2009; Luo et al,
2022).
To collect dynamic information on the continuous transition from proliferation to terminal differentiation we performed time-lapse widefield microscopy of large fields of view each containing ~3000 cells (Fig.
1C; Movie
EV1). Cultures were imaged for 23 h, starting 1.5 h after ERKi or DMSO (control) treatment. We observed that differentiation was accompanied by a decrease in cell motility, consistent with previous studies showing that myocytes are less motile than myoblasts ex vivo (Powell,
1973; Griffin et al,
2010) (Fig.
1D; Movie
EV1). The regulated reduction in cellular motility is beneficial for enhancing the cell–cell interactions that would initiate differentiation and subsequently promote fusion (Krauss et al,
2005; Buggenthin et al,
2017; Nowak et al,
2009; Luo et al,
2022). Similarly, we observed an increase in the fluorescence intensity of the F-actin marker most likely corresponding to the expression of muscle-specific actin isoforms (Fig.
1C,E). Cumulatively, these experiments demonstrate that the transition of myoblasts from the proliferative to the terminally differentiated state is accompanied by dynamic changes in motility and actin intensity.
Machine learning applied to time-series data generates quantitative single-cell differentiation trajectories
Following the association between differentiation and the population scale changes in actin intensity and motility, we hypothesized that the information encoded in single-cell migration trajectories and actin dynamics might be sufficient to computationally estimate a continuous score reflecting a myoblast’s gradual transition from an undifferentiated proliferative state to a terminally differentiated fusion-competent state. To test this hypothesis, we took a machine-learning approach: (1) extracting features from the motility/actin time series, (2) training machine-learning classification models (aka classifiers) to discriminate between the undifferentiated and differentiated states, and (3) using the confidence of these models as a quantitative measurement for cell state.
The first step in designing our machine-learning model was determining which cells and time frame can be considered as undifferentiated or differentiated. Myocytes must differentiate to become fusion-competent (Abmayr and Pavlath,
2012). Hence, we defined the cultures as differentiated for classification 2.5 h before the first fusion event was observed in the field of view. Undifferentiated cells were defined from the population grown in proliferation medium in the presence of DMSO, which continue to proliferate and remain undifferentiated, except a small fraction that begins to differentiate stochastically toward the end of the experiment due to the increase in cell density (Eigler et al,
2021). To enable continuous scoring along single-cell differentiation trajectories, we performed semi-manual single-cell tracking, where each trajectory was manually verified and corrected when necessary (Movie
EV2). Single-cell analysis confirmed our population-based results that differentiation was accompanied by a decrease in cell motility and an increase in the F-actin marker’s fluorescence intensity (Appendix Fig.
S1). We partitioned trajectories of undifferentiated and differentiated myoblasts to overlapping temporal segments of 2.5 h each, for an overall 16,636 undifferentiated and 47,819 differentiated temporal segments, extracted from 310 and 538 cells correspondingly, that were used for model training (Fig.
2A, top). From each temporal segment, we extracted the corresponding single-cell motility (dx/dt, dy/dt) and actin intensity time series. Single-cell motility/actin time-series features were extracted using the Python package “Time Series FeatuRe Extraction on the basis of Scalable Hypothesis tests” (tsfresh) that derives properties such as temporal peaks, descriptive statistics (e.g., mean, standard deviation) and autocorrelations (Christ et al,
2018). The extracted single-cell feature vectors and their corresponding undifferentiated/differentiated labels were used to train random forest classifiers (Breiman,
2001), which surpassed other machine-learning algorithms (Appendix Fig.
S2). The entire process is depicted in Fig.
2A and detailed in “Methods”.
We applied the trained motility and actin classifiers on single-cell trajectories from an experiment that was not used for training and attained a continuous quantification following the differentiation process by using overlapping temporal segments. At the population level, the single-cell-state classification performance gradually increased from an area under the receiver operating characteristic (ROC) curve (AUC) of ~0.6 to ~0.85 at 7.5–14.5 h from experimental onset (Fig.
2B,C). A classifier trained on features derived from both motility and actin time series surpassed each of the motility/actin classifiers, suggesting that motility and actin dynamics contain complementary information regarding the cells’ state (Fig.
EV3A–C). The AUC values of all models were well beyond the random value of 0.5, indicating that our classifiers can discriminate between undifferentiated and differentiated cells at the population level before appreciable cell morphological changes occur.
Next, we wondered whether we can use these classifiers to predict the differentiation state of a single cell. For a given temporal segment of a given cell, the classifier outputs a “confidence score” (i.e., differentiation score) that reflects the model’s certainty in its prediction. Low differentiation scores indicate that the cells are predicted as undifferentiated, while high scores indicate predicted differentiation. To interpret what temporal features were the most important for the model’s prediction, we applied SHapley Additive exPlanations (SHAP) (Lundberg and Lee,
2017) and used random forest’s feature importance algorithms (Breiman,
2001). Both interpretability methods highlighted temporal features related to high variance of acceleration rate or high complexity of actin intensity time series as dominant features driving the models’ prediction (Appendix Fig.
S3). We hypothesized that the differentiation score could be used as a continuous readout for the cell state. At the critical time frame of 7.5–14.5 h, at the population level, the differentiation scores of ERKi-treated cells gradually increased for the motility (Fig.
2D), the actin-based (Fig.
2E), and the combined (Fig.
EV3D,E) models while maintaining low scores for experiments of DMSO-treated cells. For the rest of the manuscript, we focused on analyzing the motility- and actin-based models, because showing that each of two independent models trained with different temporal readouts (motility, actin) can quantitatively monitor the cell differentiation process, and thus strengthening our goal of evaluating whether the confidence score of a classifier trained for a binary classification task can be used to continuously measure a biological process that evolves over time.
We conducted a single-cell analysis by measuring the Spearman correlation between single-cell differentiation score and time at the critical time interval of 7.5–14.5 h when differentiation occurs. This analysis indicated that most cells underwent a monotonic increase in differentiation scores over time (Fig.
2F). A similar gradual increase in differentiation score at 7.5–14.5 was observed when flipping the experiments used for training and testing (Appendix Fig.
S4), the differentiation score was not sensitive to the size of the temporal segment (Appendix Fig.
S5), nor to the window size used to measure actin (Appendix Fig.
S6), and was consistent across multiple independent trainings (Appendix Fig.
S7). Visualizing single-cell trajectories showed that most trajectories followed a gradual increase in their differentiation scores (Fig.
2G). Measuring the distribution of the per-cell differentiation scores’ temporal derivative (Appendix Fig.
S8A,B), their integration over time (Appendix Fig.
S8C–F), the predicted onset (Appendix Fig.
S9), and the predicted duration (Fig.
2H) of the differentiation process, suggested that the progression in single-cell differentiation is highly heterogeneous (Fig.
2H). These results suggest a heterogeneous gradual transition from an undifferentiated to a differentiated state within a typical time frame.
The motility- and actin-based classifiers’ predictions were mostly consistent at their single-cell predictions, showing high correlation beyond 0.5 for 62% of cells, and negative correlation for less than 4%, and providing further support that both models measure the continuous state transition (Appendix Fig.
S10A,B). Analysis of cells sub-groups partitioned according to the agreement between the motility and actin models, showed that lower agreement between the classifiers was associated with lower monotonicity of the differentiation scores of both models (Appendix Fig.
S10C). In most cases of disagreement between the motility and actin models, we noticed a deviation in the actin-based model when the cell entered a crowded region and/or crawled below other cells (Appendix Fig.
S11). High agreement between the two models was associated with the single cells’ motility persistence, the ratio between the direct translation (i.e., distance from start to end), and the overall distance traveled (Appendix Fig.
S10C). Qualitative and quantitative association between single cells’ motility persistence and motility-based (but not actin-based) differentiation score identified persistent motility as a functional marker for the intermediate states of myoblast differentiation (Fig.
EV4). The lack of association between the actin model and persistence suggests that the actin model encodes different dynamic properties that are linked to the differentiation. Moreover, this lack of association highlights the potential to use deviations between these models to discover mechanisms that uncouple the link between motility, actin dynamics, and myoblast differentiation. Altogether, our data suggest that machine learning can transform motility and actin dynamics to a quantitative readout characterizing the myoblast differentiation process at single-cell resolution describing a continuous myoblast state transition from an undifferentiated to the terminally differentiated states at 7.5–14.5 h post induction.
Models that discriminate between undifferentiated and differentiated states are not sufficient for the quantitative characterization of the continuous differentiation processes
Using the simplest readouts to quantify and delineate different biological conditions/states is always preferred because it provides more direct insight regarding the underlying mechanisms. Is it possible that our approach is overly complicated and exceeds what is required to quantitatively describe the differentiation process? Are straight-forward single-cell measurements sufficient to discriminate between undifferentiated and differentiated cells and follow the differentiation process? To test this possibility, we evaluated the discriminative performance of single-cell properties that are expected to deviate between the undifferentiated and differentiated cells. These included cell speed, actin intensity, the temporal derivative of actin intensity, migration persistence, and local density. The local density dramatically increased over time for cells grown in proliferation medium due to continued proliferation throughout the experiment (Fig.
3A). The mean speed and actin intensity in proliferating cultures slightly decreased and increased correspondingly over time, perhaps due to the increased density (Fig.
3B,C), and the mean temporal derivative of actin intensity fluctuated around zero for both differentiating and proliferating cells (Fig.
3D). Persistent migration of proliferating cells was lower compared to differentiating cells without a clear trend over time (Fig.
3E). Each of these four discriminative readouts, as well as their integration, could be generalized across experiments as demonstrated by using each feature to train a machine-learning model and applying this model to discriminate between the two experimental conditions in an independent experiment (Fig.
3F).
The model trained with the local density and the model trained with all four features surpassed the discrimination performance of the time-series motility and actin models (also reported in Fig.
2B,C). However, discrimination does not necessarily imply that these readouts can be used to quantitatively describe the differentiation process. Indeed, the differentiation score of each of these classifiers could not capture the differentiation process as measured by single-cell monotonic increase at the critical differentiation time interval of 7.5–14.5 h. The single-cell correlations between the differentiation score and time were low for all the single-feature classifiers, as well as for the integrated classifier (Fig.
3G). This is in contrast to our classifiers that generalized to effectively quantify the differentiation process leading to a higher correlation between the differentiation score and time (Fig.
3G—motility, actin intensity, same data as in Fig.
2G). A plausible explanation for why these effective discriminating models could not capture the continuous differentiation process is that the discriminating features captured properties attributed to the undifferentiated state. For example, the increased local cell density of proliferating cells can be used to effectively discriminate between the two populations but does not provide any information regarding the progression through differentiation. Indeed, training models that included temporal features extracted from single-cell local density dynamics showed the same or deteriorated correlation between the differentiation score and time compared to models that were not trained with local density (Appendix Fig.
S12). Motility persistence, which was the most monotonic among the single feature, and was earlier identified to be associated with differentiation (Fig.
EV4), was still far behind the integrated models in terms of monotonicity, suggesting that there is further discriminative information beyond persistence in the cells’ trajectories. These results indicate that effective discrimination between the discrete extreme states is insufficient for the quantitative characterization of continuous state transitions. Specifically, using machine learning for quantitative characterization requires extracting features that can capture the state transition and avoiding features that may confound the quantitative characterization of the process (e.g., avoiding local cell density in characterizing the differentiation process). In our case, and in agreement with other studies (Copperman et al,
2021; Wang et al,
2020; Wu et al,
2022), integration of multiple dynamic features encoding the temporal changes were necessary to continuously measure a biological process.
The transition from terminal differentiation to fusion is controlled by p38
We next aimed at harnessing our single myoblasts continuous differentiation scores to investigate the relationship between cell differentiation and fusion. We manually annotated the fusion time of 68 myoblasts that fused to 6 myofibers (Fig.
4A) and used the continuous differentiation score to determine an estimated time of the terminal differentiation state. Both the distributions of the single cells’ terminal differentiation and fusion times followed a normal-like distribution, where the variability in the predicted differentiation time was higher than that of fusion time (Fig.
4B). The time duration between terminal differentiation and fusion also followed a normal-like distribution, indicating a typical duration of ~3 h between terminal differentiation and fusion at the population scale (Fig.
4C). The mean of a ~3 h gap between predicted terminal differentiation and fusion is not trivially derived from the definition of differentiation timing during training because (1) the differentiation time at training was defined as 2.5 h before the first fusion event, and (2) the heterogeneity in fusion timing spans over ~10 h (Fig.
4B). These results suggest that cells undergo fusion within a typical time interval from their terminal differentiation. This coupling was validated by measuring a correlation between single-cell differentiation and fusion times (Fig.
4D) and was not sensitive to the threshold used to determine the terminal differentiation time (Appendix Fig.
S13). These results suggest that myoblasts must reach a differentiation checkpoint before fusion can proceed.
Previous studies have shown that co-inhibition of p38, a family of MAP kinases that play a critical role in the initiation of the differentiation program, together with a promyogenic factor, overcomes the early block in differentiation but not the later impairment of muscle cell fusion imposed by the p38 inhibitor, leading to differentiated unfused cells (Gardner et al,
2015). Following this logic, we treated primary myoblasts with the promyogenic ERKi and the p38 inhibitor BIRB 796 (p38i; 5 µM) and performed primary myoblasts live imaging experiments. There was little appreciable difference between cells treated with p38i and controls treated with DMSO, consistent with previous studies showing that p38i maintains myoblasts in a proliferative undifferentiated state (Zetser et al,
1999). Myoblasts co-treated with ERKi and p38i appeared differentiated but failed to fuse, leading to the complete absence of multinucleated myofibers, reinforcing the notion that p38 is essential for a transition from differentiation to fusion and maturation (Fig.
5A,B). Immunofluorescence staining validated that the fraction of MyoG-positive cells remained low for p38i-treated cells and increased in cultures co-treated with p38i and ERKi, indicating that co-inhibition of p38 and ERK1/2 leads to bona fide differentiation (Figs.
5C and
EV2; Movie
EV3). Moreover, cells started fusing once the p38 inhibitor was washed out, indicating that they were stalled at a “ready to fuse” state (Fig.
EV5).
However, it was not clear whether the differentiation process was altered with respect to ERKi-treated cells. Thus, we quantitatively described the differentiation process of cells co-treated with p38i and ERKi by applying our motility and actin models trained with proliferating and differentiating cells. As a control, we validated that the differentiation score profile of p38i-treated cells resembled that of proliferating cells treated with DMSO alone (Fig.
5D,E; Appendix Figs.
S14 and
S15). The motility classifier showed that the differentiation profile of p38i-ERKi-treated cells followed a trend strikingly similar to the one obtained for ERKi-treated cells and specifically included the gradual transition at the critical time of 7.5–14.5 h (Fig.
5D). The mean actin classifier score maxed at a value of ~0.5, which was not sufficient to distinguish undifferentiated from differentiated cells, suggesting that p38i alters actin dynamics (Fig.
5E). Nevertheless, the actin model predicted a monotonically increasing trend, suggesting that the cells were becoming more differentiated over time.
To test the model’s prediction that co-inhibition of p38 and ERK1/2 leads to properly differentiated cells that are ready to fuse and that actin dynamics are altered under these conditions, we acquired mass spectrometry data. We extracted proteins from primary myoblasts treated with ERKi, ERKi+p38i, p38i, and DMSO control for 24 h and used mass spectrometry-based proteomics for unbiased, high-throughput identification of differentially expressed proteins. In total 4319 proteins were identified across all experiments. Unsupervised hierarchical clustering analysis showed that p38i-treated cultures cluster with the DMSO proliferation control, consistent with previous studies showing that p38 inhibition maintains myoblasts in the proliferative undifferentiated state (Fig.
6A). The analysis also showed that ERK1/2 and p38 co-inhibited cultures cluster with the differentiating ERKi-treated cultures, suggesting that differentiation occurs under co-inhibition conditions (Fig.
6A).
The actin cytoskeleton is essential for fusion and for the subsequent formation of sarcomeres (Abmayr and Pavlath,
2012; Deng et al,
2015; Gruenbaum-Cohen et al,
2012; Onel and Renkawitz-Pohl,
2009; Richardson et al,
2008). Gene Ontology (GO) annotation enrichment analysis showed a significant upregulation of proteins involved in the process of actin cytoskeleton organization and downregulation in the process of skeletal muscle contraction and maturation in the ERKi-p38i co-inhibition cultures compared to ERKi alone, consistent with an arrest prior to fusion and maturation (Fig.
6B). These results supported our conclusion that co-inhibited cells arrest at terminal differentiation before fusion and maturation and explain the deviation in the differentiation score of the actin, but not the motility classifier, upon ERKi-p38i co-inhibition.
Strikingly, only 18 out of 504 proteins that were significantly differentially expressed between the differentiated (ERKi) and the differentiated unfused cells (ERKi-p38i) were completely missing from the differentiated unfused cells (Fig.
6C). One of these being Myomixer, which is known to be essential for fusion (Zhang et al,
2017; Quinn et al,
2017; Bi et al,
2017). These results suggest that cell differentiation and fusion can be uncoupled and that p38 is essential for the transition from terminal differentiation to fusion and maturation, supporting the notion that differentiation and fusion are coordinated processes that occur in tandem.