Bayesian Decision-Making Shapes Phenotypic Landscapes from Differentiation to Cancer

1. Introduction

Decision-making consists of identifying objectives and choosing actions that optimize outcomes among available alternatives [1]. In a similar spirit, cells can be regarded as decision-making agents that adapt their behaviour through phenotypic evolution in response to local micro-environmental cues, while operating under fixed intrinsic regulatory wiring. Cells sense information from extrinsic variables such as ligands, nutrients, cell density, or mechanical stress, which is encoded into internal states via biophysical and biochemical processes. This information is subsequently decoded into phenotypic responses, creating a feedback loop in which cellular decisions reshape the microenvironment [2,3,4,5,6].

A defining characteristic of cellular decision-making is its stochastic nature. Cells operate in noisy environments and experience fluctuations arising from low copy numbers of molecular species, environmental variability, and intrinsic stochasticity in biochemical reactions [7]. Noise is not purely detrimental; it can play a constructive role in probabilistic decision-making, generate phenotypic heterogeneity, and support bet-hedging strategies [8,9,10]. The filtering, amplification, and utilization of noise in signaling pathways have been extensively discussed in previous theoretical and experimental studies [11,12]. At the same time, phenotypic decisions are regulated through the coordinated interaction of gene regulatory networks[13,14], signaling pathways[15], and metabolic states that integrate external cues with internal dynamics [16]. Extrinsic signals sensed by receptors are transduced into intracellular signaling dynamics, which modulate transcriptional and epigenetic programs to stabilize specific phenotypic states or enable transitions between them [17]. Feedback loops, multistability, and threshold mechanisms within these networks allow cells to commit to discrete fates while retaining flexibility in changing environments[18,19].

Despite significant progress in identifying molecular mechanisms, a complete understanding of phenotypic control at a coarse-grained level remains elusive. In many biological systems, intracellular regulatory interactions are only partially characterized, and the relevant variables are high-dimensional and noisy, often it considered as a “black box” [20,21]. This raises several fundamental questions: How can cellular phenotypic control be characterized in a reduced, probabilistic framework? How do reliable decisions emerge from stochastic intracellular regulation and noisy extracellular signals? How do cells achieve robust phenotypic outcomes despite environmental uncertainty? And how does sensing accuracy—or its imperfection—shape phenotypic precision, stability, and the emergence of pathological states such as malignancy? To asnwer the latter, we postulate our main biological assumption is that differentiated cells can transition back to an undifferentiated, stem-like cell state which can explain carcinogenesis [22]. This reversion has been firstly discovered by Takahashi and Yamanaka [23] where any differentiated cell can be sent back in time to a state of pluripotency by expressing appropriate transcription factors (the Nobel Prize in Medicine 2012).

To address these questions, we adopt Bayesian learning as a candidate organizing principle for cellular phenotypic adaptation. In this framework, cells combine noisy environmental inputs with prior beliefs about their internal state to update their phenotypic configuration over time. These updates subsequently act as priors for future decisions, leading to a dynamical process of learning and adaptation. Exploiting a separation of time scales between fast environmental dynamics and slower internal phenotypic dynamics at the level of probability distributions, we formulate the evolution of phenotypic states as a Bayesian updating process augmented by replication and diffusion. This yields a coarse-grained description of phenotypic evolution that does not rely on explicit molecular mechanisms but instead captures effective statistical relationships between internal and external variables.

The remainder of the paper is organised as follows. In Section 2 we introduce the Bayesian replicator–diffusion model and derive the effective Fokker–Planck description of phenotypic dynamics. Section 3 analyses the resulting phenotypic landscapes, focusing on correlation-driven differentiation (positive or negative), sensing deterioration, and the role of proliferation in tissue homeostasis and carcinogenesis. In Section 4 we discuss the broader implications of this framework and outline directions for future work.

2. Theoretical Developments

2.1. Bayesian Inference and Cell Decision-Making

The decision of cells builds upon two important aspects i.e., $\left(\mathbf{a}\right)$ internal variables (i.e., representing genes, RNA molecules, translational proteins, metabolites, receptors, phenotypes etc.) and the $\left(\mathbf{2}\right)$ external variables (i.e., ligands, chemicals, nutrients, cellular density or stress fields etc.). Cells usually receive the information from the local microenvironment, which latter deciphered in the language of internal variables. We define the external variables as $Y\subset \mathbb{R}$ which does not evolve over time due to fast time scale and acts as an obervation. The internal variables are denoted by $X\subset \mathbb{R}$. In this paper, we focus solely on the scalar case. If the cell acts as a Bayesian decision-maker [24], then we can write the posterior probability distribution of internal variables $P\left(X\mid Y\right)$ at time t as

$$P\left(X+\delta X,t+\tau \right)={\displaystyle \frac{P\left(Y\mid X\right)P\left(X,t\right)}{P\left(Y\right)}}$$

(1)

Table 1. Table of symbols used in the Materials and Methods section.

Symbol	Description
X	Internal/phenotypic random variable
Y	External/environmental random variable
x	Internal/phenotypic states realization
y	External/environmental states realization
$P(X\mid Y)$	Posterior probability distribution of X given Y
$P(Y\mid X)$	Likelihood of perceived environmental state given internal state
$P(X,t)$	Percieved probability distribution of X at time t
$P\left(Y\right)$	Marginal percieved probability distribution of Y
$q\left(Y\right)$	Real distribution of microenvironment Y
$q(X,Y)$	True joint distribution between X and Y
$q(Y\mid X)$	True conditional distribution of Y given X
$\tau$	Relaxation time for adaptation
$\delta X$	Phenotypic displacement due to adaptation
$v\left(X\right)$	Drift velocity of the phenotypes during adaptation
D	Diffusion constant ($\sqrt{{\sigma }_X^{\prime }}$; from physical diffusion noise)
${\sigma }_X^{\prime }$	Standard deviation of noise from physical diffusion
${\sigma }_X$	Standard deviation of noise from imperfect decision-making
$\mu \left(X\right)$	Proliferation (or growth) rate
$\overline{\mu }$	Average proliferation rate $\overline{\mu }=\int \mu \left(x\right)P\left(x\right)dx$
$f\left(X\right)$	Deterministic term of environmental relation $Y=f\left(X\right)+\eta$
$\eta$	Gaussian distributed microenvironmental noise
$\rho$	Intrinsic-extrinsic $X-Y$ correlation coefficient
$\overline{Y}$	Expected Y over true distribution of Y
${\mu }_Y$	Expected Y over percieved distribution of Y
${\mu }_X$	Expected X over the distribution of X
$z_X=(X-{\mu }_X)/{\sigma }_X$	Standardized internal variable
$z_Y=(Y-{\mu }_Y)/{\sigma }_Y$	Standardized external variable
${\overline{z}}_Y$	Expected $z_Y$ under the true conditional $q(Y\mid X)$
$S_X=\langle X\rangle /{\sigma }_X={\mu }_X/{\sigma }_X$	Signal-to-noise ratio (SNR) for internal variable
$P_{\mathrm{ss}}\left(X\right)$	Steady-state probability distribution of phenotypes
Z	Normalization constant (partition function)

Table 1. Table of symbols used in the Materials and Methods section.

