{"id":"a7007eed-842d-4a3c-b690-9d54f76a993d","arxiv_id":"2607.11820","paper_version":2,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":14,"one_line_summary":"A phenotype-structured PDE model calibrated to melanoma data predicts three stable population behaviours and population-level hysteresis despite reversible single-cell phenotype switching.","lead":"This paper builds a mathematical model of how melanoma cell populations change over time, coupling each cell's MITF level to its growth, death, and movement. The model predicts three possible long-term behaviours—slow growth, faster invasive growth, and rapid oscillatory growth—and shows that a population may not return to its original state even if individual cells can switch back.","discovery_kind":"new_application","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Tri-stability and hysteresis are contingent on the uncalibrated density–stress coupling Eq. (30); robustness to alternative stress–MITF relations is untested.","rationale":"The paper does several things well: subcellular parameters are inferred from published data with MCMC; the homogenisation is carefully derived; the bifurcation analysis is systematic; the spatial model is a clear extension. The central mathematical result—that a density-dependent MITF-repressing feedback can produce bistability and limit cycles—is plausible. My concern is not that the model is internally inconsistent, but that the specific biological prediction of three behaviours and hysteresis rests on an uncalibrated functional response (Eq. 30). The authors themselves flag density as an incomplete stress proxy and dismiss the oscillatory branch as unrealistic. The robustness sweep I propose directly tests whether the claimed tri-stability is an artifact of the chosen Michaelis–Menten form; until that is done, the quantitative status of the central claim is conditional. This matches the reader's weakest-assumption diagnosis and does not alter the verdict.","tokens_in":30389,"tokens_out":11821,"duration_ms":116599,"concrete_test":"Using the provided Mathematica code, recompute the bifurcation diagram of Fig. 13 under three alternative stress couplings: (i) a(M)=exp(-M), (ii) a(M)=1/(1+M^2), (iii) a(M)=1/(1+M)+s with s a small density-independent stress offset. Keep all MAP parameters fixed. Record the existence and κ-ordering of the limit cycle, the INV/PRO steady state, and the PRO/DIF steady state, and the bistable interval. If all three regimes persist in the same κ order with a nonempty bistable window, Eq. (30)'s particular form is not load-bearing; if any regime disappears or the ordering changes, the central tri-stability claim must be treated as conditional on an unvalidated functional form.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Equations (29)-(30) couple population density M(t) to MITF transcription via a(t)=1/(1+M/M*); this feedback is the sole mechanism that converts the reversible single-cell rheostat into the three population-level behaviours and the hysteresis shown in Figs. 9, 13, 14. The spatial analogue (85) makes the same assumption with local density, and the paper concedes (Sect. 5) that the resulting universal PRO/DIF wavefront is a 'direct consequence' of this proxy. The functional form and the half-saturation constant M* are not calibrated against any measurement of the ISR–MITF response; κ is a free parameter whose only constraint is the authors' own rejection of the oscillatory branch as biologically unrealistic. A different monotone stress function (Hill coefficient ≠1, exponential, or with a density-independent stress offset) could shift or destroy the Hopf and double-fold bifurcations, and a non-monotone stress signal could eliminate the invasive steady state entirely. Since the central claim is that the rheostat *coupled to density-dependent stress* generates these behaviours, the unvalidated shape of Eq. (30) is the load-bearing premise.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper develops a multiscale, phenotype-structured PDE model of melanoma population dynamics driven by the MITF rheostat. The authors first reduce a stochastic subcellular model of MITF RNA, protein, and a downstream phenotype variable to an effective advection-diffusion phenotype flux using matched asymptotic expansions (Sect. 2, Appendix B), and then couple the average transcription rate a(t) to total cell density through Eq. (30). The spatially homogeneous population model is reported to exhibit three long-term behaviours — a limit cycle, an INV/PRO steady state, and a PRO/DIF steady state — separated by Hopf and double-fold bifurcations in the contact-inhibition parameter κ (Fig. 13), with hysteresis upon parameter variation (Fig. 14). A radially symmetric spatial extension produces travelling waves whose core behaviour follows the same κ-dependent branches (Sect. 4). Parameters are calibrated by Bayesian inference to published datasets (scRNA-seq, RNA/protein half-life, Ki-67, tumour doubling time, xenograft switching times, TUNEL index).","tokens_in":30791,"tokens_out":6232,"duration_ms":71941,"significance":"The paper has genuine strengths: the multiscale reduction is careful and the moment calculations are verified against stochastic simulations; the inference pipeline is transparent, with MCMC diagnostics and public code; and the conceptual claim that single-cell reversibility need not imply population-level reversibility is clearly demonstrated and is of broad interest. If the three-branch phenomenology is robust, the framework is a valuable contribution to phenotype-structured modelling of cancer. However, the central biological conclusions rest on an uncalibrated density-stress coupling, and the bifurcation structure is reported at one MAP parameter set without uncertainty propagation. These issues need to be addressed before the quantitative claims can be fully accepted.","major_comments":[{"comment":"The density-stress coupling a(t) = 1/(1+M/M*) is the load-bearing mechanism that generates the three long-term behaviours and the hysteresis in Figs. 9, 13, and 14. Yet M* is not constrained by the calibration data: the likelihood (Eqs. 42-45) uses only the low-density exponential phase where a≈1, and M* is scaled out of the dimensionless equations. Thus the threshold density at which ISR-driven MITF repression begins is arbitrary. Likewise κ is not identifiable from the exponential-phase data; its lower bound is imposed by rejecting the oscillatory branch as biologically unrealistic (Sect. 5). The authors should either calibrate or bound M* and κ from data, or at minimum perform a robustness analysis with alternative stress-MITF relations (e.g., Hill coefficient ≠1, exponential, density-independent offset, non-monotone signals) and show whether the double-fold and Hopf bifurcations pers","section":"§3.1, Eq. (30)"},{"comment":"The central bifurcation diagram is computed at the single MAP parameter vector (49). The MCMC posterior samples obtained in Sect. 3.3 are not propagated into the bifurcation analysis, so the reported thresholds (κ≈0.57, 4.7, 10.5) and even the existence of the three branches have no uncertainty quantification. Given the strong parameter correlations shown in Fig. 12 (e.g., ρ_max with ν, Δφ with ν), different posterior draws could shift or eliminate branches. The authors should provide credible intervals for the branch boundaries, or at least a multi-sample overlay of the bifurcation diagram, to substantiate the claim that the model 'admits' three stable behaviours rather than that one MAP parameter set does.","section":"§3.4, Fig. 13"},{"comment":"The paper lists oscillatory dynamics as one of the three stable long-term behaviours in the abstract and conclusions, but later states that the cyclic solutions are 'biologically unrealistic' and uses this to infer a lower bound κ>0.57. This is a post-hoc prior rather than a validation, and it creates tension with the central claim. If the oscillatory branch is not a viable biological prediction, the abstract and conclusions should be revised to present it as a mathematical possibility that is excluded by biological reasoning; alternatively, the authors should identify data that could test this branch.","section":"§5, Discussion"}],"minor_comments":[{"comment":"The median Ki-67 and doubling-time outputs are calibrated quantities, and the posterior predictive checks in Fig. F.3 reuse the same data used for inference. Describing these as 'predictions' is misleading; they should be termed posterior predictive checks or fitted quantities. The genuinely emergent predictions are the three long-term behaviours and hysteresis.","section":"§3.3 and Appendix F"},{"comment":"The constraint that invasive cells comprise <1% of the population during exponential growth is enforced by rejecting MCMC samples, not derived from the model. This should be stated more prominently as an assumption, as it directly shapes the inferred switching time t_s and the resulting Ki-67 prediction.","section":"§3.3, Eq. (47)"},{"comment":"The spatial model introduces two uncalibrated parameters, D_max and ζ. While the authors note this, the wave-speed results in Fig. 18 are therefore illustrative rather than quantitative. A sentence clarifying that no formal inference was attempted for the spatial model would help avoid over-interpretation.","section":"§4"},{"comment":"The right-hand side appears to contain a typo: the death term is rendered as '−ρνρ' and should presumably be '−ν' (or '−νρ'). Please correct.","section":"Eq. (29)"}],"recommendation":"major_revision","confidential_remarks":"This is a well-executed modelling paper with a sound multiscale reduction and a clear conceptual message. The main risk is that the headline three-behaviour bifurcation structure and the associated hysteresis are contingent on an uncalibrated density-stress coupling, Eq. (30), and are presented without uncertainty quantification. I do not think this is grounds for rejection — the framework is worth publishing — but the authors should be asked to address the robustness of the bifurcation structure and to soften claims of 'quantitative' prediction. The paper seems well suited to the journal's mathematical biology scope."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"The real contribution here is the upscaling: the paper starts from a stochastic RNA/protein model, homogenises it into a phenotype flux via a clean asymptotic argument, and couples that to a nonlocal population PDE. The moment calculations in Section 2 are done carefully and checked against simulations; the posterior sampling is honest with diagnostics. That part will be useful to anyone building structured population models from single-cell data. The paper also ships code and data, which is more than most.\n\nThe soft spots are exactly where the reader's report puts them. Equation (30) — density as inverse stress proxy for MITF transcription — is the mechanism that produces the three behaviours and the hysteresis, and it is not calibrated to any stress measurement. The bifurcation diagram in Fig. 13 is computed at the MAP parameters, with no uncertainty bands, even though the MCMC posterior is sitting right there. kappa is a free parameter; the only constraint on it is the authors' own decision that the oscillatory branch is biologically unrealistic. The posterior predictive checks in Fig. F.3 reuse the same Ki-67 and doubling-time data that were used for calibration, so calling them predictions is generous. And the constraint in Eq. (47) rules out invasive cells in the exponential phase, which quietly shapes the inference. These are all addressable — you could propagate the posterior through the bifurcation problem, test alternative monotone or non-monotone stress–density relations, and hold out data — but they have to be addressed before the headline biological claims carry weight.\n\nThe paper is honest about the main weakness: Section 5 admits density is an incomplete stress proxy and that the universal wavefront is a direct consequence. That candour is real. Still, the central result is contingent on an unvalidated functional form, so the biological conclusions are conditional. Where the paper is strongest is the mathematical pipeline itself, not the melanoma conclusions.\n\nWho should read it: mathematically inclined oncology modelers who want a template for deriving phenotype fluxes from subcellular stochastic models. Biologists should read it with caution. It deserves a serious referee — the derivation is novel and the authors have done the heavy lifting on calibration. A good referee should ask for robustness across the posterior, sensitivity to the stress coupling, withheld-data validation, and a proper code release with a commit hash. I'd send it out, expecting revision rather than rejection.","headline":"A careful multiscale model whose tri-stability and hysteresis are real outputs of the equations, but the load-bearing density–stress coupling is uncalibrated, so treat the biological conclusions as conditional.","tokens_in":31264,"tokens_out":1825,"would_cite":true,"duration_ms":22104,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["92D25","35Q92"],"pacs":[],"model":"deepseek-v4-flash","headline":"A phenotype-structured model of the MITF rheostat predicts three stable long-term population behaviors for melanoma, and shows that single-cell phenotype reversibility does not guarantee population-level reversibility.","keywords":["melanoma","MITF rheostat","phenotype switching","structured population model","multiscale PDE","hysteresis","cell plasticity","contact inhibition"],"falsifier":"Measure MITF transcription rate and local cell density in a growing melanoma spheroid or xenograft at multiple time points. If MITF transcription does not decrease monotonically with density, or if a population switched to an invasive phenotype returns to the proliferative/differentiated state when density is reduced (i.e., no hysteresis), the central claim is falsified.","tokens_in":1339,"feed_emoji":"🧬","tokens_out":2248,"duration_ms":50522,"temperature":0.7,"pith_summary":"The paper tries to establish that the MITF rheostat, in which MITF activity drives the switch between proliferative (PRO), invasive (INV), and differentiated (DIF) melanoma cell states, gives rise to emergent population-level dynamics that are more than the sum of single-cell behaviors. By building a multiscale model that couples subcellular MITF dynamics to a phenotypic cell population, the authors find three stable long-term outcomes depending on sensitivity to contact inhibition: a slow-growing PRO/DIF state, a faster-growing state with an invasive core, and a rapidly growing state with an oscillatory core. A key consequence is hysteresis: a population that switches from PRO/DIF to INV/PRO does not return when the original conditions are restored, even though individual cells remain reversible. The authors argue that this cautions against extrapolating single-cell plasticity directly to tumor population dynamics.","feed_headline":"Three stable fates emerge from MITF melanoma model","feed_subtitle":"Density-driven feedback makes single-cell reversibility fail at the population level, producing hysteresis in tumor phenotype.","key_machinery":"The core mechanism is the nonlocal coupling of cell density to MITF transcription, Eq. (30), which acts as a negative feedback from population stress to phenotype. This density proxy drives the bifurcations and hysteresis. The model also uses a phenotype-advection-diffusion flux derived from homogenizing fast subcellular stochastic dynamics of MITF RNA, protein, and a slower phenotype variable; this flux provides a tractable representation of single-cell plasticity within the population PDE.","core_discovery":"The central claim is that coupling the MITF rheostat to cell density—used as a proxy for microenvironmental stress—through Eq. (30), a(t) = 1/(1 + M(t)/M*), produces a bifurcation structure with three stable population behaviors. Numerical and analytical analysis of the phenotype-structured PDE shows that the contact-inhibition sensitivity κ selects between (i) a low-density PRO/DIF steady state, (ii) a higher-density INV/PRO steady state, and (iii) a limit cycle oscillating between INV and PRO phenotypes. The bistability between the two steady states generates hysteresis, demonstrated by a numerical experiment in which lowering κ and then restoring it leaves the population trapped in the IN","pith_inferences":["The qualitative tri-stability and hysteresis likely persist for other monotone decreasing functions a(M) mapping density to MITF transcription, not just the specific Michaelis-Menten form; a robustness check with alternative saturating forms would clarify the model's structural stability.","The hysteresis effect suggests that therapeutic interventions aimed at 'reversing' invasive phenotypes by removing stress (e.g., reducing mechanical confinement or increasing nutrients) may fail unless the population is pushed across the bistable threshold—possibly requiring a transient over-perturbation.","The model's predication that the wavefront is always PRO/DIF implies that spatial sampling of a melanoma edge might systematically miss invasive cells, confounding biopsy-based phenotype assessment; this is an inference not directly drawn in the paper.","The same multiscale homogenization approach could be applied to other rheostat-like transcriptional regulators, where a fast noisy gene-product pair drives a slower phenotype variable, generating population-level PDEs with calibrated parameters from single-cell data."],"forward_implications":["Restoring contact inhibition sensitivity or relieving stress after a population has switched to an invasive phenotype may not restore the original phenotype distribution, implying that treating melanoma may require crossing a hysteresis boundary rather than simply reverting conditions.","The model predicts that the invasive core of a melanoma can persist as a stable steady state even when the overall phenotype distribution is reversible at the single-cell level, offering a mechanism for phenotypic heterogeneity in tumors.","The three stable behaviors—PRO/DIF, INV/PRO, and oscillatory—correspond to different growth speeds and spatial patterns, so the model can be used to interpret radial growth phase dynamics and the transition to vertical growth.","The wavefront in the spatially resolved model is always composed of proliferative and differentiated cells, independent of the core behavior, which follows from using density as the stress proxy.","The oscillatory core, although the paper regards it as biologically unrealistic for most parameters, provides a mathematical prediction that could be tested in engineered systems with very low contact inhibition sensitivity."],"fun_headline_variants":["MITF density model yields three melanoma fates","Density-dependent MITF creates population hysteresis","Melanoma model shows irreversible phenotype switches","Three stable states emerge from MITF rheostat coupling","Phenotype reversibility fails under density feedback"],"cache_read_input_tokens":32512,"weakest_assumption_plain":"The entire feedback structure relies on Eq. (30), which assumes that total cell density M(t) is a faithful proxy for microenvironmental stress and that MITF transcription decreases monotonically with density; if the true stress signal is not monotone in density, or if other cues dominate, the predicted tri-stability and hysteresis may not arise.","fun_headline_variants_meta":{"raw":{"variants":["MITF density model yields three melanoma fates","Density-dependent MITF creates population hysteresis","Melanoma model shows irreversible phenotype switches","Three stable states emerge from MITF rheostat coupling","Phenotype reversibility fails under density feedback"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000126,"raw_usage":{"total_tokens":922,"prompt_tokens":692,"completion_tokens":230,"prompt_tokens_details":{"cached_tokens":256},"prompt_cache_hit_tokens":256,"prompt_cache_miss_tokens":436,"completion_tokens_details":{"reasoning_tokens":160}},"tokens_in":436,"tokens_out":230,"duration_ms":4300,"temperature":1.0,"reasoning_tokens":160,"cache_read_input_tokens":256,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-04T01:39:36.257211+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Measure MITF transcription rate and local cell density in a growing melanoma spheroid or xenograft at multiple time points. If MITF transcription does not decrease monotonically with density, or if a population switched to an invasive phenotype returns to the proliferative/differentiated state when density is reduced (i.e., no hysteresis), the central claim is falsified.","supporting_citations":[],"review_version":2}