Optogenetic actuator/ERK biosensor circuits identify MAPK network nodes that shape ERK dynamics

Combining single-cell measurements of ERK activity dynamics with perturbations provides insights into the MAPK network topology. We built circuits consisting of an optogenetic actuator to activate MAPK signaling and an ERK biosensor to measure single-cell ERK dynamics. This allowed us to conduct RNAi screens to investigate the role of 50 MAPK proteins in ERK dynamics. We found that the MAPK network is robust against most node perturbations. We observed that the ERK-RAF and the ERK-RSK2-SOS negative feedbacks operate simultaneously to regulate ERK dynamics. Bypassing the RSK2-mediated feedback, either by direct optogenetic activation of RAS, or by RSK2 perturbation, sensitized ERK dynamics to further perturbations. Similarly, targeting this feedback in a human ErbB2-dependent oncogenic signaling model increased the efficiency of a MEK inhibitor. The RSK2-mediated feedback is thus important for the ability of the MAPK network to produce consistent ERK outputs and its perturbation can enhance the efficiency of MAPK inhibitors.


34
The extracellular signal-regulated kinase (ERK) is part of the mitogen-activated protein  Single-cell biosensor imaging has provided new insights into MAPK signaling that 63 were not accessible with biochemical, population-averaged measurements. This 64 showed that the MAPK network can produce a wide variety of ERK dynamics such as KTR contains an Elk-1 docking domain phosphorylated by ERK (Regot et al. 2014), 139 we hypothesized that it could be negatively affected by optoFGFR-evoked Ca 2+ input 140 (Kim et al. 2014) (Appendix Figure S1E). Consistently, Ionomycin-evoked increase in 141 cytosolic Ca 2+ induced a dip in absence of light stimulation (Appendix Figure S1F). (D, mJ/cm 2 ) as their product to quantify the total energy received per illuminated area. 161 To characterize ERK dynamics, we extracted the amplitude at the maximum of the 162 peak (maxPeak), and the full width at half maximum (FWHM) of the ERK trajectories 163 ( Figure 2B). With increasing light doses, ERK peaks increased both in duration and 164 amplitude, until the latter reached saturation. Based on these observations, we 165 selected 180 mW/cm 2 and 100 ms (D = 18 mJ/cm 2 ) as the minimal light input to 166 generate an ERK transient of maximal amplitude. Using this light dose, we then 167 investigated ERK dynamics in response to multiple light pulses delivered at different 168 intervals ( Figure 2C). All stimulation regimes led to identical maximal ERK amplitude 169 ( Figure EV1A) and adaptation kinetics when optoFGFR input ceased ( Figure EV1B).     . We also included a receptor level inactivation process for EGFR, but not 241 for optoFGFR, to account for EGF-dependent regulatory mechanisms. We used a Bayesian inference approach (Mikelson and Khammash 2020) to infer the model 243 parameters from averaged ERK trajectories in response to sustained low optoFGFR 244 input with or without sustained EGFR pre-stimulation ( Figure 3E). After identification 245 of parameters that allowed the model to capture the training dataset ( Figure 3F), we 246 simulated ERK dynamics evoked by low EGFR input (adaptative, oscillatory ERK 247 dynamics), high EGFR input (adaptative ERK dynamics without oscillation) and 248 sustained high optoFGFR input (sustained ERK dynamics) ( Figure 3G). We observed 249 that our model with a NFB and EGFR inactivation was able to predict ERK dynamics  Thus, these results suggest that optoFGFR lacks receptor-dependent regulatory 269 mechanisms but allows us to investigate the intracellular MAPK feedback structure 270 shaping ERK dynamics. In our model, we used the well-established ERK-RAF NFB.

271
However, several NFBs have been mapped in the MAPK signaling cascade, whose 272 role in shaping ERK dynamics is still unknown and which could also be responsible 273 for the observed oscillatory ERK dynamics.   295 We then explored the network circuitry that shapes optoFGFR-evoked ERK dynamics 296 with an RNA interference (RNAi) screen targeting 50 MAPK signaling nodes. We        Figure EV4F) and z-scored feature values to the negative control ( Figure 5F). We observed more prominent ERK amplitude phenotypes in response to optoSOS input 434 than to optoFGFR input. Some of these phenotypes are shown in Figure 5G. Most 435 prominently, CRAF, ERK2, DUSP4 KD led to a stronger reduction in ERK amplitude 436 than observed with optoFGFR input. RSK2 KD also reduced ERK amplitude, 437 suggesting that it also regulates nodes downstream of RAS. However, RSK2 KD did 438 not decrease ERK adaptation following optoSOS input removal ( Figure EV4G), 439 suggesting that it is not involved in NFB regulation in this system. PP2A KD did not 440 induce increased ERK amplitude or baseline as observed in the optoFGFR system.           While providing the experimental throughput to perturb and analyze ERK dynamics at 641 scale, optoFGFR, that lacks an ectodomain, evoked different ERK dynamics than 642 endogenous RTKs such as FGFR and EGFR ( Figure 3A,B compared to Figure 2F).