Symbol	Description
X	Internal/phenotypic random variable
Y	External/environmental random variable
x	Internal/phenotypic states realization
y	External/environmental states realization
$P(X\mid Y)$	Posterior probability distribution of X given Y
$P(Y\mid X)$	Likelihood of perceived environmental state given internal state
$P(X,t)$	Percieved probability distribution of X at time t
$P\left(Y\right)$	Marginal percieved probability distribution of Y
$q\left(Y\right)$	Real distribution of microenvironment Y
$q(X,Y)$	True joint distribution between X and Y
$q(Y\mid X)$	True conditional distribution of Y given X
$\tau$	Relaxation time for adaptation
$\delta X$	Phenotypic displacement due to adaptation
$v\left(X\right)$	Drift velocity of the phenotypes during adaptation
D	Diffusion constant ($\sqrt{{\sigma }_X^{\prime }}$; from physical diffusion noise)
${\sigma }_X^{\prime }$	Standard deviation of noise from physical diffusion
${\sigma }_X$	Standard deviation of noise from imperfect decision-making
$\mu \left(X\right)$	Proliferation (or growth) rate
$\overline{\mu }$	Average proliferation rate $\overline{\mu }=\int \mu \left(x\right)P\left(x\right)dx$
$f\left(X\right)$	Deterministic term of environmental relation $Y=f\left(X\right)+\eta$
$\eta$	Gaussian distributed microenvironmental noise
$\rho$	Intrinsic-extrinsic $X-Y$ correlation coefficient
$\overline{Y}$	Expected Y over true distribution of Y
${\mu }_Y$	Expected Y over percieved distribution of Y
${\mu }_X$	Expected X over the distribution of X
$z_X=(X-{\mu }_X)/{\sigma }_X$	Standardized internal variable
$z_Y=(Y-{\mu }_Y)/{\sigma }_Y$	Standardized external variable
${\overline{z}}_Y$	Expected $z_Y$ under the true conditional $q(Y\mid X)$
$S_X=\langle X\rangle /{\sigma }_X={\mu }_X/{\sigma }_X$	Signal-to-noise ratio (SNR) for internal variable
$P_{\mathrm{ss}}\left(X\right)$	Steady-state probability distribution of phenotypes
Z	Normalization constant (partition function)

Here, the likelihood i.e., $P\left(Y\mid X\right)$ represents the information sampled by the cell percieved in the local microenvironment which is a percieved state or believed state sensed by the cell. Further the conditional term has been multiplied with the ratio of the probability distribution of internal variables $P\left(X,t\right)$ percieved by the cell and the probability distribution of the believed environment state without sensing $P\left(Y\right)$. Interestingly, the likelihood $P\left(Y\mid X\right)$ can be associated with biological processes of microenvironmental sensing. Cells perceive their surroundings by the engagement of different biochemical processes, like polymerizing pseudopodia, translocating receptor molecules or modifying its cytoskeleton according to mechanical signals [25,26]. In this way, they estimate the aforementioned empirical likelihood. Please note that the real distribution of the environment is defined by $q\left(Y\right)$ and the corresponding joint distribution between internal states and real environment is denoted as $q(X,Y)$ which can be further written as

$$q(X,Y,t)=P\left(X,t\right)q(Y\mid X)$$

(2)

Here we assumed that the cell state distribution is unique and define as $P\left(X,t\right)=q(X,t)$. At the same time, cellular adaptation with the microenvironment helps the cell to jointly cope with the cellular state and its corresponding microenvironment. Expanding the Equation (1)

$${\displaystyle \frac{\partial P\left(X,t\right)}{\partial t}}={\displaystyle \frac{1}{\tau }}\left({\displaystyle \frac{P\left(Y\mid X\right)}{P\left(Y\right)}}-1\right)P\left(X,t\right)+\left(\mu \left(X\right)-\overline{\mu }\right)P\left(X,t\right)-v\left(X\right){\displaystyle \frac{\partial P\left(X,t\right)}{\partial X}}+D{\displaystyle \frac{{\partial }^{2}P\left(X,t\right)}{\partial X^2}}$$

(3)

where $\tau$ is the relaxation time which dictates how fast the the phenotype density will evolve over time t and $\delta X$ is the phenotypic displacement due to adaptation. The $v\left(X\right)={\displaystyle \frac{\delta X}{\tau }}$ incorporates the drift of phenotypes making and D is the diffusion constant as $\sqrt{{\sigma }_X^{\prime }}$. The ${\sigma }_X^{\prime }$ is the standard deviation of the noise coming from physical diffusion which is different as the ${\sigma }_X$ is defined as the standard deviation of the noise due to imperfect decision making which has a different origin. The proliferation rate of phenotype density is defined by $\mu \left(X\right)$ and average proliferation rate is given by $\overline{\mu }$. In Equation (3), we additively combine a replicator term and the corresponding diffusion process. Moreover, for simplicity, we assume that intrinsic variable depends on the phenotypic variable in the following non-linear quasi-steady state way:

$$Y=f\left(X\right)+\eta$$

(4)

where the deterministic part is described as $f\left(X\right)$ and $\eta$ is a Gaussian distributed noise which has mean 0 and standard deviation ${\sigma }^2$. This formulation implies the existence of a time scale seperation between Y (i.e, fast in time) and X (i.e, slow over time). In the Gaussian approximation employed below, the correlation coefficient $\rho$ between X and Y is determined by the local slope $f^{\prime }\left(X\right)$, so that $\rho >0$ corresponds to reinforcing intrinsic–extrinsic coupling, whereas $\rho <0$ represents antagonistic or compensatory feedback between phenotype and environment.

2.1.1. Case-I: Cell State Adaptation Limit ( $\rho \ll 1$)

As there exists a time-scale seperation at the distribution level where X can influence Y while not the other way around. Moreover, we assume that the joint distribution $P(X,Y)$ behaves like a bivariate Gauusian distribution, where we can calculate the corresponding marginals $P\left(X\right)$ and $P\left(Y\right)$. The latter will be both Gaussians with means ${\mu }_X,{\mu }_Y$ and standard deviations ${\sigma }_X,{\sigma }_Y$, respectively,The dependence ratio can be defined as

$$r(X,Y):={\displaystyle \frac{P(X,Y)}{P\left(X\right)P\left(Y\right)}}={\displaystyle \frac{P(Y\mid X)}{P\left(Y\right)}}.$$

(5)

Assuming that the cell is closely to its adaptation limit with respect to the environment, implies the limit $\rho \ll 1$. In turn, we accordingly expanded the quantity $r(X,Y)-1$ as

$$r(X,Y)-1={\displaystyle \frac{P\left(Y\right|X)}{P\left(Y\right)}}-1\approx \rho z_X{z}_Y+{\displaystyle \frac{{\rho }^2}{2}}(z_X^2-1)(z_Y^2-1)+\mathcal{O}\left({\rho }^3\right),$$

(6)

where the $z_X={\displaystyle \frac{X-{\mu }_X}{{\sigma }_X}}$ and $z_Y={\displaystyle \frac{Y-{\mu }_Y}{{\sigma }_Y}}$ the standarized intrinsic and extrinsic random variables. Averaging over the real distribution of the environment $q(Y\mid X)$ in Equation (3), one can write the first term in RHS as

$$\int \left({\displaystyle \frac{P\left(Y\right|X)}{P\left(Y\right)}}-1\right)q(Y\mid X)dY=\rho z_X{\overline{z}}_Y+{\displaystyle \frac{{\rho }^2}{2}}(z_X^2-1)({\overline{z}}_Y^2-1),$$

(7)

where we did use Mehler-Hermite (Gaussian Lancaster) expansion series expansion [27,28,29] until 2nd order with respect to X and Y. Other terms after averaging over $q(Y\mid X)$ in the Equation (3) other terms independent of Y so they vanish under the centered Gaussian average, while higher-order cumulants are neglected due to weak-correaltion limit (i.e., $\rho \ll 1$) as they have functional dependecy on X. We define ${\overline{z}}_Y={\displaystyle \frac{\overline{Y}-{\mu }_Y}{{\sigma }_Y}}$ and $\overline{Y}=\int yq(Y=y|X)dy$ the true environmental average. The term ${\overline{z}}_Y$ represents the sensing bias term as it quantifies the difference between true average on the environment and the cell-percieved distributioin average on the environment. Now, the term in Equation (7) can be written as a second-order polynomial in X as

