Dynamic Regulation of JAK-STAT Signaling Through the Prolactin Receptor Predicted by Computational Modeling

Abstract

Introduction: The expansion of insulin-producing beta cells during pregnancy is critical to maintain glucose homeostasis in the face of increasing insulin resistance. Prolactin receptor (PRLR) signaling is one of the primary mediators of beta cell expansion during pregnancy, and loss of PRLR signaling results in reduced beta cell mass and gestational diabetes. Harnessing the proliferative potential of prolactin signaling to expand beta cell mass outside of the context of pregnancy requires quantitative understanding of the signaling at the molecular level.

Methods: A mechanistic computational model was constructed to describe prolactin-mediated JAK-STAT signaling in pancreatic beta cells. The effect of different regulatory modules was explored through ensemble modeling. A Bayesian approach for likelihood estimation was used to fit the model to experimental data from the literature.

Results: Including receptor upregulation, with either inhibition by SOCS proteins, receptor internalization, or both, allowed the model to match experimental results for INS-1 cells treated with prolactin. The model predicts that faster dimerization and nuclear import rates of STAT5B compared to STAT5A can explain the higher STAT5B nuclear translocation. The model was used to predict the dose response of STAT5B translocation in rat primary beta cells treated with prolactin and reveal possible strategies to modulate STAT5 signaling.

Conclusions: JAK-STAT signaling must be tightly controlled to obtain the biphasic response in STAT5 activation seen experimentally. Receptor up-regulation, combined with SOCS inhibition, receptor internalization, or both is required to match experimental data. Modulating reactions upstream in the signaling can enhance STAT5 activation to increase beta cell survival.

Keywords: Beta cell biology; Ensemble modeling; Feedback control; Intracellular signaling.

Figure 1

Model schematic of JAK-STAT signaling in pancreatic beta cells. PRL binds to the PRLR:JAK2 complex (RJ), which induces receptor dimerization and activation by JAK2 kinase activity. The activated receptor PRL:RJ2* phosphorylates STAT5, which dimerizes and transports into the nucleus, where it promotes transcription of target genes. Phosphatases attenuate the signaling at the membrane (SHP-2), in the cytosol (PPX), and in the nucleus (PPN). Signaling modules for ensemble modeling include (a) STAT5-induced SOCS negative feedback, (b) STAT5-induced receptor up-regulation, and (c) ligand-induced receptor internalization. Green indicates positive feedback; red indicates inhibition of signaling. ECM, extracellular matrix.

Figure 2

Ensemble Modeling Predicts the Number of Peaks in STAT5 Phosphorylation. (a) Simulated time courses were classified into three shapes based on the number of peaks in STAT5 phosphorylation over 6 hours of PRL stimulation. (b) Bar chart shows the percentage of Monte Carlo simulations from each model structure that were classified into each shape shown in panel A. Row labels correspond to the inclusion or exclusion of regulatory modules shown in Fig. 1. The data labels in red show the number of simulations that were classified as “Multiple Peaks” for each structure. n = 100,000 simulations per structure (800,000 total). MP Multiple Peaks.

Figure 3

Classification of simulations into qualitative shapes. Simulated time course of STAT5 phosphorylation for each shape shows the mean (solid line) and 95% confidence interval (shaded area) of all Monte Carlo simulations (800,000 total) classified into that shape. All shapes are mutually exclusive, that is, all simulations were uniquely assigned to one shape (see Fig. S1 for decision tree). Simulations that did not reach a threshold level of 1% of STAT5 phosphorylated were labeled as “weak activation” and filtered out, n = 436,731.

Figure 4

Breakdown of simulations matching desired shape by structure. The y-axis shows the eight model structures defined by the inclusion or exclusion of the regulatory modules. Horizontal bars show the percentage contribution of each model structure to the 101 simulations that matched the desired shape shown in Fig. 3h.

Figure 5

Parameters Correlated with STAT5 phosphorylation. Pearson correlation between each kinetic parameter or initial value and five quantitative characteristics of the STAT5 phosphorylation time course. (a) Illustration of five characteristics. (b) Activation strength. (c) Negative feedback strength. (d) Positive feedback strength. (e) Time of attenuation. (f) Time of reactivation. Only parameters with statistically significant (p < 0 .05) correlations are shown in the waterfall plots. The five parameters most highly correlated with each characteristic are labeled. RJ initial value of PRLR:JAK2 complex, k6 phosphorylation rate of STAT5, k5 activation rate of JAK2, k4 dimerization rate of PRLR:JAK2 complexes, deg_ratio ratio of degradation rate of ligand-bound receptor complexes to unbound complexes, k2 ligand binding on rate, k12 rate of dephosphorylation of pSTAT5 by cytoplasmic phosphatase, PPX initial value of cytoplasmic phosphatase, k11 binding rate of cytoplasmic phosphatase to pSTAT5, k_3 receptor complex dimerization off rate, k3 receptor complex dimerization on rate. The full list of correlated parameters and their Pearson correlation values are given in Supplemental File S1.

Figure 6

Model calibration. Model predictions for (a) Phosphorylated JAK2, normalized to the 10 min time point, (b) Phosphorylated STAT5A, normalized to the 30 min time point, (c) Phosphorylated STAT5B, normalized to the 30-min time point, (d) Ratio of nuclear to cytosolic STAT5A and STAT5B, and (d) Fold change of Bcl-xL. Lines show mean value of model predictions with shading indicating the standard deviation across the 1000 parameter sets from the posterior distribution. Squares show experimental data points from Brelje et al. for panels A, B, C, and D or from Fujinaka et al. for panel E. Error bars are included for experimental data points that had error bars shown in the previously published work. All experimental data are for INS-1 cells treated with PRL at 200 ng/mL. Thirty-three parameters were fit simultaneously to the six data sets using a Bayesian likelihood estimation approach. Dark blue, STAT5A; light blue, STAT5B in (d).

Figure 7

Dose Response predictions. (a) Model predicted time course of STAT5B import into the nucleus under various concentrations of PRL ligand, simulated for 60 minutes. The red dotted line emphasizes the values at the 30 minute time point, which is plotted in the bar chart in panel B. (b) Model predicted dose response data for 30-min timepoint (blue) compared to experimental data from Brelje et al. treating rat primary beta cells with PRL (grey). Values are normalized to the amount of STAT5B in the nucleus with no PRL stimulation (0 ng/mL dose). Error bars for model predictions show standard deviation of predictions across the 1,000 posterior parameter sets.

Figure 8

Model perturbations. (a) The effect of varying the initial ligand-binding rate k2 and the cytosolic phosphatase dephosphorylation rate k12 between 0.1-fold and 10-fold of the fitted parameter values. (b) Varying the initial values of the receptor-JAK2 complex RJ and the cytosolic phosphatase PPX between 0.1- and 10-fold of the fitted values. Coloring of the heat map indicates the initial peak in the STAT5B cytoplasm to nucleus ratio averaged across the 1000 posterior parameter sets.