649
OptoFGFR therefore provides a simplified system that allowed us to focus on 650 intracellular feedback structures, without confounding receptor level regulations. Our 651 Bayesian inference modeling approach, that is parameter agnostic, could provide 652 simple intuitions about the receptor-level and negative feedback structures that shape 653 ERK dynamics in response to optoFGFR and EGFR inputs. However, even if we had 654 access to many ERK dynamics phenotypes, our modeling approach did not allow us 655 to explore more sophisticated MAPK network topologies such as the presence of two 656 NFBs or multiple node isoforms. We interpreted our data using some of the feedback  The RSK2-mediated feedback can be targeted to potently inhibit oncogenic 716 ErbB2 signaling 717 Our data suggest that the RSK2-mediated NFB is important for MAPK signaling 718 robustness downstream of our prototypic optoFGFR RTK (Figure 6). We found that 719 the RSK2-mediated NFB likely also operates downstream of EGFR and oncogenic 720 ErbB2 signaling in MCF10A cells (Figure 7). In response to EGF stimulation, or ErbB2 721 overexpression, a subset of RSK-inhibited cells displayed wider ERK pulses, 722 suggesting that the RSK2 NFB is also involved in ERK adaptation in this system 723 ( Figure 7G cluster 3, Figure 7K cluster 1 and 2). Further, RSK inhibition led to a high inputs were delivered at specific frequencies and intensities (see below). MCF10A 826 image acquisition was performed at 5-minute time resolution. Growth factor 827 stimulations were done by manually pipetting EGF and bFGF during the experiment. 828 We used mCitrine intensity to quantify the expression level of the optogenetic 829 constructs. However, as mCitrine excitation leads to optoFGFR or optoSOS activation, 830 we acquired one frame with the ERK-KTR-mRuby2, the H2B-miRFP703 and the 831 mCitrine-tagged optoFGFR or optoSOS only at the end of each NIH3T3 experiments.

832
All experiments were carried on at 37°C with 5% CO2.

834
Optogenetic stimulation 835 Light stimulations were delivered with a 470 nm LED light source that was hardware-  Lumencor) (Appendix Figure S1B). All experiments were carried on at 37°C with 5% 872 CO2.

874
Image processing pipeline 875 Nuclear segmentation was done in  the binarized mask obtained with Ilastik using Fiji (Appendix Figure S1C).

903
Quantification of ERK activity 904 We wrote a set of custom R scripts to automatically calculate the ERK-KTR where the variance of the measurement noise ! was estimated from the data.

932
Appendix Table S1 shows all modelled species, their notation used for the equation, 933 as well as the initial values. We assume that in the beginning of the experiment, all 934 species are in the inactive form, reflecting the fact that the cells have been starved.

935
The total concentrations of all species have been normalized to 1. The model 936 equations are shown in Appendix Table S2. The phosphorylation events are modelled 937 with Michaelis-Menten kinetics. The NFB is modelled through the modeling species 938 and its "active" version * which affects the dephosphorylation rate of 939 linearly. The activation, endocytosis, and recycling of the EGF receptor is modelled 940 linearly. The model parameters are described in Appendix Table S3. For the modeling 941 of the two smaller models (without feedback ( Figure EV1E) or without endocytosis 942 ( Figure EV1H)), we set the corresponding parameters ( #)* and !,, ) to zero.

943
For the parameter inference, we used a Nested Sampling algorithm as described in 944 (Mikelson and Khammash 2020). The inference was performed on the ETH High-945 performance Cluster Euler and was done using the parallel implementation on 48 946 cores. The algorithm was run for 24 hours or until the algorithm stopped because the 947 termination criterion -./0 (see (Mikelson and Khammash 2020) for details) was −∞.

948
As prior distributions, we chose for all parameters non-informative log-uniform priors 949 between 10 -5 and 10 5 , except for "#"$ for which we chose a uniform prior on the 950 interval [0, 1] and for for which we chose a log-uniform prior between 10 -5 and 1.

951
Predictive distributions can be found on Figure 3F,G, EV1F,G,I,J.

953
RNAi perturbation screen 954 We used Ingenuity Pathway Analysis (IPA, Qiagen) to select proteins directly  Table S4). We then imported this protein list in STRING (Jensen et 958 al. 2009) to generate an interaction network with a minimum interaction score of 0.4.

959
The final interactome was manually modified to display the protein names to facilitate 960 the readout ( Figure 4A). We targeted these selected proteins with RNA interference,