$$h\left(X\right)=a+bX+c{X}^2$$

(8)

with the coeffecients as

$$a={\displaystyle \frac{\rho {\overline{z}}_Y\langle X\rangle }{{\sigma }_X}}+{\displaystyle \frac{{\rho }^2({\overline{z}}_Y^2-1)}{2}}\left({\displaystyle \frac{{\langle X\rangle }^2}{{\sigma }_X^2}}+1\right)=\rho {\overline{z}}_Y{S}_X+{\displaystyle \frac{{\rho }^2({\overline{z}}_Y^2-1)}{2}}\left(1+{\left(S_X\right)}^2\right)$$

(9)

$$b={\displaystyle \frac{\rho {\overline{z}}_Y}{{\sigma }_X}}-{\displaystyle \frac{{\rho }^2({\overline{z}}_Y^2-1)\langle X\rangle }{{\sigma }_X^2}}={\displaystyle \frac{\rho {\overline{z}}_Y}{{\sigma }_X}}-{\displaystyle \frac{{\rho }^2({\overline{z}}_Y^2-1)S_X}{{\sigma }_X}}$$

(10)

$$c={\displaystyle \frac{{\rho }^2({\overline{z}}_Y^2-1)}{2{\sigma }_X^2}}.$$

(11)

We define Signal-to-noise ratio or SNR as ${\displaystyle \frac{\langle X\rangle }{{\sigma }_X}}$ and $\langle X\rangle ={\mu }_X$.

2.1.2. Case-II: Operating Around the Expected Phenotypic State Limit ($z_X\to 0$)

In this case, we focus on cell internal states characterized around an adapted fixed point, where the corrsponding fluctuations around the mean are small. This is realized by the limit of the standarized intrnsic variable $z_X=(X-{\mu }_X)/{\sigma }_X\to 0$. At a biological level, this is more of a homeostatic and environmentally-sensitive system, where the cell reacts to external fluctuations by trying to maintain a stable internal state where the external signals do not prompt any major change in the internal state of the cell but instead trigger some corrective measures that correct the state in a subtle manner. For this kind of cellular event, one can write dependence ratio for the $z_X\to 0$ limit and intrinsic-extrinsic correlation $\rho \in \left(-1,1\right)$, similar to Equation (7), as

$$\begin{array}{cc}\hfill & r(X,Y)-1={\displaystyle \frac{P\left(Y\right|X)}{P\left(Y\right)}}-1\approx r_0\left(z_Y\right)(1+\left[A{S}_X+B{S}_X^2\right]\hfill \\ \hfill & +\left[{\displaystyle \frac{A}{{\sigma }_x}}-{\displaystyle \frac{2B{S}_X}{{\sigma }_x}}\right]X+\left[{\displaystyle \frac{B}{{\sigma }_x^2}}\right]{\left(X\right)}^2)-1+\mathcal{O}\left({\left(\frac{x}{{\sigma }_x}-S_X\right)}^3\right)\hfill \end{array}$$

(12)

After averaging over $q(Y\mid X)$ in the Equation (12), we obtain:

$$\int \left({\displaystyle \frac{P\left(Y\right|X)}{P\left(Y\right)}}-1\right)q(Y\mid X)dY=\tilde{r}+r_0\left(\left[A{S}_X+B{S}_X^2\right]+\left[{\displaystyle \frac{A}{{\sigma }_x}}-{\displaystyle \frac{2B{S}_X}{{\sigma }_x}}\right]X+\left[{\displaystyle \frac{B}{{\sigma }_x^2}}\right]{\left(X\right)}^2\right).$$

(13)

Here the constant A and B are defined by

$$A={\displaystyle \frac{\rho }{1-{\rho }^2}}{\overline{z}}_Y\&B={\displaystyle \frac{{\rho }^2}{2{(1-{\rho }^2)}^2}}({\overline{z}}_Y^2-{(1-{\rho }^2)}^2)$$

(14)

Where $r_0(z_Y,\rho )={\displaystyle \frac{1}{\sqrt{1-{\rho }^2}}}exp(-{\displaystyle \frac{{\rho }^2}{2(1-{\rho }^2)}}{z}_Y^2)\ge 0$ and $\tilde{r}=\int (r_0\left(z_Y\right)-1)q(Y=y\mid X)dy\approx r_0({\overline{z}}_Y,\rho )-1$. The latter can be positive for large $\rho$ and small ${\overline{z}}_Y$. From Equation (13) we can actually write a second-order polynomial similar to the case-I as

$$\tilde{h}\left(X\right)=\tilde{r}+\tilde{A}+\tilde{B}X+\tilde{C}{X}^2$$

(15)

where we absorbed the constants in terms of coeffecients as

$$\begin{array}{cc}\hfill & \tilde{A}=A{S}_X{r}_0+B{S}_X^2{r}_0,\hfill \\ \hfill & \tilde{B}={\displaystyle \frac{A}{{\sigma }_x}}{r}_0-{\displaystyle \frac{2B{S}_X}{{\sigma }_x}}{r}_0,\hfill \\ \hfill & \tilde{C}={\displaystyle \frac{B}{{\sigma }_x^2}}{r}_0.\hfill \end{array}$$

(16)

So, at the limit of perfect sensing (i.e.,${\overline{z}}_Y\approx 0$), all the three coeffecients in the polynomial term becomes

$$\begin{array}{cc}\hfill & \tilde{A}=-{\displaystyle \frac{{\rho }^2}{2}}{S}_X^2{r}_0<0,\hfill \\ \hfill & \tilde{B}=-{\displaystyle \frac{{\rho }^2}{{\sigma }_x}}{S}_X{r}_0,\hfill \\ \hfill & \tilde{C}=-{\displaystyle \frac{{\rho }^2}{2{\sigma }_x^2}}{r}_0<0.\hfill \end{array}$$

(17)

The replicator term encapsulates the population level selection acting on internal cellular states, implemented as a bias of the distribution toward phenotypes with higher proliferation or survival rates. This contribution, unlike the information acquisition represented by the Bayesian update term, reflects a Darwinian amplification of advantageous internal configurations. As we define $(\mu \left(X\right)-\overline{\mu })P(X,t)$ represents differential proliferation or survival of cells as a function of their internal state X. One can expand this term around ${\overline{x}}^t$ using Taylor series as upto 2^nd order as

$$\mu \left(X\right)-\overline{\mu }={\mu }^{\prime }(X-\overline{X})+{\displaystyle \frac{1}{2}}{\mu }^{″}\left({(X-\overline{X})}^2\right).$$

(18)

Interstingly, we did a Mehler–Hermite (Gaussian Lancaster) expansion [27,28,29] for Bayesian adaptation term and Taylor series expansion for proliferation term there exists a common universal 2^nd order polynomial. Now, we can combine the terms together coming from Bayesian adaptation with the proliferation term and can define a kind of effective fitness term

$${\displaystyle \frac{1}{\tau }}{h}_{eff}\left(X\right)+(\mu \left(X\right)-\overline{\mu })=\gamma X^2+\alpha X+\beta .$$

For the weak correlation function case (case-I) and around the mean value limit (case-II), one can rewrite the new coeffecients with replication term as

$$\begin{array}{cc}\mathbf{Case-I}& \mathbf{Case-II}\\ \begin{array}{cc}\hfill h_{eff}& =h\left(x\right)\hfill \\ \hfill \alpha & ={\displaystyle \frac{b}{\tau }}+{\mu }^{\prime }-{\mu }^{\prime \prime }\hfill \\ \hfill \beta & ={\displaystyle \frac{a}{\tau }}-{\mu }^{\prime }\overline{\mu }+{\displaystyle \frac{1}{2}}{\mu }^{\prime \prime }{\overline{\mu }}^2\hfill \\ \hfill \gamma & ={\displaystyle \frac{c}{\tau }}+{\displaystyle \frac{1}{2}}{\mu }^{\prime \prime }\hfill \end{array}& \begin{array}{cc}\hfill h_{eff}& =\tilde{h}\left(x\right)\hfill \\ \hfill \alpha & ={\displaystyle \frac{\tilde{B}}{\tau }}+{\mu }^{\prime }-{\mu }^{\prime \prime }\hfill \\ \hfill \beta & ={\displaystyle \frac{\tilde{r}}{\tau }}+{\displaystyle \frac{\tilde{A}}{\tau }}-{\mu }^{\prime }\overline{\mu }+{\displaystyle \frac{1}{2}}{\mu }^{\prime \prime }{\overline{\mu }}^2\hfill \\ \hfill \gamma & ={\displaystyle \frac{\tilde{C}}{\tau }}+{\displaystyle \frac{1}{2}}{\mu }^{\prime \prime }\hfill \end{array}\end{array}$$

(19)

2.2. Effective Fokker-Planck Equation for Phenotypic Dynamics

Now we calculate an effective force acting on the phenotype which emerges from Bayesian adaptation, proliferation and diffusion processes in the cotext of Fokker-Planck equation (FKE). The FKE helps us to statistically evaluate the stability, variability, and prevalence of each phenotype within a population. Robust/dominant phenotypes are identified from steady states and landscape minima, while transient probability flows report on phenotypic switching, heterogeneity, and responsiveness due to noise or environmental perturbations. Now, let’s write the FKE as

$${\displaystyle \frac{\partial P\left(X,t\right)}{\partial t}}=\left(\gamma X^2+\alpha X+\beta \right)P\left(X,t\right)-v\left(X\right){\displaystyle \frac{\partial P\left(X,t\right)}{\partial X}}+D{\displaystyle \frac{{\partial }^{2}P\left(X\right)}{\partial X^2}}$$

(20)

Interestingly, the modified equation above can be written as the Fokker-Planck equation [30] with a negative sign in front of the drift term as Equation (20) is the fitness equation with advection term and diffusion term. So, from the negative sign in the force tells us about the energy minimzation (or maximimsing fitness) approach in the modified Fokker-Planck equation setting. for the effective fitness drift ( considering conservative nature of the force ) where

$$\begin{array}{cc}\hfill & {\displaystyle \frac{\partial K\left(X\right)}{\partial X}}=\gamma X^2+\alpha X+\beta ,\hfill \\ \hfill & U\left(X\right)=-\int K\left(X\right)dX+C\hfill \end{array}$$

(21)

So, this tell us precisely that effective force $K\left(X\right)$ and effective potential function $U\left(X\right)$ in the system follows as

$$\begin{array}{cc}\hfill & K\left(X\right)={\displaystyle \frac{\gamma }{3}}{X}^3+{\displaystyle \frac{\alpha }{2}}{X}^2+\beta X+\tilde{C},\hfill \\ \hfill & U\left(X\right)=-\left({\displaystyle \frac{\gamma }{12}}{X}^4+{\displaystyle \frac{\alpha }{6}}{X}^3+\beta {\left(X\right)}^2+\tilde{C}X\right)+C.\hfill \end{array}$$

(22)

Please note that, the potential landscape $U\left(X\right)$ whose curvature and asymmetry are determined by $\alpha$ and $\gamma$. The cubic term $(\propto \alpha )$ introduces asymmetry while the quartic term $(\propto \gamma )$ sets the confinement strength. Although, for stability requirement, one typically requires $(\gamma >0)$ so that $U\left(X\right)\to +\infty \mathrm{as}\left|X\right|\to \infty$. Now, the steady state distribution of the system follows for conservative force ($\tilde{C}=0$)as

$$\begin{array}{cc}\hfill & P_{\mathrm{ss}}\left(x\right)={\displaystyle \frac{1}{Z}}exp\left[{\displaystyle \frac{1}{D}}\left({\displaystyle \frac{\gamma }{12}}{X}^4+{\displaystyle \frac{\alpha }{6}}{X}^3+{\displaystyle \frac{\beta }{2}}{X}^2\right)\right],\hfill \\ \hfill & Z={\int }_{-\infty }^{\infty }{e}^{-U\left(X\right)/D}dX.\hfill \end{array}$$

(23)

So, the effecitive potential function is $U_{eff}$ which is defined as

$$U_{eff}=-{\displaystyle \frac{1}{D}}\left({\displaystyle \frac{\gamma }{12}}{X}^4+{\displaystyle \frac{\alpha }{6}}{X}^3+{\displaystyle \frac{\beta }{2}}{X}^2\right)$$

(24)

From our Fokker-Planck equation of effective fitness (with both conservative and non-conservative part) one can write an overdamped Langevin’s dynamics of the phenotype as

$${\displaystyle \frac{dX}{dt}}={\displaystyle \frac{\gamma }{3}}{X}^3+{\displaystyle \frac{\alpha }{2}}{X}^2+\beta X+\tilde{C}+\xi \left(t\right),$$

(25)

where $\xi \left(t\right)$ is a white or Gaussian noise where the mean of the noise is 0 while the standard deviation is defined by ${\sigma }^2$.

3. Results

3.1. Bayesian Adaptation Shapes Phenotypic Fitness Landscapes

We begin by analysing the deterministic structure of the phenotypic dynamics generated by the effective fitness function derived in Section 2. In the absence of stochastic fluctuations, the overdamped Langevin equation governing the evolution of the phenotypic variable X reads

$${\displaystyle \frac{dX}{dt}}=K\left(X\right)={\displaystyle \frac{\gamma }{3}}{X}^3+{\displaystyle \frac{\alpha }{2}}{X}^2+\beta X,$$

(26)

where the coefficients $(\alpha ,\beta ,\gamma )$ arise from the combined effects of Bayesian adaptation and proliferation (Equation (19)). Equation (26) defines an effective force that can be written as the gradient of an emergent potential $U_{\mathrm{eff}}\left(X\right)$,

$$K\left(X\right)=-{\displaystyle \frac{d{U}_{\mathrm{eff}}}{dX}},U_{\mathrm{eff}}\left(X\right)=-\left({\displaystyle \frac{\gamma }{12}}{X}^4+{\displaystyle \frac{\alpha }{6}}{X}^3+{\displaystyle \frac{\beta }{2}}{X}^2\right),$$

(27)

up to an additive constant. The structure of $U_{\mathrm{eff}}\left(X\right)$ determines the stability, multistability, and large-amplitude behaviour of phenotypic dynamics.

Fixed points and stability.

The stationary phenotypic states are given by the roots of $K\left(X\right)=0$, which yields

$$X^{*}=0,X_{\pm }^{*}={\displaystyle \frac{-3\alpha \pm \sqrt{\Delta }}{4\gamma }},\Delta =9{\alpha }^2-48\gamma \beta .$$

(28)

The linear stability of a fixed point $X^{*}$ is determined by the sign of

$$K^{\prime }\left(X^{*}\right)=\gamma {\left(X^{*}\right)}^2+\alpha X^{*}+\beta .$$

(29)

A fixed point is stable if $K^{\prime }\left(X^{*}\right)<0$ and unstable if $K^{\prime }\left(X^{*}\right)>0$. The sign of $\beta$ controls the local stability of the origin, while $\gamma$ determines whether the potential is confining at large $\left|X\right|$.

Classification of phenotypic landscapes.

Depending on the signs of $(\beta ,\gamma )$ and on the discriminant $\Delta$, the effective potential admits a small number of qualitatively distinct configurations. We identify five canonical phenotypic landscapes, illustrated in Figure 1:

• Phenotypic fixation (homeostatic monostability):$\beta <0$ and $\gamma <0$. The potential has a single deep minimum at $X=0$ and all trajectories relax toward a stable homeostatic phenotype.
• Phenotypic switch (bistable landscape):$\gamma <0$ and $\Delta >0$ with $\beta <0$. Two stable minima coexist and are separated by an unstable fixed point, enabling decision-like switching between discrete phenotypic states.
• Critical switch (marginal stability):$\gamma \approx 0$ or $\Delta \approx 0$. One stable and one unstable fixed point coexist, corresponding to a near-critical landscape at the boundary between bistability and monostability.
• Phenotypic relaxation (canalizing monostability):$\gamma >0$, $\beta >0$, and $\Delta >0$. The landscape possesses a single stable minimum, but nearby unstable fixed points generate slow relaxation and enhanced sensitivity to noise.
• Phenotypic explosion (runaway landscape):$\beta >0$ and $\gamma >0$ with $\Delta <0$. The origin is unstable and the potential fails to confine the dynamics at large amplitudes, leading to runaway phenotypic amplification.

These regimes exhaust all possible qualitative behaviours of the cubic force in Equation (26). Importantly, they arise without assuming any specific molecular regulatory mechanisms: the landscape structure is fully determined by coarse-grained parameters encoding information flow, noise, and proliferation. The Table 2 summarizes the above.

Role of the coefficients.

The three coefficients $(\alpha ,\beta ,\gamma )$ admit clear mathematical and biological interpretations. The linear term $\beta$ controls the local stability of the reference phenotype $X=0$. The quadratic term $\alpha$ introduces asymmetry and bias in phenotypic space, reflecting directional effects of sensing or regulation. The cubic term $\gamma$ determines nonlinear saturation at large amplitude $\left|X\right|$ and thus governs whether the landscape is confining or runaway. Subsequent sections relate these coefficients explicitly to correlation, sensing bias, phenotypic fidelity, and proliferation, and use this classification to interpret differentiation, tissue homeostasis, and carcinogenesis.

3.2. Correlation-Driven Emergence of Cell Fate Decisions: From Stem to Differentiated Cells

In the Gaussian setting considered here, the correlation coefficient $\rho$ between the internal state X and the environmental state Y is directly related to the mutual information,

$$I(X;Y)=-{\displaystyle \frac{1}{2}}ln\left(1-{\rho }^2\right),$$

(30)

so that $\rho$ can be interpreted as an information-theoretic measure of how strongly intrinsic and extrinsic variables are coupled. Large $\rho$ corresponds to high mutual information: microenvironmental fluctuations are tightly tracked by the internal phenotype. Conversely, small $\rho$ implies weak information flow from the environment to the phenotype and, therefore, a more decoupled and internally stabilized state.

To connect this information-theoretic control parameter to cell fate decisions, we construct phase diagrams in the plane spanned by the sensing bias ${\overline{z}}_Y$ and the phenotypic signal-to-noise ratio $S_X$ for different fixed values of $\rho$ (Figure 2). Each point in this plane is classified according to the sign structure of $(\alpha ,\beta ,\gamma ,\Delta )$ and thus belongs to one of the four canonical phenotypic regimes defined in Table 2: phenotypic fixation (homeostatic monostable landscape), phenotypic switch (bistable decision-making landscape), phenotypic relaxation, or phenotypic explosion (runaway landscape).

For correlations close to one ($\rho \approx 0.9$), the diagram is dominated by the explosive regime (red region in Figure 2, left panel). In this limit, internal and external states are almost perfectly correlated and small changes in the microenvironment induce large phenotypic fluctuations. The effective potential becomes shallow or even unconfined, and trajectories are easily driven away from any homeostatic configuration. We interpret this regime as a coarse-grained description of highly plastic stem-like cells, whose phenotype can be strongly driven by the microenvironment.

As $\rho$ decreases (from $0.8$ to $0.65$ in Figure 2), the phase diagram reorganizes. The parameter region associated with phenotypic explosion shrinks, while domains corresponding to phenotypic fixation (homeostatic monostable wells) and phenotypic switch (bistable landscapes) emerges. In these intermediate-correlation regimes, phenotypic regulation of cells with inefficient sensing can still respond sensitively to environmental cues, but improved sensing leads to effective potentials with one or more well-defined minima. Biologically, this corresponds to partially committed progenitor states: cells begin to acquire more stable identities, yet retain a degree of plasticity through bistable decision-making between alternative phenotypic fates.

For low correlations ($\rho \approx 0.3$, rightmost panel in Figure 2), the diagram is dominated by homeostatic and bistable regimes, whereas the explosive regime becomes confined to a small region at high bias. In this limit, the mutual information $I(X;Y)$ is small, the effective potential wells are deep, and phenotypic fluctuations induced by the microenvironment are strongly suppressed. We associate this behaviour with differentiated cells that exhibit robust fate decisions: they either maintain a single, stable phenotype (phenotypic fixation) or switch in a controlled, plastic manner between a discrete set of phenotypic states (bistable decision-making).

Another interesting point is that along the differentiation course, lowering of the correlation $\rho$, not only the potential wells emerges but also deepen. In Figure 3, we plot the potential function with respect to the phenotype in terms of $\rho ,S_X,{\sigma }_X$ and we clearly observe the deepening of the stability regions. To generalize the result, we can analytically estimate the depth of the potential function from the conservative part of the force ($\tilde{C}=0$) as

$$\begin{array}{cc}\hfill \Delta U_{\mathrm{eff}}& =U_{\mathrm{eff}}\left(X_u\right)-U_{\mathrm{eff}}\left(X_s\right)={\displaystyle \frac{{\Delta }^{3/2}}{192{\gamma }^2}}={\displaystyle \frac{{\left(9{\alpha }^2-48\beta \gamma \right)}^{3/2}}{192{\gamma }^2}}\hfill \\ \hfill \Delta U_{\mathrm{eff}}& ={\displaystyle \frac{\sqrt{3}{(1-{\rho }^2)}^5{\sigma }_X^4{\tau }^{2}exp\left(\frac{{\rho }^2{\overline{z}}_Y^2}{1-{\rho }^2}\right)}{16{\rho }^4{\left({({\rho }^2-1)}^2-{\overline{z}}_Y^2\right)}^2}}\hfill \\ \hfill & \times [-{\displaystyle \frac{1}{{(1-{\rho }^2)}^5{\sigma }_X^2{\tau }^2}}({\rho }^{2}exp\left({\displaystyle \frac{{\rho }^2{\overline{z}}_Y^2}{{\rho }^2-1}}\right)[{\rho }^2{S}_X^2{\overline{z}}_Y^4-{({\rho }^2-1)}^2{\overline{z}}_Y^2(2{\rho }^2{S}_X^2-5)\hfill \\ \hfill & +{({\rho }^2-1)}^4({\rho }^2{S}_X^2-8)-14\rho ({\rho }^2-1)S_X{\overline{z}}_Y^3+14\rho {({\rho }^2-1)}^3{S}_X{\overline{z}}_Y]\hfill \\ \hfill & +8{\rho }^2{(1-{\rho }^2)}^{5/2}exp\left({\displaystyle \frac{{\rho }^2{\overline{z}}_Y^2}{2({\rho }^2-1)}}\right)\left({({\rho }^2-1)}^2-{\overline{z}}_Y^2\right)){]}^{3/2}.\hfill \end{array}$$

(31)

For large correlation $\rho \to 1$ , the coupling between X and Y becomes strong and the potential barrier $\Delta U_{\mathrm{eff}}$ decreases, leading to a shallow effective potential. This facilitates frequent state switching and continuous environmental tracking without strong phenotypic preference. On the contrary, for small correlation $\rho$ between the internal state X and the environment Y, indicating weak coupling between sensing and phenotypic redulation. In this regime, the potential barrier $\Delta U_{\mathrm{eff}}$ is large, corresponding to a deep potential well. This reflects a stable phenotype or homeostasis associated with adaptation consolidation and strong memory. This implies that stem cells acquire new identities/phenotypes and these phenotypes become robust in nature.

Taken together, we propose cell differentiation theory as a correlation-driven bifurcation process. High-correlation states represent stem-like cells with large microenvironmentally driven phenotypic fluctuations and explosive landscapes, whereas progressive reduction of $\rho$—and thus of mutual information between X and Y—leads to the emergence and stabilization of homeostatic and bistable landscapes that encode differentiated cell fates. In this view, development and its dysregulation can be understood as deformations of a Bayesian phenotypic fitness landscape controlled by the information flow between intrinsic and extrinsic degrees of freedom.

3.3. Microenvironmental Sensing Deterioration as Pathway to Cancer

In the previous subsection, we interpreted the correlation $\rho$ between internal and external states as an information-theoretic control parameter that drives the emergence of stem-like versus differentiated cell fates. Here, we focus on a different axis of control: the quality of microenvironmental sensing at fixed $\rho$, as summarized in the sensing–fidelity phase diagram of Figure 4. The horizontal axis represents the phenotypic signal-to-noise ratio $S_X$ (phenotypic fidelity), while the vertical axis captures the sensing bias ${\overline{z}}_Y$, which quantifies the mismatch between the true environmental statistics and the distribution perceived by the cell. For $\rho =0.3$, each point in this plane again falls into one of the four canonical regimes of Table 2.

Normal, differentiated cells are expected to operate in a region of high phenotypic fidelity and accurate sensing, corresponding to large $S_X$ and small bias ${\overline{z}}_Y$ (green and yellow regions in Figure 4). In this part of the diagram, the effective fitness landscape is either homeostatic monostable (phenotypic fixation) or bistable (phenotypic switch). Both regimes are consistent with specialized cell types: they maintain a stable phenotype under perturbations, yet can undergo controlled, switch-like changes of state when appropriate signals are received. In terms of the hallmarks of cancer, such cells remain dependent on external growth cues and are sensitive to growth-suppressive signals, since the shape of the landscape is still strongly constrained by the sensed microenvironment [31].

As the sensing bias ${\overline{z}}_Y$ increases, the perceived environment becomes progressively distorted. This deterioration of microenvironmental sensing can be interpreted as the cumulative effect of receptor loss, pathway desensitization, or persistent autocrine signaling [32]. In the phase diagram, increasing ${\overline{z}}_Y$ at fixed, high $S_X$ moves the system vertically upward. Beyond a critical threshold, the corresponding point crosses from the homeostatic or bistable region into the explosive regime, even though $S_X$ remains large. In this regime the effective potential loses its confining minimum and phenotypic trajectories exhibit runaway amplification: the cell no longer settles into a stable fate but instead experiences uncontrolled dynamics in phenotype space.

This transition provides a coarse-grained interpretation of several classical cancer hallmarks. The loss of environmental control over the landscape may correspond to sustained proliferative signalling and insensitivity to growth-suppressive cues: once in the explosive regime, the cell’s phenotypic dynamics are dominated by internal amplification rather than by external regulation. Importantly, our framework shows that this malignant-like behaviour does not require low phenotypic fidelity. Even cells with high $S_X$ - that is, cells capable of expressing a well-defined, reproducible phenotype - can become effectively autonomous if their environmental sensing becomes sufficiently biased. In other words, cancer can be viewed as a disease of biased microenvironmental perception and not a disease of "broken parts".

3.4. Proliferation, Tissue Homeostasis and Carcinogenesis

So far, our analysis has focused on the phenotypic landscape generated by Bayesian adaptation alone. However, in a growing tissue, phenotypes are also subject to differential proliferation and death, which can either stabilise homeostatic states or promote runaway, cancer-like dynamics. In our framework, these effects are captured by the proliferation rate $\mu \left(X\right)$ entering the effective fitness polynomial

$${\displaystyle \frac{1}{\tau }}{h}_{\mathrm{eff}}\left(X\right)+\left(\mu \left(X\right)-\overline{\mu }\right)=\gamma X^2+\alpha X+\beta ,$$

(32)

where $h_{\mathrm{eff}}\left(X\right)$ arises from Bayesian updating and $\overline{\mu }$ is the average proliferation rate. For the environmentally coupled regime (case II), the coefficients including proliferation are given by Equation (19):

$$\begin{array}{ccc}\alpha ={\displaystyle \frac{\tilde{B}}{\tau }}+{\mu }^{\prime }-{\mu }^{\prime \prime },\hfill & \beta ={\displaystyle \frac{\tilde{r}+\tilde{A}}{\tau }}-{\mu }^{\prime }\overline{\mu }+{\displaystyle \frac{1}{2}}{\mu }^{\prime \prime }{\overline{\mu }}^2,\hfill & \gamma ={\displaystyle \frac{\tilde{C}}{\tau }}+{\displaystyle \frac{1}{2}}{\mu }^{\prime \prime },\hfill \end{array}$$

(33)

where ${\mu }^{\prime }={\partial }_X\mu$ and ${\mu }^{\prime \prime }={\partial }_X^2\mu$ are evaluated around the reference phenotype, and $(\tilde{A},\tilde{B},\tilde{C})$ encode the Bayesian adaptation part through $\rho$, ${\overline{z}}_Y$ and $S_X$ (Eqs. (14)–(16)).

The deterministic phenotypic dynamics follow from the effective force

$${\displaystyle \frac{dX}{dt}}=K\left(X\right)={\displaystyle \frac{\gamma }{3}}{X}^3+{\displaystyle \frac{\alpha }{2}}{X}^2+\beta X,$$

(34)

so that the homeostatic state $X^{*}=0$ is locally stable if $K^{\prime }\left(0\right)=\beta <0$ and unstable if $\beta >0$, while the sign of $\gamma$ controls the large-amplitude behaviour of the landscape (Table 2). In particular, the explosive, cancer-like regime is characterised by $\gamma >0$ and $\beta >0$, leading to a single unstable fixed point and a runaway trajectory in phenotype space.

Linear proliferation: the case ${\mu }^{″}=0$.

We first consider the simplest case of a locally linear proliferation profile around the reference phenotype,

$$\mu \left(X\right)\approx \overline{\mu }+{\mu }^{\prime }(X-\overline{X}),{\mu }^{″}=0.$$

(35)

In this case Equation (33) simplifies to

$$\alpha ={\alpha }_{\mathrm{adapt}}+{\mu }^{\prime },\beta ={\beta }_{\mathrm{adapt}}-{\mu }^{\prime }\overline{\mu },\gamma ={\gamma }_{\mathrm{adapt}},$$

(36)

where we have defined the adaptation-only coefficients

$${\alpha }_{\mathrm{adapt}}={\displaystyle \frac{\tilde{B}}{\tau }},{\beta }_{\mathrm{adapt}}={\displaystyle \frac{\tilde{r}+\tilde{A}}{\tau }},{\gamma }_{\mathrm{adapt}}={\displaystyle \frac{\tilde{C}}{\tau }}.$$

(37)

Thus, linear proliferation does not modify the large-amplitude saturation term $\gamma$ (and therefore cannot by itself create or remove the explosive regime), but it tilts the landscape by shifting the linear and quadratic contributions.

The local stability of the homeostatic state is determined by

$$K^{\prime }\left(0\right)=\beta ={\beta }_{\mathrm{adapt}}-{\mu }^{\prime }\overline{\mu }.$$

(38)

For a fixed adaptation landscape, the sign of $\beta$ is controlled by the slope ${\mu }^{\prime }$:

• If ${\beta }_{\mathrm{adapt}}<0$, adaptation alone already stabilises $X^{*}=0$. A range of values of ${\mu }^{\prime }>0$ preserves $\beta <0$, corresponding to robust homeostasis.
• If ${\beta }_{\mathrm{adapt}}>0$, adaptation alone would make $X^{*}=0$ unstable. A suitable choice of ${\mu }^{\prime }<0$ can further enhance instability (make $\beta$ more positive).

Biologically, ${\mu }^{\prime }$ quantifies how proliferation responds to small deviations in phenotype. If cells deviating from the homeostatic phenotype proliferate faster than average (${\mu }^{\prime }>0$ with $\overline{\mu }>0$), the term $-{\mu }^{\prime }\overline{\mu }$ tends to decrease $\beta$ and can partially compensate an unstable adaptation contribution. Conversely, if cells deviating from homeostasis proliferate slower (${\mu }^{\prime }<0$), $\beta$ is increased and the origin can be destabilised. Importantly, because $\gamma$ remains equal to ${\gamma }_{\mathrm{adapt}}$, the presence or absence of a confining well at large $\left|X\right|$ is still determined entirely by sensing properties (correlation $\rho$, bias ${\overline{z}}_Y$ and SNR $S_X$). In the linear case ${\mu }^{″}=0$, proliferation therefore primarily modulates the local stability and bias of the landscape rather than its global confinement.

Healthy tissue homeostasis corresponds, in this language, to a regime where the adaptation part and the proliferation profile combine to yield $\beta <0$ and $\gamma <0$: perturbations of the phenotype decay and the effective potential possesses a single deep minimum (phenotypic fixation) or a controlled bistable structure (phenotypic switch). Carcinogenic transformation in the ${\mu }^{″}=0$ limit requires that the adaptation parameters first drive the system into a regime with ${\gamma }_{\mathrm{adapt}}>0$ (loss of saturation at large $\left|X\right|$), and that the combined effect of adaptation and proliferation yields $\beta >0$. In that case the only fixed point at $X=0$ becomes repelling and the phenotype undergoes unbounded amplification, consistent with sustained proliferative behaviour that has escaped microenvironmental control.

Beyond linear proliferation: curvature ${\mu }^{″}\ne 0$ and growth control.

When ${\mu }^{″}\ne 0$, proliferation no longer acts only as a tilt but also reshapes the large-amplitude behaviour of the landscape through $\gamma$:

$$\gamma ={\gamma }_{\mathrm{adapt}}+{\displaystyle \frac{1}{2}}{\mu }^{″}.$$

(39)

A negative curvature ${\mu }^{″}<0$ corresponds to saturating growth (proliferation is reduced for extreme phenotypes) and makes $\gamma$ more negative, deepening the confining well and protecting against explosive dynamics. This is consistent with growth control mechanisms in healthy tissues, such as contact inhibition [33] and resource limitation [34], which penalise phenotypes that stray too far from homeostasis.

Conversely, a positive curvature ${\mu }^{″}>0$ corresponds to superlinear growth (proliferation increases for extreme phenotypes) and makes $\gamma$ more positive. If the adaptation part already brings the system close to the boundary between confinement and explosion, even a moderate positive ${\mu }^{\prime \prime }$ can flip $\gamma$ to positive and push the landscape into the explosive regime. In this situation, proliferation no longer counteracts deviations from homeostasis but amplifies them, in line with the hallmarks of cancer such as sustained proliferative signalling and insensitivity to growth-suppressive cues.

In summary, proliferation shapes phenotypic landscapes at two complementary levels. With ${\mu }^{″}=0$, linear dependence of proliferation on phenotype tilts the landscape and controls the local stability of homeostatic states through $\beta$, modulating whether a given adaptive landscape supports fixation or instability. When ${\mu }^{″}\ne 0$, the curvature of the proliferation profile directly alters the nonlinear saturation encoded in $\gamma$, either reinforcing tissue homeostasis (for ${\mu }^{″}<0$) or facilitating carcinogenesis (for ${\mu }^{″}>0$) by converting marginally confined landscapes into explosive, cancer-like ones.

3.5. Negative Correlation-Driven Bifurcation: The Robustness-Plasticity Trade-Off

So far, our analysis has focused primarily on positive correlations between intrinsic phenotypic states and extrinsic microenvironmental variables. In this regime, the phase diagrams in the $({\overline{z}}_Y,S_X)$ plane display extended regions of bistability and phenotypic explosion, enabling both phenotypic plasticity and malignant runaway. Here, we investigate the complementary case of negative correlations, motivated by the numerical phase diagrams shown in Figure 5.

Our simulations reveal two robust qualitative effects of negative correlation. First, the region associated with phenotypic explosion (red area) is strongly reduced compared to the positive correlation case. Second, and more subtly, the bistable decision-making region (green area) is also significantly diminished. As a consequence, negative correlations promote global stability and suppress runaway dynamics, but at the cost of reducing the parameter space supporting bistable phenotypic switching.

This behaviour can be understood by revisiting the origin and interpretation of the correlation parameter $\rho$. In our model, environmental variable is described by the stochastic relation

$$Y=f\left(X\right)+\eta \approx f\left(X_0\right)+f^{\prime }\left(X_0\right)(X-X_0)+\eta ,$$

(40)

where $f\left(X\right)$ captures the deterministic coupling between phenotype and environment, and $\eta$ is Gaussian noise. Here, the correlation coefficient between X and Y is therefore proportional to the local slope $\rho \propto f^{\prime }\left(X_0\right)$. Positive correlation corresponds to $f^{\prime }\left(X_0\right)>0$, meaning that increases in the phenotype are reinforced by environmental cues, whereas negative correlation corresponds to $f^{\prime }\left(X_0\right)<0$, implying an antagonistic or compensatory coupling between intrinsic and extrinsic variables.

Negative correlation thus acts as an effective negative feedback between phenotype and environment. Such feedback suppresses large phenotypic excursions and strongly stabilises the effective potential, explaining the reduced extent of the explosive regime observed in the phase diagrams. However, the same stabilising mechanism also suppresses bistability. Bistable decision-making relies on the coexistence of two competing stable attractors separated by an energy barrier, which in turn requires a sufficiently strong positive coupling between intrinsic state and environmental signal. When $\rho <0$, this positive reinforcement is weakened or reversed, narrowing the parameter range over which bistable phenotypic landscapes can exist.

From a biological perspective, this can alternative differentiation course via negative intrinsic–extrinsic correlations bifurcation cannot be the dominant regime in physiological tissues. Many well-documented cellular processes, such as epithelial–mesenchymal transition (EMT/MET) [35] or immune activation [36], rely on bistable or multistable dynamics that enable controlled phenotypic plasticity. The suppression of bistability observed for $\rho <0$ would therefore imply overly rigid phenotypic behaviour, characterizing a subset of physiological tissues. Candidates for such tissues can be the ones that cancer occurances are rare or non-existence, such as heart, skeletal muscle, mature neurons and more [37].

At the same time, the stabilising effect of negative correlation has a potentially important therapeutic implication. In the cancer-like regime, phenotypic explosion corresponds to a loss of confinement in the effective landscape, driven by biased sensing and deregulated proliferation. Introducing or enhancing negative feedback between intrinsic phenotypes and microenvironmental signals—effectively pushing the system toward $\rho <0$—can strongly reduce the explosive region and drive the system back toward a stable, homeostatic phenotype. In this sense, negative correlation provides a theoretical mechanism for phenotypic re-stabilisation rather than phenotypic plasticity.

4. Discussion

In this work, we developed a coarse-grained, information-theoretic framework in which cellular phenotypic adaptation is described as Bayesian decision-making embedded in a replicator–diffusion process. Bayesian updating of noisy environmental cues generates an effective drift term, while replication and diffusion complete an emergent Fokker–Planck description of phenotype dynamics. The resulting effective potential defines a phenotypic fitness landscape that is not assumed a priori but derived from sensing fidelity, environmental bias, correlation between intrinsic and extrinsic variables, and proliferation. Within this landscape, we identified homeostatic fixation, bistable switching, and explosive regimes, and used them to interpret differentiation, tissue homeostasis, and carcinogenesis.

Within the Gaussian setting considered here, the correlation coefficient $\rho$ between intrinsic and extrinsic states plays a central conceptual role. Because $\rho$ directly encodes the mutual information $I(X;Y)$, it acts as an information-theoretic control parameter for the shape of the phenotypic landscape: high correlation yields shallow or unconfined landscapes with high plasticity, whereas decreasing $\rho$ deepens potential wells and stabilises discrete phenotypic fates. In this view, differentiation corresponds to a correlation-driven bifurcation process in which stem-like, highly plastic cells move from an explosive regime towards robustly differentiated states as mutual information decreases. Proliferation modulates this picture at two levels: linear responses in ${\mu }^{\prime }\left(X\right)$ tilt the landscape and change local stability, while curvature ${\mu }^{\prime \prime }\left(X\right)$ alters global confinement, so that superlinear growth can convert marginally stable landscapes into fully explosive ones.

Our numerical exploration of negative intrinsic–extrinsic correlations ($\rho <0$) reveals an interesting asymmetry. When $\rho$ is negative, the region of explosive dynamics in the $({\overline{z}}_Y,S_X)$ phase space becomes strongly suppressed, but the bistable decision-making region is also substantially reduced. Thus, negative correlations enhance overall robustness at the cost of strongly limiting phenotypic plasticity. This is consistent with the interpretation of $\rho$ as the local slope of the environment-phenotype coupling in Equation (4),since $\rho >0$ corresponds to reinforcing feedback between phenotype and environment, whereas $\rho <0$ encodes an antagonistic, compensatory coupling. While such antagonistic coupling may be too rigid to account for the rich bistable phenotypic dynamics observed in physiological processes such as EMT/MET, it suggests a potential strategy for re-stabilising pathological states.

These observations have direct implications for cancer. In our framework, carcinogenesis emerges when deteriorating microenvironmental sensing (increasing ${\overline{z}}_Y$) and loss of growth control (positive curvature ${\mu }^{\prime \prime }$) jointly drive the system into an explosive regime in which the effective potential loses its confining minimum. In this regime, cells become effectively autonomous from growth-regulating cues: their phenotypic trajectories are dominated by internal amplification rather than external regulation, consistent with hallmarks such as sustained proliferative signalling and insensitivity to growth-suppressive cues [31]. Crucially, this malignant-like behaviour does not require low phenotypic fidelity: even cells with high signal-to-noise ratio $S_X$—and therefore precise, reproducible phenotypes—can undergo phenotypic explosion if their environmental sensing becomes sufficiently biased. Therapeutically, our results suggest two complementary routes to restore control: (i) reducing effective sensing bias or introducing antagonistic feedback to shift the effective $\rho$ towards less explosive regimes, and (ii) enforcing saturating or sublinear proliferation (negative ${\mu }^{\prime \prime }$) to deepen confining wells and penalise extreme phenotypes. The former bears similarities with the ideas of differentiation therapy [38] and feedback-restoring therapies [39]. Both approaches aim to push the system back from the explosive region of the phase diagram into homeostatic or controlled bistable regimes.

The framework also provides a natural lens on exploration–exploitation trade-offs in cell fate decision-making [7]. High-correlation, weakly confining landscapes support large phenotypic excursions and multiple accessible fates: stem-like or progenitor cells in this regime “explore” phenotype space, using microenvironmental cues to sample a broad range of states. As correlation decreases and the landscape develops deeper wells, cells progressively “exploit” a reduced subset of phenotypes, committing to robust differentiated fates. Bistable landscapes occupy an intermediate regime, balancing exploration (via occasional switching) with exploitation (via stable wells) and thereby capturing controlled phenotypic plasticity. This interpretation is consistent with observations that multipotent lineages exhibit stronger drift and diffusion in high-dimensional state space than fully committed cells, and that state-dependent stochasticity is essential for maintaining non-equilibrium developmental dynamics.

Beyond biological interpretation, our results suggest concrete ways to inform learning algorithms for cell dynamics, particularly neural stochastic differential equation frameworks such as scDiffEq [40]. In those models, drift and diffusion are represented by neural networks and trained from single-cell data using optimal-transport losses. Our theory provides a family of analytically tractable drift fields—cubic forces derived from Bayesian adaptation and proliferation—that can serve as inductive biases or regularisation targets. Conversely, neural SDE models can be used to infer empirical drift and diffusion fields from high-dimensional data and then project them onto low-dimensional phenotypic coordinates, allowing direct estimation of effective coefficients $(\alpha ,\beta ,\gamma )$ and of the inferred $\rho ,{\overline{z}}_Y,S_X,\mu \left(X\right)$ from data. In this way, our analytic framework and data-driven neural differential equation approaches become mutually constraining: theory guides model architecture and regularisation, while learned dynamics provide empirical tests of the predicted landscape regimes. A similar endeavor, inspired by the ideas of adaptation and learning, was the development of geodesic learning [41].

Several limitations of the present framework point to important directions for future work. First, we assumed Gaussian statistics between intrinsic and extrinsic variables, which allowed us to identify correlation with mutual information in closed form and to derive a simple cubic force. Non-Gaussian or multimodal statistics would generate more rugged, hierarchical landscapes and may be essential to capture complex fate trees and oscillatory dynamics. Second, we restricted ourselves to one-dimensional phenotypic and environmental variables; real single-cell state spaces are high-dimensional, and extending the theory to higher dimensions will be necessary to connect more directly to transcriptomic or epigenetic data and to exploit the full power of neural SDE learning frameworks. Third, we neglected explicit cell–cell interactions and spatial mechanics, which are known to feed back on phenotypes and may qualitatively alter the landscape at the tissue level. Finally, our therapeutic and algorithmic implications are conceptual: translating effective parameters such as $\rho$, ${\overline{z}}_Y$ or ${\mu }^{\prime \prime }$ into concrete molecular interventions or into practical regularisation schemes for learning algorithms will require careful calibration against experimental data. Despite these limitations, the present work offers a compact, interpretable scaffold on which more detailed mechanistic models and data-driven neural differential equation approaches can be systematically built.