# On a Stochastic PDE Model for Epigenetic Dynamics

Pablo Padilla-Longoria  IIMAS, Universidad Nacional Autonoma de Mexico (UNAM)
Ciudad Universitaria, Mexico City, Mexico
<EMAIL_ADDRESS>Jesus Sierra  CIMAT
A.P. 402, 3600, Guanajuato. Gto, Mexico
<EMAIL_ADDRESS>
((Day Month Year); (Day Month Year))

###### Abstract

We model the process of cell fate determination and mutation of the flower
$Arabidopsis-thaliana$ through a stochastic reaction-diffusion system governed
by a potential field and additive noise. The potential mimics the flower’s
epigenetic landscape as defined by Waddington, and the noise represents
environmental fluctuations. We show through numerical simulations that the
system eventually exits the local minima, traversing the epigenetic landscape
in the spatial order that, in many of the realizations, corresponds to the
correct architecture of the flower, that is, following the observed
geometrical features of the meristem. We use the theory of large deviations to
estimate the exit time, characterize the associated invariant measure, and
discuss its phenotypic implications from an evolutionary perspective. We also
investigate an optimal control problem for our system to study the reversal of
epigenetic mutations. This setting can be interpreted not only in the context
of epigenetic aspects of evolution and development (eco-evo-devo), but also
has possible applications in developing epi-drugs for diseases like cancer and
Alzheimer’s.

###### keywords:

Epigenetics; Mutations; Large deviations; Invariant measure; Stochastic PDE;
Stochastic control

(xxxxxxxxxx)

AMS Subject Classification: 35R60, 60H15, 92D10, 37N25, 49N90

##  1 Introduction

In this paper, we study floral architecture development as a prototypical case
of a system in which regulatory networks are responsible for the emergence of
structures and patterns. We are particularly interested in the influence of
epigenetic factors, especially environmental fluctuations, on this process.

There are approximately 250,000 species of flowering plants (Angiosperms). The
organs of the flower in most of them (the only known exception being the
flower $Lacandonia-schismatica$) are organized in four concentric rings
(whorls): sepals, petals, stamens, and carpels (from the outer rim to the
center).

We work with $Arabidopsis$ $thaliana$, the first plant whose complete genome
was sequenced and has been extensively studied.[12] In Ref. alvarez2008floral,
using experimental data, the authors obtained the gene regulatory network
(GRN) that determines the fate of floral organ cells in $Arabidopsis$
$thaliana$. Based on this model, Cortés-Poza and Padilla-Longoria constructed
a system of reaction–diffusion equations governed by a potential field
corresponding to the epigenetic landscape of the flower’s organ formation.[7]
In this work, we want to introduce perturbations due to epigenetic factors,
particularly environmental fluctuations, into the system as additive noise;
our goal is to study this additional mechanism in connection with mutations
observed in $Arabidopsis$ $thaliana$.

We consider an energy landscape with isolated minima. Furthermore, the
deterministic part of the dynamics drives the system, by steepest descent, to
the vicinity of one of these minima, where it remains for a very long time. We
will see that the random perturbations push the system significantly up and
away from this minimum. After some time, the system will manage to escape the
basin of attraction of the minimum it is currently in and find its way toward
the location of another minimum. The cycle can then repeat.

We will use the theory of large deviations (see Ref. freidlin1998random) to
estimate the time of escape from the location of one minimum to another, which
is exponentially long relative to the height of the energy barrier between
these minima measured in units of the amplitude of the random perturbation.
The theory of large deviations also provides information about the paths of
maximum likelihood by which the transitions between minima occur.

On the other hand, we believe our results can be applied to study the
connection between environmental fluctuations and diseases like Alzheimer’s
and cancer. As mentioned in Ref. cancerepig, it is now recognized that
epigenetics plays a role in the development of cancer (carcinogenesis); see
also Ref. esteller2008epigenetics,sharma2010epigenetics. For example, abnormal
epigenetic modifications in specific oncogenes and tumor suppressors genes can
result in uncontrolled cell growth and division that can cause cancer.
Besides, epigenetic alterations in regions of DNA outside of genes can also
give rise to cancer. It is now accepted that the environment and human
behavior are the principal causes of abnormal epigenetic modifications.

Unlike genetic mutations, epimutations are reversible, which gives the
possibility to reverse the epimutations in cancer cells through the so-called
epi-drugs.[20] The targets of epigenetic therapy are the enzymes involved in
epigenetic modifications. This motivates the introduction of the theory of
control into our system: our model can help to find and characterize the
elements (controls) needed for such reversal (see, e.g., Ref.
lu2021mathematical).

The rest of the paper is organized as follows. Section 2 introduces the
stochastic epigenetic model. In Section 3, we present our main results. We
review the necessary theory and notation in Section 4 and demonstrate our
results in Section 5. Finally, Section 6 shows our numerical simulations.

##  2 Stochastic Epigenetic Model

We model the epigenetic landscape of the flower as a potential field with four
different basins of attraction, each corresponding to a different flower organ
(sepals, petals, stamens, and carpels). For the size of each basin, we will
use the reciprocal of the number of (dynamically observed) initial conditions
that land in each steady state of the dynamical system (see Ref.
cortes2022variational); this guarantees that equilibrium points that are
reached more often will have larger basins, and conversely, equilibrium points
that are reached fewer times will have a smaller basin. The centers of the
basins are located in $\mathbb{R}^{2}$; see Ref. cortes2022variational.

From an evolutionary perspective, our model provides significant information.
For instance, the fact that the basins of attraction of carpels and stamens
are the largest ones, as well as the most likely end points of realizations of
the solutions, has a natural evolutionary interpretation; that is, in terms of
fitness, the most essential cell fates are those of reproductive organs,
namely, stamens and carpels. Furthermore, it will follow from the theory of
large deviations that the time spent in a basin of attraction is considerably
(in fact exponentially) longer than that corresponding to transitions between
these basins, which is consistent with the theory of punctuated equilibrium
(see Ref. gould2002structure) that accounts for the discontinuity observed in
the fossil records and the ubiquity of “missing links” in the evolutionary
process.

We now define the potential field on the plane $\left(u,v\right)$ determined
by the epigenetic landscape in the following way:

$\displaystyle F\left(u,v\right)=\min$
$\displaystyle\left\\{a_{1}\left[\left(u-u_{1}\right)^{2}+\left(v-v_{1}\right)^{2}\right],a_{2}\left[\left(u-u_{2}\right)^{2}+\left(v-v_{2}\right)^{2}\right],\right.$
$\displaystyle\left.a_{3}\left[\left(u-u_{3}\right)^{2}+\left(v-v_{3}\right)^{2}\right],a_{4}\left[\left(u-u_{4}\right)^{2}+\left(v-v_{4}\right)^{2}\right]\right\\},$
(1)

where $\left(u_{1},v_{1}\right),$ $\left(u_{2},v_{2}\right),$
$\left(u_{3},v_{3}\right),$ and $\left(u_{4},v_{4}\right)$ are the centers of
the basins and $a_{1},a_{2},a_{3},a_{4}$ define the size of each basin. As in
Ref. cortes2022variational, we consider the smooth (mollified) version of $F$.
Moreover, we use a Gaussian filter to smooth out $F$ in our numerical
simulations.

Define $\tilde{f}=-\partial F\left(u,v\right)/\partial u$,
$\tilde{g}=-\partial F\left(u,v\right)/\partial v$. We study the following
stochastic reaction-diffusion system with periodic boundary conditions on
$\mathcal{O}=\left[0,R\right],\ R>0$, and with homogeneous Neumann boundary
conditions on $\mathcal{O}=\left(0,R\right)$:

$\displaystyle\frac{\partial u}{\partial t}\left(t,x\right)=$ $\displaystyle
d_{1}\frac{\partial^{2}u}{\partial
x^{2}}\left(t,x\right)+\tilde{f}\left(u\left(x,t\right),v\left(x,t\right)\right)+\sigma_{1}w_{1}\left(x,t\right),$
$\displaystyle\frac{\partial v}{\partial t}\left(t,x\right)=$ $\displaystyle
d_{2}\frac{\partial^{2}v}{\partial
x^{2}}\left(t,x\right)+\tilde{g}\left(u\left(x,t\right),v\left(x,t\right)\right)+\sigma_{2}w_{2}\left(x,t\right),~{}x\in\mathcal{O},~{}t>0,$
(2)

where $d_{1},d_{2}>0$ are diffusion constants and $\sigma_{1}$, $\sigma_{2}>0$
represent the magnitude of the $noise$. Furthermore, we consider a white noise
perturbation:
$w_{k}\left(t,x\right)=\partial^{2}\hat{w}_{k}\left(t,x\right)/\partial
t\partial x$, $k=1,2$, where $\hat{w}_{k}\left(t,x\right)$ are independent
Brownian sheets. We will focus on the case $\sigma_{1}=\sigma_{2}=\sigma>0$
and write
$u^{\sigma}\left(t,x\right)=\left(u\left(t,x\right),v\left(t,x\right)\right)$
for the corresponding solution of (2).

We call $M\subset\mathbb{R}^{2}$ a phenotype (or cellular type) region of the
epigenetic landscape if $M$ contains a unique center
$\left(u_{k},v_{k}\right)$, $k\in\left\\{1,\ldots,4\right\\}$, $M$ is closed
and convex (hence simply connected), and $\left(\tilde{f},\tilde{g}\right)$
points strictly into $M$ on $\partial M$.

As in Definition 14.5,[24] a closed subset, $\Sigma\subset\mathbb{R}^{2}$, is
an invariant region for the deterministic epigenetic system if any solution
$u^{\sigma=0}\left(t,x\right)$ having all of its boundary and initial values
in $\Sigma$, satisfies $u^{\sigma=0}\left(t,x\right)\in\Sigma$ for all
$x\in\mathcal{O}$ and $t>0$.

Later, we will show (using the theory of large deviations) that our epigenetic
model can generate a mutation process. Our goal afterward is to study a
mechanism (control) capable of reversing such mutation; this control’s
objective will be to return the system from the mutated state caused by
environmental fluctuations. Hence, we look for controls capable of moving the
system state from a mutated (neighborhood of a) basin to a specific non-
mutated one in the epigenetic potential. This approach can produce significant
information for developing therapies and the so-called epi-drugs to treat
diseases like cancer and Alzheimer’s.

From our previous discussion, we have to consider an optimal control problem
with endpoint/state constraints. Unfortunately, this problem is not well
understood up to now. Hence, we will start with a closely related optimal
control problem for which we can show a Pontriagyn-type maximum principle.
Afterward, we will study the case with endpoint/state constraints using set-
valued analysis and additional restrictions.

For our study of optimal control problems for the stochastic epigenetic model,
it is convenient to adopt the semigroup theory approach of Da Prato and
Zabczyk for stochastic PDE.[9] We note that, for the Brownian sheets
$\hat{w}_{k}\left(t,x\right)$, $k=1,2$ in (2), the distributional derivative
$\partial\hat{w}_{k}\left(t,x\right)/dx$ can be identified, up to a constant,
with a cylindrical Wiener process in $L^{2}\left(\mathbb{T}\right)$; see
Section 4.1.5 in Ref. da2014stochastic.

Let $H$ and $V$ be two separable Hilbert spaces, and denote by
$\mathcal{L}_{2}^{0}=\mathcal{L}_{2}^{0}\left(V;H\right)$ the space of
Hilbert-Schmidt operators from $V$ into $H$. Let
$\left\\{W\left(t\right)\right\\}_{t\in\left[0,T\right]}$ be a $V-$valued,
$\mathbf{F}-$adapted cylindrical Wiener process on the filtered probability
space $\left(\Omega,\mathcal{F},\mathbf{F},P\right)$,
$\mathbf{F}=\left\\{\mathcal{F}_{t}\right\\}_{t\in\left[0,T\right]}$ with
standard conditions. Let $\mathbb{F}$ be the progressive $\sigma-$field (in
$\left[0,T\right]\times\Omega$) with respect to $\mathbf{F}$.

Following Ref. lu2021mathematical, we study the controlled stochastic PDE

$\displaystyle dx\left(t\right)=$
$\displaystyle\left(Ax\left(t\right)+a\left(t,x\left(t\right),\alpha\left(t\right)\right)\right)dt$
$\displaystyle+b\left(t,x\left(t\right),\alpha\left(t\right)\right)dW\left(t\right)\quad\textrm{in
}\left(0,T\right],$ (3) $\displaystyle x\left(0\right)=$ $\displaystyle
x_{0},$

where $A$ generates a $C_{0}-$semigroup,
$\left\\{S\left(t\right)\right\\}_{t\geq 0}$, on $H$,
$a\left(\cdot\right):\left[0,T\right]\times\Omega\times H\times U\rightarrow
H$ and $b\left(\cdot\right):\left[0,T\right]\times\Omega\times H\times
U\rightarrow\mathcal{L}_{2}^{0}$; $U$ is a separable Hilbert space. Moreover,
$x\left(\cdot\right)$ is the state variable (valued in $H$) and
$\alpha\left(\cdot\right)$ is the control variable (valued in $U$).

For $\xi\in\mathcal{O}\subset\mathbb{R}$, $d_{1},d_{2}>0$, and appropriate
spaces $H,V,U$ on $\mathcal{O}$ (e.g.,
$H^{1}\left(\mathcal{O};\mathbb{R}^{2}\right)$, plus boundary conditions), we
will set

$\displaystyle
x\left(t,\xi\right)=\left(\begin{array}[]{c}u\left(t,\xi\right)\\\
v\left(t,\xi\right)\end{array}\right),\textrm{
}A=\left(\begin{array}[]{c}d_{1}\frac{\partial^{2}}{\partial\xi^{2}}\\\
d_{2}\frac{\partial^{2}}{\partial\xi^{2}}\end{array}\right),\textrm{
}W\left(t\right)=\left(\begin{array}[]{c}W_{1}\left(t\right)\\\
W_{2}\left(t\right)\end{array}\right),\textrm{
}\alpha\left(t,\xi\right)=\left(\begin{array}[]{c}\alpha_{1}\left(t,\xi\right)\\\
\alpha_{2}\left(t,\xi\right)\end{array}\right),$ (12) $\displaystyle
e\left(t,x\left(t\right),\alpha\left(t\right)\right)\left(\xi\right)=\left(\begin{array}[]{c}e_{1}\left(t,u\left(t,\xi\right),v\left(t,\xi\right),\alpha_{1}\left(t,\xi\right),\alpha_{2}\left(t,\xi\right)\right)\\\
e_{2}\left(t,u\left(t,\xi\right),v\left(t,\xi\right),\alpha_{1}\left(t,\xi\right),\alpha_{2}\left(t,\xi\right)\right)\end{array}\right),$
(15) $\displaystyle
a\left(t,x\left(t\right),\alpha\left(t\right)\right)\left(\xi\right)=\left(\begin{array}[]{c}\tilde{f}\left(u\left(t,\xi\right),v\left(t,\xi\right)\right)\\\
\tilde{g}\left(u\left(t,\xi\right),v\left(t,\xi\right)\right)\end{array}\right)+e\left(t,x\left(t\right),\alpha\left(t\right)\right)\left(\xi\right),$
(18)
$\displaystyle\left(b\left(t,x\left(t\right),\alpha\left(t\right)\right)w\right)\left(\xi\right)=\sigma\left(\begin{array}[]{c}b_{1}\left(t,u\left(t,\xi\right),v\left(t,\xi\right),\alpha_{1}\left(t,\xi\right),\alpha_{2}\left(t,\xi\right)\right)w_{1}\\\
b_{2}\left(t,u\left(t,\xi\right),v\left(t,\xi\right),\alpha_{1}\left(t,\xi\right),\alpha_{2}\left(t,\xi\right)\right)w_{2}\end{array}\right),$
(21) $\displaystyle
w\left(\cdot\right)=\left(\begin{array}[]{c}w_{1}\left(\cdot\right)\\\
w_{2}\left(\cdot\right)\end{array}\right),$ (24)

where $\tilde{f}=-\partial F\left(u,v\right)/\partial u$, $\tilde{g}=-\partial
F\left(u,v\right)/\partial v$, $F$ is given by (1), $W_{1},W_{2}$ are
independent cylindrical Wiener processes, $\sigma>0$, and
$e_{1},e_{2},b_{1,}b_{2}$ have suitable regularity and growth (see below).

We start with the set of controls given by

$\mathcal{U}\left[0,T\right]=\left\\{\alpha:\left[0,T\right]\times\Omega\rightarrow
U:\alpha\left(\cdot\right)\textrm{ is }\mathbf{F}\textrm{-adapted}\right\\},$

which indicates that our controls are at least nonanticipative; we will see
bellow that, under our assumptions, $\mathcal{U}\left[0,T\right]$ is the set
of admissible controls. In Ref. lu2021mathematical, $U$ can be a separable
metric space; however, we only deal with the Hilbert case.

Define the cost functional $\mathcal{J}\left(\alpha\left(\cdot\right)\right)$
associated to (3) by

$\mathcal{J}\left(\alpha\left(\cdot\right)\right)=E\left[\int_{0}^{T}g\left(t,x\left(t\right),\alpha\left(t\right)\right)dt+h\left(x\left(T\right)\right)\right],\textrm{
}\forall\alpha\left(\cdot\right)\in\mathcal{U}\left[0,T\right].$ (25)

Our first objective is to study the following optimal control problem for the
controlled equation (3) with the cost functional (25):

Find $\bar{\alpha}\left(\cdot\right)\in\mathcal{U}\left[0,T\right]$ such that

$\mathcal{J}\left(\bar{\alpha}\left(\cdot\right)\right)=\underset{\alpha\left(\cdot\right)\in\mathcal{U}\left[0,T\right]}{\inf}\mathcal{J}\left(\alpha\left(\cdot\right)\right).$
(26)

We call any $\bar{\alpha}\left(\cdot\right)$ satisfying the last expression an
optimal control; the corresponding state $\bar{x}\left(\cdot\right)$ is an
optimal state and
$\left(\bar{x}\left(\cdot\right),\bar{\alpha}\left(\cdot\right)\right)$ is an
optimal pair.

As mentioned in Section 12.1,[22] to establish Pontryagin-type necessary
conditions for an optimal pair
$\left(\bar{x}\left(\cdot\right),\bar{\alpha}\left(\cdot\right)\right)$ of
Problem (26), we have to consider (along with appropriate assumptions) the
following $H-$valued backward stochastic differential equation:

$\displaystyle dy\left(t\right)=$
$\displaystyle-A^{*}y\left(t\right)dt-\left(a_{x}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right)^{*}y\left(t\right)+b_{x}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right)^{*}Y\left(t\right)\right.$
$\displaystyle\left.-g_{x}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right)\right)dt+Y\left(t\right)dW\left(t\right)\textrm{
in }\left[0,T\right),$ (27) $\displaystyle y\left(T\right)=$ $\displaystyle-
h_{x}\left(\bar{x}\left(T\right)\right).$

Next, we deal with the (more technical) optimal control problem with
endpoint/state constraints. Consider the controlled stochastic differential
equation (3) with $\alpha\in\mathcal{U}_{2}$,

$\mathcal{U}_{2}=\left\\{\alpha\left(\cdot\right):\left[0,T\right]\rightarrow
U:\alpha\left(\cdot\right)\in
L_{\mathbb{F}}^{2}\left(0,T;H_{1}\right)\right\\},$

where $U$ is a nonempty closed subset of the separable Hilbert space $H_{1}$.
Let $\mathcal{K}_{a}$ be a nonempty closed subset of $H$, and $h:\Omega\times
H\rightarrow\mathbb{R}$, $g^{j}:H\rightarrow\mathbb{R}$ ($j=0,\ldots,n$). We
associate to the control system (3) a Mayer cost functional,
$\mathcal{J}_{M}\left(\cdot\right)$, given by

$\mathcal{J}_{M}\left(\alpha\left(\cdot\right),x_{0}\right)=E\left[h\left(x\left(T\right)\right)\right],$
(28)

along with the state constraint

$E\left[g^{0}\left(x\left(t\right)\right)\right]\leq 0,\textrm{ for all
}t\in\left[0,T\right],$ (29)

and the initial-final states constraints

$x_{0}\in\mathcal{K}_{a},\quad
E\left[g^{j}\left(x\left(T\right)\right)\right]\leq 0,\quad j=1,\ldots,n.$
(30)

In this case, the set of admissible controls (with initial datum $x_{0}$) is
given by

$\mathcal{U}_{ad}^{x_{0}}=\left\\{\alpha\in\mathcal{U}_{2}:\textrm{ the
corresponding solution }x\left(t\right)\textrm{ of (\ref{eq:cont_sys})
satisfies (\ref{eq:state_const}) and (\ref{eq:ini_fin_const})}\right\\}.$

We study the following optimal control problem: find
$\left(\bar{x}_{0},\bar{\alpha}\left(\cdot\right)\right)\in\mathcal{K}_{a}\times\mathcal{U}_{ad}^{x_{0}}$
such that

$\mathcal{J}_{M}\left(\bar{x}_{0},\bar{\alpha}\left(\cdot\right)\right)=\underset{\left(x_{0},\alpha\left(\cdot\right)\right)\in\mathcal{K}_{a}\times\mathcal{U}_{ad}^{x_{0}}}{\inf}\mathcal{J}_{M}\left(x_{0},\alpha\left(\cdot\right)\right).$
(31)

We note that it is possible to study the more general Bolza problem using the
previous formulation (see Ref. frankowska2020first Section 1).

Let

$\varphi_{1}\left[t\right]=\varphi_{x}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right),\quad\varphi_{2}\left[t\right]=\varphi_{u}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right),$

where $\varphi$ can be either $a,b,f,g,$ or $h$(with appropriate regularity).
To analyze Problem (31), we will need the auxiliary linearized stochastic
control system

$\displaystyle dx_{1}\left(t\right)=$
$\displaystyle\left(Ax_{1}\left(t\right)+a_{1}\left[t\right]x_{1}\left(t\right)+a_{2}\left[t\right]\alpha_{1}\left(t\right)\right)dt$
$\displaystyle+\left(b_{1}\left[t\right]x_{1}\left(t\right)+b_{2}\left[t\right]\alpha_{1}\left(t\right)\right)dW\left(t\right)\textrm{
in }\left(0,T\right],$ $\displaystyle x_{1}\left(0\right)=$ $\displaystyle
x_{1},$ (32)

along with its first-order adjoint equation

$\displaystyle dy\left(t\right)=$
$\displaystyle-\left(A^{*}y\left(t\right)+a_{1}\left[t\right]^{*}y\left(t\right)+b_{1}\left[t\right]^{*}Y\left(t\right)\right)dt$
$\displaystyle+d\psi\left(t\right)+Y\left(t\right)dW\left(t\right)\textrm{ in
}\left(0,T\right],$ $\displaystyle y\left(T\right)=$ $\displaystyle y_{T},$
(33)

where $y_{T}\in L_{\mathcal{F}_{T}}^{2}\left(\Omega;H\right)$ and $\psi\in
L_{\mathbb{F}}^{2}\left(\Omega;BV_{0}\left(\left[0,T\right];H\right)\right)$.

From the biological perspective, it is natural to ask if there is a canonical
cost function associated with the epigenetic landscape. There are at least two
options which seem realistic: to consider an energy function or a minimal time
problem. The former requires to define an energy function corresponding to the
given system. In general, this function will depend on the detailed structure
of the underlying genetic regulatory network and the specific system. An
alternative, which relies on the assumption that the energetic cost spent in
traversing from a state (genetic expression profile), $A$, to another state,
$B$, is a convex function of the distance in the state space, e.g., the
distance squared. As for the second option, it is compatible with the
hypothesis that the energy spent is a monotone function of the time required
to reach state $B$ starting at $A$. Both problems are interesting from the
mathematical perspective; they arise naturally, especially the optimal time
problem, in the context of epigenetic therapy, where designing a control
policy for the disease is crucial.

##  3 Main Results

######  Theorem 3.1.

Let $\sigma_{1}=\sigma_{2}=\sigma>0$ and
$u\left(0,x\right),v\left(0,x\right)\in
C\left(\mathbb{T};\mathbb{R}\right)=C\left(\mathbb{T}\right)$. The stochastic
epigenetic system (2) has a unique generalized solution,
$u^{\sigma}\left(x,t\right)$, and the random process
$u^{\sigma}\left(t\right)=u^{\sigma}\left(t,\cdot\right)$ in the state space
$C\left(\mathbb{T}\right)$ is a Markov-Feller process. Furthermore,

  1. (i)

The process $u^{\sigma}\left(t\right)$ has a unique normed stationary measure,
$\nu^{\sigma}$, in $C\left(\mathbb{T};\mathbb{R}^{2}\right)$ such that, for
any Borel set $\Gamma\in C\left(\mathbb{T};\mathbb{R}^{2}\right)$ and any
$u_{0}\in C\left(\mathbb{T};\mathbb{R}^{2}\right)$,

$P_{u_{0}}\left\\{\underset{T\rightarrow\infty}{\lim}\frac{1}{T}\int_{0}^{T}\chi_{\Gamma}\left(u^{\sigma}\left(t\right)\right)dt=\nu^{\sigma}\left(\Gamma\right)\right\\}=1,$

where $\chi_{\Gamma}\left(u\right)$ is the indicator of the set $\Gamma\subset
C\left(\mathbb{T};\mathbb{R}^{2}\right)$.

  2. (ii)

Let
$\hat{\varphi}^{\left(1\right)},\ldots,\hat{\varphi}^{\left(4\right)}\in\mathbb{R}^{2}$
be the points where the epigenetic potential $F$ achieves its absolute
minimum, i.e., the points $\left(u_{1},v_{1}\right),$
$\left(u_{2},v_{2}\right),$ $\left(u_{3},v_{3}\right),$ and
$\left(u_{4},v_{4}\right)$. Let
$\triangle_{k}=\det\left(\mathbf{H}\left(F\right)\left(\hat{\varphi}^{k}\right)\right)=4a_{k}^{2}>0$
, $k=1,\ldots 4$, where $\mathbf{H}$ is the Hessian matrix. Then, the measure
$\nu^{\sigma}$ weakly converges as $\sigma\rightarrow 0$ to the measure
$\nu^{0}$ concentrated at the m points
$\hat{\varphi}^{\left(1\right)},\ldots,\hat{\varphi}^{\left(4\right)}\in\mathbb{R}^{2}$
and

$\nu^{0}\left(\hat{\varphi}^{\left(k\right)}\right)=\triangle_{k}^{-1}\left(\sum_{j=1}^{m}\triangle_{j}^{-1}\right)^{-1}.$

Since for arabidopsis thaliana $a_{k}=1/c_{k}$, $k=1,\ldots,4$ (see Section 3
Ref. cortes2022variational), our previous result suggests
$a_{k}=1/\left(2\sqrt{c_{k}}\right)$.

  3. (iii)

$\varphi_{0}\left(x\right)\equiv\hat{\varphi}^{\left(k\right)}$,
$k\in\left\\{1,\ldots,4\right\\}$ is an asymptotically stable equilibrium
point of the deterministic epigenetic system, i.e., (2) with
$\sigma_{1}=\sigma_{2}=0$. In addition, the phenotype (or cellular type)
regions associated to each $\hat{\varphi}^{\left(k\right)}$ are invariant.
Assume that $D$ is a regular region (see Subsection 4.1 below) in
$C\left(\mathbb{T};\mathbb{R}^{2}\right)$ which contains $\varphi_{0}$ such
that all the trajectories (of the deterministic system) starting from $g\in
D\cup\partial D$ tend to $\varphi_{0}$ without leaving $D$. Let

$U\left(\varphi\right)=\int_{\mathbb{T}}\left[\frac{1}{2}\sum_{k=1}^{2}d_{k}\left(\frac{d\varphi_{k}}{dx}\right)^{2}+F\left(\varphi\left(x\right)\right)\right]dx,\\\
\varphi\left(x\right)=\left(\varphi_{1}\left(x\right),\varphi_{2}\left(x\right)\right)\in
C\left(\mathbb{T};\mathbb{R}^{2}\right).$

Let $\tau^{\sigma}=\tau_{D}^{\sigma}=\inf\left\\{t:\textrm{
}u^{\sigma}\left(t\right)\notin D\right\\}$ be the first exist time of
$u^{\sigma}\left(t\right)$ from $D$. Then

$\underset{\sigma\rightarrow 0^{+}}{\lim}\textrm{ }\sigma^{2}\ln
E_{g}\tau^{\sigma}=2\underset{\varphi\in\partial
D}{\min}\left(U\left(\varphi\right)-U\left(\varphi_{0}\right)\right),\textrm{
}g\in D.$

Now, we have the following set of assumptions:

###### Assumption 3.1

Let $e\left(\cdot,\cdot,\cdot\right):\left[0,T\right]\times\Omega\times
H\times U\rightarrow H$ and
$b\left(\cdot,\cdot,\cdot\right):\left[0,T\right]\times\Omega\times H\times
U\rightarrow\mathcal{L}_{2}^{0}$ satisfy:

  1. (i)

For any $\left(x,\alpha\right)\in H\times U$, the functions
$e\left(\cdot,x,\alpha\right):\left[0,T\right]\times\Omega\rightarrow H$ and
$b\left(\cdot,x,\alpha\right)\rightarrow\mathcal{L}_{2}^{0}$ are
$\mathbb{F}-$measurable,

  2. (ii)

For any $x\in H$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$, the functions
$e\left(t,x,\cdot\right):U\rightarrow H$ and
$b\left(t,x,\cdot\right):U\rightarrow\mathcal{L}_{2}^{0}$ are continuous,

  3. (iii)

For any $\left(x_{1},x_{2},\alpha\right)\in H\times H\times U$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$,

$\left\\{\begin{array}[]{l}\left|e\left(t,x_{1},\alpha\right)-e\left(t,x_{2},\alpha\right)\right|_{H}+\left|b\left(t,x_{1},\alpha\right)-b\left(t,x_{2},\alpha\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C\left|x_{1}-x_{2}\right|_{H},\\\
\left|e\left(t,0,\alpha\right)\right|_{H}+\left|b\left(t,0,\alpha\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C.\end{array}\right.$

###### Assumption 3.2

Let $g\left(\cdot,\cdot,\cdot\right):\left[0,T\right]\times\Omega\times
H\times U\rightarrow\mathbb{R}$ and $h\left(\cdot\right):\Omega\times
H\rightarrow\mathbb{R}$ be two functions satisfying:

  1. (i)

For any $\left(x,\alpha\right)\in H\times U$,
$g\left(\cdot,x,\alpha\right):\left[0,T\right]\times\Omega\rightarrow\mathbb{R}$
is $\mathbb{F}-$measurable and $h\left(x\right):\Omega\rightarrow\mathbb{R}$
is $\mathcal{F}_{T}-$measurable,

  2. (ii)

For any $x\in H$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$,
$g\left(t,x,\cdot\right):U\rightarrow\mathbb{R}$ is continuous,

  3. (iii)

For any $\left(x_{1},x_{2},\alpha\right)\in H\times H\times U$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$,

$\left\\{\begin{array}[]{l}\left|g\left(t,x_{1},\alpha\right)-g\left(t,x_{2},\alpha\right)\right|+\left|h\left(x_{1}\right)-h\left(x_{2}\right)\right|\leq
C\left|x_{1}-x_{2}\right|_{H},\\\
\left|g\left(t,0,\alpha\right)\right|+\left|h\left(0\right)\right|\leq
C.\end{array}\right.$

###### Assumption 3.3

The control region $U$ is a convex subset of a separable Hilbert space,
$\tilde{H}$, and the metric of $U$ is induced by the norm of $\tilde{H}$, that
is,
$d\left(\alpha_{1},\alpha_{2}\right)=\left|\alpha_{1}-\alpha_{2}\right|_{\tilde{H}}$.

###### Assumption 3.4

For a.e. $\left(t,\omega\right)\in\left(0,T\right)\times\Omega$, the functions
$e\left(t,\cdot,\cdot\right):H\times U\rightarrow H$,
$b\left(t,\cdot,\cdot\right):H\times U\rightarrow\mathcal{L}_{2}^{0}$,
$g\left(t,\cdot,\cdot\right):H\times U\rightarrow\mathbb{R}$, and
$h\left(\cdot\right):H\rightarrow\mathbb{R}$ are $C^{1}$. Furthermore, for any
$\left(x,\alpha\right)\in H\times U$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$, we have

$\left\\{\begin{array}[]{l}\left|e_{x}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(H\right)}+\left|b_{x}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(H;\mathcal{L}_{2}^{0}\right)}+\left|g_{x}\left(t,x,\alpha\right)\right|_{H}+\left|h_{x}\left(x\right)\right|_{H}\leq
C,\\\
\left|e_{\alpha}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(\tilde{H};H\right)}+\left|b_{\alpha}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(\tilde{H};\mathcal{L}_{2}^{0}\right)}+\left|g_{\alpha}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(\tilde{H}\right)}\leq
C\end{array}\right.$

######  Theorem 3.2.

Consider the controlled stochastic epigenetic system (3)-(24). Let Assumptions
3.1-3.4 hold. Then,

  1. (i)

For any $x_{0}\in L_{\mathcal{F}_{0}}^{p_{0}}\left(\Omega;H\right)$,
$p_{0}\geq 2$, and $\alpha\left(\cdot\right)\in\mathcal{U}\left[0,T\right]$,
system (3) has a unique mild solution, $x\left(\cdot\right)\equiv
x\left(\cdot;x_{0},\alpha\right)\in
C_{\mathbb{F}}\left(\left[0,T\right];L^{p_{0}}\left(\Omega;H\right)\right)$,
such that

$\left|x\left(\cdot\right)\right|_{C_{\mathbb{F}}\left(\left[0,T\right];L^{p_{0}}\left(\Omega;H\right)\right)}\leq
C\left(1+\left|x_{0}\right|_{L_{\mathcal{F}_{0}}^{p_{0}}\left(\Omega;H\right)}\right).$

Moreover, equation (27) is well-posed in the sense of transposition solution
(see Definition 4.13 Ref. lu2021mathematical).

  2. (ii)

Let $\left(\bar{x}\left(\cdot\right),\bar{\alpha}\left(\cdot\right)\right)$ be
an optimal pair for Problem (26) with $x_{0}\in
L_{\mathcal{F}_{0}}^{2}\left(\Omega;H\right)$. Then,

$\textrm{Re}\left\langle
a_{u}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right)^{*}y\left(t\right)+b_{u}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right)^{*}Y\left(t\right)-g_{u}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right),\right.\\\
\left.\alpha-\bar{\alpha}\left(t\right)\right\rangle_{\tilde{H}}\leq 0,$

a.e. $\left(t,\omega\right)\in\left[0,T\right]\times\Omega$, $\forall\alpha\in
U$, where $\left(y\left(\cdot\right),Y\left(\cdot\right)\right)$ is the
transposition solution of (27).

Next, we need some notation and ideas from set-valued analysis; see Subsection
4.2 for details. In addition, we have the following assumptions:

###### Assumption 3.5

$e\left(\cdot,\cdot,\cdot,\cdot\right):\left[0,T\right]\times H\times
H_{1}\times\Omega\rightarrow H$ and
$b\left(\cdot,\cdot,\cdot,\cdot\right):\left[0,T\right]\times H\times
H_{1}\times\Omega\rightarrow\mathcal{L}_{2}^{0}$ are two maps such that

  1. (i)

For any $\left(x,\alpha\right)\in H\times H_{1}$,
$e\left(\cdot,x,\alpha,\cdot\right):\left[0,T\right]\times\Omega\rightarrow H$
and
$b\left(\cdot,x,\alpha,\cdot\right):\left[0,T\right]\times\Omega\rightarrow\mathcal{L}_{2}^{0}$
are $\mathcal{B}\left(\left[0,T\right]\right)\times\mathcal{F}$ measurable and
$\mathbb{F}-$adapted,

  2. (ii)

For any $\left(t,x,\omega\right)\in\left[0,T\right]\times H\times\Omega$,
$e\left(t,x,\cdot,\omega\right):H_{1}\rightarrow H$ and
$b\left(t,x,\cdot,\omega\right):H_{1}\rightarrow\mathcal{L}_{2}^{0}$ are
continuous and

$\begin{cases}\left|e\left(t,x_{1},\alpha,\omega\right)-e\left(t,x_{2},\alpha,\omega\right)\right|_{H}+\left|b\left(t,x_{1},\alpha,\omega\right)-b\left(t,x_{2},\alpha,\omega\right)\right|_{\mathcal{L}_{2}^{0}}\leq\\\
\qquad
C\left|x_{1}-x_{2}\right|_{H}\quad\forall\left(t,x_{1},x_{2},\alpha,\omega\right)\in\left[0,T\right]\times
H\times H\times H_{1}\times\Omega\\\
\left|e\left(t,0,\alpha,\omega\right)\right|_{H}+\left|b\left(t,0,\alpha,\omega\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C,\quad\forall\left(t,\alpha,\omega\right)\in\left[0,T\right]\times
H_{1}\times\Omega.\end{cases}$

###### Assumption 3.6

For a.e. $\left(t,\omega\right)\in\left[0,T\right]\times\Omega$, the functions
$e\left(t,\cdot,\cdot,\omega\right):H\times H_{1}\rightarrow H$ and
$b\left(t,\cdot,\cdot,\omega\right):H\times
H_{1}\rightarrow\mathcal{L}_{2}^{0}$ are differentiable, and
$\left(e_{x}\left(t,x,\alpha,\omega\right),e_{\alpha}\left(t,x,\alpha,\omega\right)\right)$
and
$\left(b_{x}\left(t,x,\alpha,\omega\right),b_{\alpha}\left(t,x,\alpha,\omega\right)\right)$
are uniformly continuous with respect to $x\in H$ and $\alpha\in U$ (Fréchet
differentiability). There exists a nonnegative $\eta\in
L_{\mathbb{F}}^{2}\left(0,T;\mathbb{R}\right)$ such that for a.e.
$\left(t,\omega\right)\in\left[0,T\right]\times\Omega$ and for all $x\in H$
and $\alpha\in H_{1}$,

$\begin{cases}\left|e\left(t,0,\alpha,\omega\right)\right|_{H}+\left|b\left(t,0,\alpha,\omega\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C\left(\eta\left(t,\omega\right)+\left|\alpha\right|_{H_{1}}\right),\\\
\left|e_{x}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H\right)}+\left|e_{\alpha}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H_{1};H\right)}+\left|b_{x}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H;\mathcal{L}_{2}^{0}\right)}\\\
\hfill+\left|b_{\alpha}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H_{1};\mathcal{L}_{2}^{0}\right)}\leq
C.\end{cases}$

###### Assumption 3.7

The functional $h\left(\cdot,\omega\right):H\rightarrow\mathbb{R}$ is
differentiable $\mathbb{P}-$a.s., and there exists an $\eta\in
L_{\mathcal{F}_{T}}^{2}\left(\Omega\right)$ such that for any $x,\tilde{x}\in
H$,

$\begin{cases}h\left(x,\omega\right)\leq
C\left(\eta\left(\omega\right)^{2}+\left|x\right|_{H}^{2}\right),\quad\left|h_{x}\left(0,\omega\right)\right|_{H}\leq
C\eta\left(\omega\right),\quad\mathbb{P}-a.s.,\\\
\left|h_{x}\left(x,\omega\right)-h_{x}\left(\tilde{x},\omega\right)\right|_{H}\leq
C\left|x-\tilde{x}\right|_{H},\quad\mathbb{P}-a.s.\end{cases}$

###### Assumption 3.8

For $j=0,\ldots,n$, the functional $g^{j}:H\rightarrow\mathbb{R}$ is
differentiable and for any $x,\tilde{x}\in H$,

$\left|g^{j}\left(x\right)\right|\leq
C\left(1+\left|x\right|_{H}^{2}\right),\quad\left|g_{x}^{j}\left(x\right)-g_{x}^{j}\left(\tilde{x}\right)\right|_{H}\leq
C\left|x-\tilde{x}\right|_{H}.$

Define the Hamiltonian

$\mathbb{H}^{epig}\left(t,x,\alpha,p,q\right)=\left\langle
p,a\left(t,x,\alpha\right)\right\rangle_{H}+\left\langle
q,b\left(t,x,\alpha\right)\right\rangle_{\mathcal{L}_{2}^{0}},$

where $\left(t,x,\alpha,p,q\right)\in\left[0,T\right]\times H\times
H_{1}\times H\times\mathcal{L}_{2}^{0}$ , with $a\left(\cdot\right)$,
$b\left(\cdot\right)$ given by (24).

######  Theorem 3.3.

Consider the controlled stochastic epigenetic system (3)-(24). Let Assumptions
3.5-3.8 hold. Then

  1. (i)

For any $x_{0}\in H$ and $\alpha\left(\cdot\right)\in\mathcal{U}_{2}$, system
(3) has a unique mild solution, $x\left(\cdot\right)\equiv
x\left(\cdot;x_{0},\alpha\right)\in
L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$, such
that

$\left|x\left(\cdot\right)\right|_{L^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)}\leq
C\left(1+\left|x_{0}\right|_{H}\right).$

Moreover, for any $\alpha_{1}\in\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)$
and $x_{1}\in T_{\mathcal{K}_{a}}^{b}\left(\bar{x}_{0}\right)$, (32) has a
unique solution, $x_{1}\left(\cdot\right)\in
L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$, and
for $\psi\in
L_{\mathbb{F}}^{2}\left(\Omega;BV_{0}\left(\left[0,T\right];H\right)\right)$,
(33) has a unique transposition solution $\left(y,Y\right)\in
D_{\mathbb{F}}\left(\left[0,T\right];L^{2}\left(\Omega;H\right)\right)\times
L_{\mathbb{F}}^{2}\left(0,T;\mathcal{L}_{2}^{0}\right)$.

  2. (ii)

Let
$\left(\bar{x}\left(\cdot\right),\bar{\alpha}\left(\cdot\right),\bar{x}_{0}\right)$
be an optimal triple of Problem (31). If
$E\left|g_{x}^{0}\left(\bar{x}\left(t\right)\right)\right|_{H}\neq 0$ for any
$t\in\mathcal{I}_{0}\left(\bar{x}\right)$, then there exists
$\lambda_{0}\in\left\\{0,1\right\\}$, $\lambda_{j}\geq 0$ for
$j\in\mathcal{I}\left(\bar{x}\right)$ and
$\psi\in\left(\mathcal{Q}^{\left(1\right)}\right)^{-}$ with
$\psi\left(0\right)=0$ verifying

$\lambda_{0}+\sum_{j\in\mathcal{I}\left(\bar{x}\right)}\lambda_{j}+\left|\psi\right|_{L_{\mathbb{F}}^{2}\left(\Omega;BV\left(0,T;H\right)\right)}\neq
0,$

such that the corresponding transposition solution
$\left(y\left(\cdot\right),Y\left(\cdot\right)\right)$ of the first order
adjoint equation (33) with
$y\left(T\right)=-\lambda_{0}h_{x}\left(\bar{x}\left(T\right)\right)-\sum_{j\in\mathcal{I}\left(\bar{x}\right)}\lambda_{j}g_{x}^{j}\left(\bar{x}\left(T\right)\right)$
satisfies the variational inequality

$E\left\langle
y\left(0\right),\nu\right\rangle_{H}+E\int_{0}^{T}\left\langle\mathbb{H}_{\alpha}^{epig}\left[t\right],v\left(t\right)\right\rangle_{H_{1}}dt\leq
0,\textrm{
}\forall\nu\in\mathcal{T}_{\mathcal{K}_{a}}\left(\bar{x}_{0}\right),\textrm{
}\forall v\left(\cdot\right)\in\mathcal{T}_{\Phi}\left(\bar{\alpha}\right),$

where
$\mathbb{H}_{\alpha}^{epig}\left[t\right]=\mathbb{H}_{\alpha}^{epig}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right),y\left(t\right),Y\left(t\right)\right)$.
Furthermore, if
$\mathcal{G}^{\left(1\right)}\cap\mathcal{Q}^{\left(1\right)}\cap\mathcal{E}^{\left(1\right)}\neq\emptyset$,
the above holds with $\lambda_{0}=1$.

##  4 Preliminaries

For a thorough review of ergodic theory for stochastic PDE, see, e.g., Ref.
da1996ergodicity,hairer2008ergodic. For the the theory of Large deviations,
see Ref. freidlin1998random,varadhan2016large. See Ref.
aastrom2012introduction,fleming2012deterministic,lu2021mathematical for a
review of stochastic optimal control.

###  4.1 Random perturbations and large deviations

In this section, we focus on the results presented in Ref. freidlin1988random
for the system

$\displaystyle\frac{\partial u_{k}^{\varepsilon}\left(t,x\right)}{\partial
t}=$ $\displaystyle D_{k}\frac{\partial^{2}u_{k}^{\varepsilon}}{\partial
x^{2}}+f_{k}\left(x,u_{1}^{\varepsilon},\ldots,u_{n}^{\varepsilon}\right)+\varepsilon\hat{\varsigma}_{k}\left(t,x\right),$
(34) $\displaystyle u_{k}^{\varepsilon}\left(0,x\right)=$ $\displaystyle
g_{k}\left(x\right),$

$t>0$, $x\in\mathbb{T}$, $k=1,\ldots,n$. The perturbations
$\hat{\varsigma}_{k}\left(t,x\right)$, $k=1,\ldots,n$, are Gaussian fields
which have independent values for different $t$. In what follows, we consider
a white noise perturbation:
$\hat{\varsigma}_{k}\left(t,x\right)=\partial^{2}\varsigma_{k}\left(t,x\right)/\partial
t\partial x$, where $\varsigma_{k}\left(t,x\right)$ are independent Brownian
sheets for different $k$. We assume that the functions
$f_{k}\left(x,u\right)$, $x\in\mathbb{T}$, $u\in\mathbb{R}^{n}$, are Lipschitz
continuous. Moreover, $g_{k}\left(x\right)\in C\left(\mathbb{T}\right)$ and
$D_{k}>0$, $k=1,\ldots,n$.

A generalized solution of (34) is a measurable function
$u^{\varepsilon}\left(t,x\right)=\left(u_{1}^{\varepsilon}\left(t,x\right),\ldots,u_{n}^{\varepsilon}\left(t,x\right)\right)$
such that, with probability 1,

$\int_{\mathbb{T}}u_{k}^{\varepsilon}\left(t,x\right)\varphi\left(x\right)dx-\int_{\mathbb{T}}g_{k}\left(x\right)\varphi\left(x\right)dx\\\
=\int_{0}^{t}\int_{\mathbb{T}}\left[u_{k}^{\varepsilon}\left(s,x\right)D_{k}\varphi^{\prime\prime}\left(x\right)-f_{k}\left(x,u_{k}^{\varepsilon}\left(s,x\right)\right)\varphi\left(x\right)\right]dsdx+\varepsilon\int_{\mathbb{T}}\varphi^{\prime}\left(x\right)\varsigma_{k}\left(t,x\right)dx,$

for any $\varphi\in C^{\infty}\left(\mathbb{T}\right)$, $k=1,\ldots,n$, and
$t>0$.

Under our assumptions, Theorem 1 in Ref. freidlin1988random ensures that (34)
has a unique generalized solution. Furthermore, the random process
$u^{\varepsilon}\left(t\right)=u^{\varepsilon}\left(t,\cdot\right)$ in the
state space $C\left(\mathbb{T}\right)$ is a Markov-Feller process.

We need the following auxiliary (linear) equation:

$\frac{\partial v^{\varepsilon}\left(t,x\right)}{\partial
t}=D\frac{\partial^{2}v^{\varepsilon}}{\partial x^{2}}-\alpha
v^{\varepsilon}+\varepsilon\frac{\partial^{2}\varsigma}{\partial t\partial
x},\textrm{ }v^{\varepsilon}\left(0,x\right)=g\left(x\right)\in
C\left(\mathbb{T};\mathbb{R}\right),$ (35)

where $\varsigma$ is a Brownian sheet, $\varepsilon,D,\alpha>0$. (35) has a
unique generalized solution and
$v^{\varepsilon}\left(t\right)=v^{\varepsilon}\left(t,\cdot\right)$ is a
Markov process in the state space $C\left(\mathbb{T};\mathbb{R}\right)$.
Moreover, by Lemma 1 in Ref. freidlin1988random,
$v^{\varepsilon}\left(t\right)$ has a unique stationary distribution,
$\mu^{\varepsilon}=\mu_{\alpha}^{\varepsilon}$, which is mean zero Gaussian
with correlation function

$\displaystyle B\left(x,y\right)=$
$\displaystyle\frac{\varepsilon^{2}}{2\pi}\sum_{k=1}^{\infty}\frac{1}{\lambda_{k}}\cos
k\left(x-y\right)\textrm{ }x,y\in\mathbb{T},$ $\displaystyle\lambda_{k}=$
$\displaystyle Dk^{2}+\alpha.$

Let
$f\left(x,u\right)=\left(f_{1}\left(x,u\right),\ldots,f_{n}\left(x,u\right)\right)$,
where $u\in\mathbb{R}^{n}$ and $x\in\mathbb{T}$ is a parameter. We call $f$ a
potential field provided there exists a function, $F\left(x,u\right)$,
continuously differentiable in the variables $u\in\mathbb{R}^{n}$ and such
that $f_{k}\left(x,u\right)=-\partial F\left(x,u\right)/\partial u_{k}$,
$x\in\mathbb{T}$, $u\in\mathbb{R}^{n}$, $k=1,\ldots,n$. Our main interest is
system (34) with the potential field $f\left(x,u\right)=-\nabla
F\left(x,u\right)$.

Consider
$B\left(x,u\right)=\left(B_{1}\left(x,u\right),\ldots,B_{n}\left(x,u\right)\right),$
$x\in\mathbb{T}$, $u\in C\left(\mathbb{T}\right)$, where

$B_{k}\left(x,u\right)=D_{k}\frac{d^{2}u}{dx^{2}}+f_{k}\left(x,u\right).$

With $\varepsilon=0$, equation (34) defines the semiflow $u\left(t\right)$ in
$C\left(\mathbb{T}\right)$ given by
$du\left(t\right)/dt=B\left(u\left(t\right)\right)$. Let

$U\left(\varphi\right)=\int_{\mathbb{T}}\left[\frac{1}{2}\sum_{k=1}^{n}D_{k}\left(\frac{d\varphi_{k}}{dx}\right)^{2}+F\left(x,\varphi\left(x\right)\right)\right]dx,\\\
\varphi\left(x\right)=\left(\varphi_{1}\left(x\right),\ldots,\varphi_{n}\left(x\right)\right)\in
C\left(\mathbb{T};\mathbb{R}^{n}\right).$ (36)

One can verify that the variational derivative of $U\left(\varphi\right)$
(with negative sign) is

$-\frac{\delta
U\left(\varphi\right)}{\delta\varphi_{k}}=D_{k}\frac{d\varphi_{k}}{dx}+f_{k}\left(x,\varphi\left(x\right)\right)=B_{k}\left(x,\varphi\right),$

which leads us to consider $U\left(\varphi\right)$ as the potential of the
field $B\left(\varphi\right)$ (see Section 3 of Ref. freidlin1988random for
details). Our goal is to define the stationary measure of
$u^{\varepsilon}\left(t\right)$ in $C\left(\mathbb{T};\mathbb{R}^{n}\right)$
and the action functional for the family of fields
$u^{\varepsilon}\left(t,x\right)$ in
$C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right)$.

Consider the Gaussian measure
$\mu^{\varepsilon}=\mu_{\alpha_{1},\ldots,\alpha_{n}}^{\varepsilon}$,
$\alpha=\left(\alpha_{1},\ldots,\alpha_{k}\right)$, $\alpha_{k}>0$,
$k=1,\ldots,n$, defined as the direct product of the measures
$\mu^{\varepsilon}=\mu_{\alpha_{k}}^{\varepsilon}$ in
$C\left(\mathbb{T};\mathbb{R}\right)$, where $\mu_{\alpha_{k}}^{\varepsilon}$
is the stationary distribution of the process $v^{\varepsilon}\left(t\right)$
with $\alpha=\alpha_{k}$ in (35). Denote by $E^{\alpha}$ the expectation with
respect to the measure $\mu_{\alpha_{1},\ldots,\alpha_{n}}^{\varepsilon}$:

$E^{\alpha}G\left(\varphi\right)=\int_{C\left(\mathbb{T};\mathbb{R}^{n}\right)}G\left(\varphi\right)\mu_{\alpha_{1},\ldots,\alpha_{n}}^{\varepsilon}\left(d\varphi\right).$

Assume that
$f\left(x,u\right)=\left(f_{1}\left(x,u\right),\ldots,f_{n}\left(x,u\right)\right)$
has a potential, $F\left(x,u\right)$, $x\in\mathbb{T}$, $u\in\mathbb{R}^{n}$,
and let for some $\alpha=\left(\alpha_{1},\ldots,\alpha_{k}\right)$,
$\alpha_{k}>0$, $k=1,\ldots,n$,

$A_{\varepsilon}=E^{\alpha}\exp\left\\{-\frac{2}{\varepsilon^{2}}\int_{\mathbb{T}}\left[F\left(x,\varphi\left(x\right)\right)-\frac{1}{2}\sum_{k=1}^{n}\alpha_{k}\varphi_{k}^{2}\left(x\right)\right]dx\right\\}<\infty.$

Let $\nu^{\varepsilon}$ be the measure on
$C\left(\mathbb{T};\mathbb{R}^{n}\right)$ such that

$\frac{d\nu^{\varepsilon}}{d\mu_{\alpha_{1},\ldots,\alpha_{n}}^{\varepsilon}}\left(\varphi\right)=A_{\varepsilon}^{-1}\exp\left\\{-\frac{2}{\varepsilon^{2}}\int_{\mathbb{T}}\left[F\left(x,\varphi\left(x\right)\right)-\frac{1}{2}\sum_{k=1}^{n}\alpha_{k}\varphi_{k}^{2}\left(x\right)\right]dx\right\\}.$
(37)

Then, by Theorem 2 in Ref. freidlin1988random, $\nu^{\varepsilon}$ is the
unique normed stationary measure of the process
$u^{\varepsilon}\left(t\right)$ in $C\left(\mathbb{T};\mathbb{R}^{n}\right)$
defined by (34). Furthermore, for any Borel set $\Gamma\in
C\left(\mathbb{T};\mathbb{R}^{n}\right)$ and any $u_{0}\in
C\left(\mathbb{T};\mathbb{R}^{n}\right)$,

$P_{u_{0}}\left\\{\underset{T\rightarrow\infty}{\lim}\frac{1}{T}\int_{0}^{T}\chi_{\Gamma}\left(u^{\varepsilon}\left(t\right)\right)dt=\nu^{\varepsilon}\left(\Gamma\right)\right\\}=1,$

where $\chi_{\Gamma}\left(u\right)$ is the indicator of the set $\Gamma\subset
C\left(\mathbb{T};\mathbb{R}^{n}\right)$.

Now consider a potential field not depending on $x$, i.e.,
$f\left(x,u\right)=f\left(u\right)$. In this case, the potential

$U\left(\varphi\right)=\int_{\mathbb{T}}\left[\frac{1}{2}\sum_{k=1}^{n}D_{k}\left(\frac{d\varphi_{k}}{dx}\right)^{2}+F\left(\varphi\left(x\right)\right)\right]dx$

takes its absolute minimum at the functions
$\hat{\varphi}=\left(\hat{\varphi}_{1},\ldots,\hat{\varphi}_{n}\right)$ with
the components independent of $x$. The vector $\hat{\varphi}$ is defined as
the point at which the absolute minimum of the function $F\left(z\right)$,
$z\in\mathbb{R}^{n}$, is attained. If the absolute minimum of
$F\left(z\right)$ is achieved at several points
$\hat{\varphi}^{\left(1\right)},\ldots,\hat{\varphi}^{\left(m\right)}\in\mathbb{R}^{n}$
, then the limit measure is distributed over these points. More precisely,
assume that $f\left(u\right)=-\nabla F\left(u\right)$. Furthermore, assume
that $F\left(u\right)$ is three times continuously differentiable, satisfies
the inequality $F\left(u\right)>\alpha\left|u\right|+\beta$ for some
$\alpha>0$, $\beta\in\left(-\infty,\infty\right)$ and attains its absolute
minimum at $m$ points
$\hat{\varphi}^{\left(1\right)},\ldots,\hat{\varphi}^{\left(m\right)}\in\mathbb{R}^{n}$.
Moreover, let the critical points be nondegenerate, that is,
$\triangle_{k}=\det\left(\partial^{2}F\left(\hat{\varphi}^{k}\right)/\partial
u_{i}\partial u_{j}\right)\neq 0$ for $k=1,\ldots m$. Then, by Theorem 5 in
Ref. freidlin1988random, the measure $\nu^{\varepsilon}$ weakly converges as
$\varepsilon\rightarrow 0$ to the measure $\nu^{0}$ concentrated at $m$ points
$\hat{\varphi}^{\left(1\right)},\ldots,\hat{\varphi}^{\left(m\right)}\in\mathbb{R}^{n}$
and

$\nu^{0}\left(\hat{\varphi}^{\left(k\right)}\right)=\triangle_{k}^{-1}\left(\sum_{j=1}^{m}\triangle_{j}^{-1}\right)^{-1}.$

To follow the exposition in Ref. freidlin1988random, for the rest of this
subsection we assume that the field $f\left(x,u\right)$ is not necessarily
potential. However, we only study the first exit time of
$u^{\varepsilon}\left(t\right)$ from a $regular$ region.

Let $H^{1,2}$ be the Sobolev space of functions of $t\in\left[0,T\right]$ and
$x\in\mathbb{T}$ with values in $\mathbb{R}^{n}$ which have square integrable
generalized derivatives of first order in $t$ and second order in $x$. We have
$H^{1,2}\subset
C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right)$, with
continuous embedding.

By Theorem 6 in Ref. freidlin1988random, the action functional for the family
of fields $u^{\varepsilon}\left(t,x\right)$, $0\leq t\leq T$,
$x\in\mathbb{T}$, in
$C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right)$ as
$\varepsilon\rightarrow 0^{+}$ has the form
$\varepsilon^{-2}S^{u}\left(\varphi\right)$ with

$S^{u}\left(\varphi\right)=\left\\{\begin{array}[]{l}\frac{1}{2}\int_{\mathbb{T}}\int_{0}^{T}\sum_{k=1}^{n}\left|\frac{\partial\varphi_{k}}{\partial
t}-D_{k}\frac{\partial^{2}\varphi_{k}}{\partial
x^{2}}-f_{k}\left(x,\varphi\left(t,x\right)\right)\right|^{2}dtdx,\textrm{
}\varphi\in H^{1,2}\\\ +\infty\textrm{ if }\varphi\in
C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right)\backslash
H^{1,2}.\end{array}\right.$

Moreover, the functional $S^{u}\left(\varphi\right)$ is lower semicontinuous
on $C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right)$ and for
every $s<\infty$, $g\in C\left(\mathbb{T};\mathbb{R}^{n}\right)$, the set
$\left\\{\varphi\in
C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right):\varphi\left(0,x\right)=g\left(x\right),S^{u}\left(\varphi\right)\leq
s\right\\}$ is compact in
$C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right)$.

We need the functional

$V\left(g,h\right)=\inf\left\\{S^{u}\left(\varphi\right):\varphi\in
C\left(\left[0,T\right]\times\mathbb{T};\mathbb{R}^{n}\right),\varphi\left(0,x\right)=g\left(x\right),\right.\\\
\left.\varphi\left(T,x\right)=h\left(x\right),\textrm{ }T\geq 0\right\\},\quad
g,h\in C\left(\mathbb{T};\mathbb{R}^{n}\right).$

Assume that $\varphi_{0}\in C\left(\mathbb{T};\mathbb{R}^{n}\right)$ is an
asymptotically stable equilibrium point of (34) with $\varepsilon=0$ and let
$D$ be a bounded open region in $C\left(\mathbb{T};\mathbb{R}^{n}\right)$
containing $\varphi_{0}$. The region $D\subset
C\left(\mathbb{T};\mathbb{R}^{n}\right)$ is called regular if for every
$\varphi\in\partial D$ there is a twice continuously differentiable function,
$h=h_{\varphi}\in C\left(\mathbb{T};\mathbb{R}^{n}\right)$, such that
$\varphi+th$ is an interior point of the complement of $D\cup\partial D$ for
all $t\geq 0$ small enough.

Let
$\tau^{\varepsilon}=\tau_{D}^{\varepsilon}=\inf\left\\{t:u^{\varepsilon}\left(t\right)\notin
D\right\\}$ be the first exit time of $u\left(t\right)$ from $D$ and
$V_{0}=\inf\left\\{V\left(\varphi_{0},\varphi\right):\varphi\in\partial
D\right\\}$. Assume, that $D\in C\left(\mathbb{T};\mathbb{R}^{n}\right)$ is
regular and that $\varphi_{0}\in D$ is an asymptotically stable point of (34)
with $\varepsilon=0$. Furthermore, assume that every trajectory of (34) with
$\varepsilon=0$ starting at a point $g\in D\cup\partial D$ does not leave $D$
for $t>0$ and tends to $\varphi_{0}$ as $t\rightarrow\infty$. Then, by Theorem
8 in Ref. freidlin1988random, for any $g\in D$

$\underset{\varepsilon\rightarrow 0^{+}}{\lim}\textrm{ }\varepsilon^{2}\ln
E_{g}\tau^{\varepsilon}=V_{0}.$

Moreover, if there is a unique $\varphi^{*}\in\partial D$ for which
$V\left(\varphi_{0},\varphi^{*}\right)=V_{0}$ , then the process
$u^{\varepsilon}\left(t\right)$ exists $D$ for the first time near
$\varphi^{*}$, that is, for any $\delta>0$ and $\varphi\in D$

$\underset{\varepsilon\rightarrow
0^{+}}{\lim}P_{\varphi}\left\\{\underset{x\in\mathbb{T}}{\sup}\left|u_{\tau^{\varepsilon}}^{\varepsilon}\left(x\right)-\varphi^{*}\left(x\right)\right|>\delta\right\\}=0.$

Let $H^{1}$ be the Sobolev space of functions on $\mathbb{T}$ with values in
$\mathbb{R}^{n}$ with square integrable first-order generalized derivatives.
Now consider the functional $\mathscr{U}\left(\varphi\right)$ on
$C\left(\mathbb{T};\mathbb{R}^{n}\right)$ taking finite values on $H^{1}$ and
$+\infty$ on $C\left(\mathbb{T};\mathbb{R}^{n}\right)\backslash H^{1}$. The
functional $\mathscr{U}\left(\varphi\right)$ is called regular if it is lower
semicontinuous on $C\left(\mathbb{T};\mathbb{R}^{n}\right)$ equipped with the
uniform convergence topology and the sets $\left\\{\varphi\in
C\left(\mathbb{T};\mathbb{R}^{n}\right):\left\|\varphi\right\|\leq
b,\mathscr{U}\left(\varphi\right)\leq a\right\\}$ are compact in
$C\left(\mathbb{T};\mathbb{R}^{n}\right)$ for any
$a,b\in\left(0,\infty\right)$.

We state Theorem 9 of Ref. freidlin1988random in the following

######  Proposition 4.1.

Assume that $\varphi_{0}\in C\left(\mathbb{T};\mathbb{R}^{n}\right)$ is an
asymptotically stable equilibrium point of (34) with $\varepsilon=0$. Let a
regular region, $D\subset C\left(\mathbb{T};\mathbb{R}^{n}\right)$, be such
that $\varphi_{0}\in D$ and every trajectory of (34) with $\varepsilon=0$
starting at a point $g\in D\cup\partial D$ does not leave $D$ for $t>0$ and
tends to $\varphi_{0}$ as $t\rightarrow\infty$. Furthermore, assume that there
is a regular functional, $\mathscr{U}\left(\varphi\right)$, and an operator,
$L\left(\varphi\right)=\left(L_{1}\left(\varphi\right),\ldots,L_{n}\left(\varphi\right)\right)$,
$\varphi\in H^{1}$, such that

  1. (i)

For $\varphi\in H^{2}$ the variational derivatives
$\delta\mathscr{U}\left(\varphi\right)/\delta\varphi_{k}$, $k=1,\ldots,n$, are
defined and

$\left(\nabla\mathscr{U}\left(\varphi\right),L\left(\varphi\right)\right)=\int_{\mathbb{T}}\sum_{k=1}^{n}\frac{\delta\mathscr{U}}{\delta\varphi_{k}}\left(\varphi\left(x\right)\right)L_{k}\left(\varphi\left(x\right)\right)dx=0,\textrm{
}\varphi\in H^{2}.$

  2. (ii)

For the field
$B\left(\varphi\right)=\left(B_{1}\left(\varphi\right),\ldots,B_{n}\left(\varphi\right)\right)$
we have

$B\left(\varphi\right)=-\nabla\mathscr{U}\left(\varphi\right)+L\left(\varphi\right),\textrm{
}\varphi\in H^{2}.$

  3. (iii)

For any $g\in H^{1}\cap\left(D\cup\partial D\right)$, there exists a function
$v\left(t,x\right)=\left(v_{1}\left(t,x\right),\ldots,v_{n}\left(t,x\right)\right)$,
$t>0$, $x\in\mathbb{T}$, such that

$\frac{\partial v_{k}\left(t,\cdot\right)}{\partial
t}=-\frac{\delta\mathscr{U}\left(v\left(t,\cdot\right)\right)}{\delta
v_{k}}-L_{k}\left(v_{k}\left(t,\cdot\right)\right),\textrm{
}t>0,k=1,\ldots,n,$ $v\left(0,x\right)=g\left(x\right),\textrm{
}\underset{t\rightarrow\infty}{\lim}\underset{x\in\mathbb{T}}{\sup}\left|v\left(t,x\right)-\varphi_{0}\left(x\right)\right|=0.$

Then for $g\in H^{1}\cap\left(D\cup\partial D\right)$

$\inf\left\\{S^{u}\left(\varphi\right),\varphi\left(0,x\right)=\varphi_{0}\left(x\right),\varphi\left(T,x\right)=g\left(x\right),T>0\right\\}=2\left(\mathscr{U}\left(g\right)-\mathscr{U}\left(\varphi_{0}\right)\right),$

and for any $g\in D$

$\underset{\varepsilon\rightarrow 0^{+}}{\lim}\textrm{ }\varepsilon^{2}\ln
E_{g}\tau^{\varepsilon}=2\underset{g\in\partial
D}{\min}\left(\mathscr{U}\left(g\right)-\mathscr{U}\left(\varphi_{0}\right)\right).$

###  4.2 Control theory for stochastic PDE

As already stated, we must consider an optimal control problem with
endpoint/state constraints. Unfortunately, this problem is not well understood
up to now. Hence, we will start with a closely related optimal control problem
for which we can show a Pontriagyn-type maximum principle. Afterward, we will
study the case with endpoint/state constraints using set-valued analysis and
additional restrictions.

We have, in addition to Assumptions 3.2-3.3, the following two assumptions for
the first case (Problem (26)):

###### Assumption 4.1

Let $a\left(\cdot,\cdot,\cdot\right):\left[0,T\right]\times\Omega\times
H\times U\rightarrow H$ and
$b\left(\cdot,\cdot,\cdot\right):\left[0,T\right]\times\Omega\times H\times
U\rightarrow\mathcal{L}_{2}^{0}$ satisfy:

  1. (i)

For any $\left(x,\alpha\right)\in H\times U$, the functions
$a\left(\cdot,x,\alpha\right):\left[0,T\right]\times\Omega\rightarrow H$ and
$b\left(\cdot,x,\alpha\right)\rightarrow\mathcal{L}_{2}^{0}$ are
$\mathbb{F}-$measurable,

  2. (ii)

For any $x\in H$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$, the functions
$a\left(t,x,\cdot\right):U\rightarrow H$ and
$b\left(t,x,\cdot\right):U\rightarrow\mathcal{L}_{2}^{0}$ are continuous,

  3. (iii)

For any $\left(x_{1},x_{2},\alpha\right)\in H\times H\times U$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$,

$\left\\{\begin{array}[]{l}\left|a\left(t,x_{1},\alpha\right)-a\left(t,x_{2},\alpha\right)\right|_{H}+\left|b\left(t,x_{1},\alpha\right)-b\left(t,x_{2},\alpha\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C\left|x_{1}-x_{2}\right|_{H},\\\
\left|a\left(t,0,\alpha\right)\right|_{H}+\left|b\left(t,0,\alpha\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C.\end{array}\right.$

###### Assumption 4.2

For a.e. $\left(t,\omega\right)\in\left(0,T\right)\times\Omega$, the functions
$a\left(t,\cdot,\cdot\right):H\times U\rightarrow H$,
$b\left(t,\cdot,\cdot\right):H\times U\rightarrow\mathcal{L}_{2}^{0}$,
$g\left(t,\cdot,\cdot\right):H\times U\rightarrow\mathbb{R}$, and
$h\left(\cdot\right):H\rightarrow\mathbb{R}$ are $C^{1}$. Furthermore, for any
$\left(x,\alpha\right)\in H\times U$ and a.e.
$\left(t,\omega\right)\in\left(0,T\right)\times\Omega$, we have

$\left\\{\begin{array}[]{l}\left|a_{x}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(H\right)}+\left|b_{x}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(H;\mathcal{L}_{2}^{0}\right)}+\left|g_{x}\left(t,x,\alpha\right)\right|_{H}+\left|h_{x}\left(x\right)\right|_{H}\leq
C,\\\
\left|a_{\alpha}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(\tilde{H};H\right)}+\left|b_{\alpha}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(\tilde{H};\mathcal{L}_{2}^{0}\right)}+\left|g_{\alpha}\left(t,x,\alpha\right)\right|_{\mathcal{L}\left(\tilde{H}\right)}\leq
C\end{array}\right.$

Under Assumption 4.1, Proposition 12.1 of Ref. lu2021mathematical ensures
that, for any $x_{0}\in L_{\mathcal{F}_{0}}^{p_{0}}\left(\Omega;H\right)$,
$p_{0}\geq 2$, and $\alpha\left(\cdot\right)\in\mathcal{U}\left[0,T\right]$,
system (3) has a unique mild solution, $x\left(\cdot\right)\equiv
x\left(\cdot;x_{0},\alpha\right)\in
C_{\mathbb{F}}\left(\left[0,T\right];L^{p_{0}}\left(\Omega;H\right)\right)$,
such that

$\left|x\left(\cdot\right)\right|_{C_{\mathbb{F}}\left(\left[0,T\right];L^{p_{0}}\left(\Omega;H\right)\right)}\leq
C\left(1+\left|x_{0}\right|_{L_{\mathcal{F}_{0}}^{p_{0}}\left(\Omega;H\right)}\right).$

Let Assumptions 3.2 and 4.1 hold. Then, by Theorem 4.19 of Ref.
lu2021mathematical, equation (27) is well-posed in the sense of transposition
solution (see Definition 4.13 of Ref. lu2021mathematical); if we consider that
the filtration $\mathbf{F}$ is the natural one and $y_{T}\in
L_{\mathcal{F}_{T}}^{p}$, $p\in\left(1,2\right]$, then the solution is mild
(see Section 4.2.1 of Ref. lu2021mathematical for the notions of solutions and
Section 4.2.2 of Ref. lu2021mathematical for the case of natural filtration).

We state Theorem 12.4 of Ref. lu2021mathematical in the following

######  Proposition 4.2.

Let Assumptions 3.2-3.3 and 4.1-4.2 hold. Let
$\left(\bar{x}\left(\cdot\right),\bar{\alpha}\left(\cdot\right)\right)$ be an
optimal pair for Problem (26) with $x_{0}\in
L_{\mathcal{F}_{0}}^{2}\left(\Omega;H\right)$. Then,

$\mathrm{Re}\left\langle
a_{u}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right)^{*}y\left(t\right)+b_{u}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right)^{*}Y\left(t\right)-g_{u}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right)\right),\right.\\\
\left.\alpha-\bar{\alpha}\left(t\right)\right\rangle_{\tilde{H}}\leq 0,$

a.e. $\left(t,\omega\right)\in\left[0,T\right]\times\Omega$, $\forall\alpha\in
U$, where $\left(y\left(\cdot\right),Y\left(\cdot\right)\right)$ is the
transposition solution of (27).

Now, we turn to the case of necessary optimality conditions for controlled
stochastic differential equations with control and state constraints. Here, we
follow closely the results presented in Ref. frankowska2020first.

We have, in addition to Assumptions 3.7-3.8, the following two assumptions for
Problem (31):

###### Assumption 4.3

$a\left(\cdot,\cdot,\cdot,\cdot\right):\left[0,T\right]\times H\times
H_{1}\times\Omega\rightarrow H$ and
$b\left(\cdot,\cdot,\cdot,\cdot\right):\left[0,T\right]\times H\times
H_{1}\times\Omega\rightarrow\mathcal{L}_{2}^{0}$ are two maps such that

  1. (i)

For any $\left(x,\alpha\right)\in H\times H_{1}$,
$a\left(\cdot,x,\alpha,\cdot\right):\left[0,T\right]\times\Omega\rightarrow H$
and
$b\left(\cdot,x,\alpha,\cdot\right):\left[0,T\right]\times\Omega\rightarrow\mathcal{L}_{2}^{0}$
are $\mathcal{B}\left(\left[0,T\right]\right)\times\mathcal{F}$ measurable and
$\mathbb{F}-$adapted,

  2. (ii)

For any $\left(t,x,\omega\right)\in\left[0,T\right]\times H\times\Omega$,
$a\left(t,x,\cdot,\omega\right):H_{1}\rightarrow H$ and
$b\left(t,x,\cdot,\omega\right):H_{1}\rightarrow\mathcal{L}_{2}^{0}$ are
continuous and

$\begin{cases}\left|a\left(t,x_{1},\alpha,\omega\right)-a\left(t,x_{2},\alpha,\omega\right)\right|_{H}+\left|b\left(t,x_{1},\alpha,\omega\right)-b\left(t,x_{2},\alpha,\omega\right)\right|_{\mathcal{L}_{2}^{0}}\leq\\\
\hfill
C\left|x_{1}-x_{2}\right|_{H}\quad\forall\left(t,x_{1},x_{2},\alpha,\omega\right)\in\left[0,T\right]\times
H\times H\times H_{1}\times\Omega\\\
\left|a\left(t,0,\alpha,\omega\right)\right|_{H}+\left|b\left(t,0,\alpha,\omega\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C,\quad\forall\left(t,\alpha,\omega\right)\in\left[0,T\right]\times
H_{1}\times\Omega.\end{cases}$

###### Assumption 4.4

For a.e. $\left(t,\omega\right)\in\left[0,T\right]\times\Omega$, the functions
$a\left(t,\cdot,\cdot,\omega\right):H\times H_{1}\rightarrow H$ and
$b\left(t,\cdot,\cdot,\omega\right):H\times
H_{1}\rightarrow\mathcal{L}_{2}^{0}$ are differentiable, and
$\left(a_{x}\left(t,x,\alpha,\omega\right),a_{\alpha}\left(t,x,\alpha,\omega\right)\right)$
and
$\left(b_{x}\left(t,x,\alpha,\omega\right),a_{\alpha}\left(t,x,\alpha,\omega\right)\right)$
are uniformly continuous with respect to $x\in H$ and $\alpha\in U$. There
exists a nonnegative $\eta\in L_{\mathbb{F}}^{2}\left(0,T;\mathbb{R}\right)$
such that for a.e. $\left(t,\omega\right)\in\left[0,T\right]\times\Omega$ and
for all $x\in H$ and $\alpha\in H_{1}$,

$\begin{cases}\left|a\left(t,0,\alpha,\omega\right)\right|_{H}+\left|b\left(t,0,\alpha,\omega\right)\right|_{\mathcal{L}_{2}^{0}}\leq
C\left(\eta\left(t,\omega\right)+\left|\alpha\right|_{H_{1}}\right),\\\
\left|a_{x}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H\right)}+\left|a_{\alpha}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H_{1};H\right)}+\left|b_{x}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H;\mathcal{L}_{2}^{0}\right)}\\\
\hfill+\left|b_{\alpha}\left(t,x,\alpha,\omega\right)\right|_{\mathcal{L}\left(H_{1};\mathcal{L}_{2}^{0}\right)}\leq
C.\end{cases}$

Under Assumption 4.3, Lemma 2.1 of Ref. frankowska2020first ensures that, for
any $x_{0}\in H$ and $\alpha\left(\cdot\right)\in\mathcal{U}_{2}$, system (3)
has a unique mild solution, $x\left(\cdot\right)\equiv
x\left(\cdot;x_{0},\alpha\right)\in
L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$, such
that

$\left|x\left(\cdot\right)\right|_{L^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)}\leq
C\left(1+\left|x_{0}\right|_{H}\right).$

We review some basic ideas from set-valued analysis needed to study (31); see
Ref. aubin2009set for details.

Let $Z$ be a Banach space and consider any subset $\mathcal{K}\subset Z$. For
$z\in\mathcal{K}$, the Clarke tangent cone
$\mathcal{C}_{\mathcal{K}}\left(z\right)$ to $\mathcal{K}$ at $z$ is defined
by

$\mathcal{C}_{\mathcal{K}}\left(z\right)=\left\\{v\in
Z:\underset{\varepsilon\rightarrow 0,y\in\mathcal{K},y\rightarrow
z}{\lim}\frac{dist\left(y+\varepsilon
v,\mathcal{K}\right)}{\varepsilon}=0\right\\},$

where
$dist\left(w,\mathcal{K}\right)=\underset{y\in\mathcal{K}}{\inf}\left|y-w\right|_{Z},$
$w\in Z$. Moreover, the adjacent cone $T_{\mathcal{K}}^{b}\left(z\right)$ to
$\mathcal{K}$ at $z$ is given by

$T_{\mathcal{K}}^{b}\left(z\right)=\left\\{v\in
Z:\underset{\varepsilon\rightarrow 0^{+}}{\lim}\frac{dist\left(y+\varepsilon
v,\mathcal{K}\right)}{\varepsilon}=0\right\\}.$

$\mathcal{C}_{\mathcal{K}}\left(z\right)$ is a closed convex cone in $Z$ and
$\mathcal{C}_{\mathcal{K}}\left(z\right)\subset
T_{\mathcal{K}}^{b}\left(z\right)$. If $\mathcal{K}$ is convex, then
$\mathcal{C}_{\mathcal{K}}\left(z\right)=T_{\mathcal{K}}^{b}\left(z\right)=\textrm{cl}\left\\{\alpha\left(\hat{z}-z\right)\textrm{
}:\textrm{ }z\geq 0,\textrm{ }\hat{z}\in\mathcal{K}\right\\}$, where
$\textrm{cl}\left\\{\cdot\right\\}$ represents the closure of the set.

For a cone $\mathcal{K}$ in $Z$, the closed convex cone
$\mathcal{K}^{-}=\left\\{\xi\in Z^{*}\textrm{ }:\textrm{
}\xi\left(z\right)\leq 0\textrm{ for all }z\in\mathcal{K}\right\\}$ is called
the dual cone of $\mathcal{K}$.

Let $\left(\varTheta,\varSigma\right)$ be a measurable space and
$F:\Theta\rightsquigarrow Z$ be a set-valued map. The domain of $F$ is
$\textrm{Dom}\left(F\right)=\left\\{\xi\in\Theta\textrm{ }:\textrm{
}F\left(\xi\right)\cap B\neq\emptyset\right\\}$. $F$ is measurable if
$F^{-1}\left(B\right)=\left\\{\xi\in\Theta\textrm{ }:\textrm{
}F\left(\xi\right)\cap B\neq\emptyset\right\\}\in\Sigma$ for any Borel set
$B\in\mathcal{B}\left(Z\right)$.

Assume that $\left(\Theta,\Sigma,\mu\right)$ is a complete $\sigma-$finite
measure space and $F$ is a set-valued map from $\Theta$ to the separable
Banach space $\tilde{Z}$ with nonempty closed image. Then, by Lemma 2.3,[14]
$F$ is measurable if and only if its graph belongs to
$\Sigma\otimes\mathcal{B}\left(\tilde{Z}\right)$.

A map, $\zeta:\left(\Omega,\mathcal{F}\right)\rightsquigarrow Z$, is a set
valued random variable if it is measurable. A map,
$\Psi:\left[0,T\right]\times\Omega\rightsquigarrow Z$, is a measurable set-
valued stochastic process if $\Psi$ is
$\mathcal{B}\left[0,T\right]\otimes\mathcal{F}-$measurable; $\Psi$ is adapted
if $\Psi\left(t,\cdot\right)$ is $\mathcal{F}_{t}-$measurable for all
$t\in\left[0,T\right]$.

Let

$\mathcal{G}=\left\\{B\in\mathcal{B}\left(\left[0,T\right]\right)\otimes\mathcal{F}\textrm{
}:\textrm{ }B_{t}\in\mathcal{F}_{t},\textrm{ }\forall
t\in\left[0,T\right]\right\\},$

where $B_{t}=\left\\{\omega\in\Omega\textrm{ }:\textrm{
}\left(t,\omega\right)\in B\right\\}$. Let m be the Lebesgue measure on
$\left[0,T\right]$. We consider the completion of the measure space
$\left(\left[0,T\right]\times\Omega,\mathcal{G},\mu=\textrm{m}\times
P\right)$, and due to Lemma 2.4,[14] we use the same notation for the
completion.

Let $H$ be a separable Hilbert space. By Lemma 2.5,[14] a set-valued
stochastic process, $F:\left[0,T\right]\times\Omega\rightsquigarrow H$, is
$\mathcal{B}\left(\left[0,T\right]\right)\otimes\mathcal{F}-$measurable and
$\mathbb{F}-$adapted if and only if $F$ is $\mathcal{G}-$measurable.

Let $\Phi$ be a set-valued stochastic process such that

  1. (i)

$\Phi$ in
$\mathcal{B}\left(\left[0,T\right]\right)\otimes\mathcal{F}-$measurable and
$\mathbb{F}-$adapted,

  2. (ii)

for a.e. $\left(t,\omega\right)\in\left[0,T\right]\times\Omega$,
$\Phi\left(t,\omega\right)$ is a nonempty closed convex cone in $H_{1}$.

  3. (iii)

$\Phi\left(t,\omega\right)\subset
T_{U}^{b}\left(\bar{\alpha}\left(t,\omega\right)\right)$, for a.e.
$\left(t,\omega\right)\in\left[0,T\right]\times\Omega$.

Define

$\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)=\left\\{\alpha\left(\cdot\right)\in
L_{\mathbb{F}}^{2}\left(0,T;H_{1}\right)\textrm{ }:\textrm{
}\alpha\left(t,\omega\right)\in\Phi\left(t,\omega\right),\textrm{ a.e.
}\left(t,\omega\right)\in\left[0,T\right]\times\Omega\right\\}.$

$\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)$ is a closed convex cone in
$L_{\mathbb{F}}^{2}\left(0,T;H_{1}\right)$. Moreover, since
$0\in\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)$,
$\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)$ is nonempty.

Under Assumption 4.3, for any
$\alpha_{1}\in\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)$ and $x_{1}\in
T_{\mathcal{K}_{a}}^{b}\left(\bar{x}_{0}\right)$, (32) has a unique solution,
$x_{1}\left(\cdot\right)\in
L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$; see
Section 3 of Ref. frankowska2020first. Furthermore, under Assumptions 4.3-4.4
and $\psi\in
L_{\mathbb{F}}^{2}\left(\Omega;BV_{0}\left(\left[0,T\right];H\right)\right)$,
Lemma 3.5 of Ref. frankowska2020first ensures that (33) has a unique
transposition solution $\left(y,Y\right)\in
D_{\mathbb{F}}\left(\left[0,T\right];L^{2}\left(\Omega;H\right)\right)\times
L_{\mathbb{F}}^{2}\left(0,T;\mathcal{L}_{2}^{0}\right)$; see Definition
3.1.[14]

Let $T_{\mathcal{K}_{a}}\left(\bar{x}_{0}\right)$ be a nonempty closed convex
cone contained in $T_{\mathcal{K}_{a}}^{b}\left(\bar{x}_{0}\right)$. Define

$\mathcal{G}^{\left(1\right)}=\left\\{x_{1}\left(\cdot\right)\in
L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)\textrm{
}:\textrm{ }\right.\\\ \left.x_{1}\left(\cdot\right)\textrm{ solves
(\ref{eq:Lin_sys}) with
}\alpha_{1}\in\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)\textrm{ and
}x_{1}\in T_{\mathcal{K}_{a}}\left(\bar{x}_{0}\right)\right\\},$
$\mathcal{I}_{0}\left(\bar{x}\right)=\left\\{t\in\left[0,T\right]\textrm{
}:\textrm{ }E\left[g^{0}\left(\bar{x}\left(t\right)\right)\right]=0\right\\},$
$\mathcal{I}\left(\bar{x}\right)=\left\\{j\in\left\\{1,\ldots,n\right\\}\textrm{
}:\textrm{ }E\left[g^{j}\left(\bar{x}\left(T\right)\right)=0\right]\right\\},$
$\mathcal{Q}^{\left(1\right)}=\left\\{z\left(\cdot\right)\in
L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)\textrm{
}:\textrm{ }E\left\langle
g_{x}^{0}\left(\bar{x}\left(t\right)\right),z\left(t\right)\right\rangle_{H}<0,\textrm{
}\forall t\in\mathcal{I}_{0}\left(\bar{x}\right)\right\\},$
$\mathcal{E}^{\left(1,j\right)}=\left\\{z\left(\cdot\right)\in
L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)\textrm{
}:\textrm{ }E\left\langle
g_{x}^{j}\left(\bar{x}\left(T\right)\right),z\left(T\right)\right\rangle_{H}<0\right\\},\textrm{
}\forall j\in\mathcal{I}\left(\bar{x}\right),$
$\mathcal{E}^{\left(1\right)}=\bigcap_{j\in\mathcal{I}\left(\bar{x}\right)}\mathcal{E}^{\left(1,j\right)}.$

Since $\mathcal{T}_{\Phi}\left(\bar{\alpha}\right)$ and
$\mathcal{T}_{\mathcal{K}_{a}}\left(\bar{x}_{0}\right)$ are nonempty convex
cones, $\mathcal{G}^{\left(1\right)}$ is a nonempty convex cone in
$L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$.
Moreover, if $\mathcal{I}_{0}\left(\bar{x}\right)=\emptyset$ (resp.
$\mathcal{I}\left(\bar{x}\right)=\emptyset$), then
$\mathcal{Q}^{\left(1\right)}=L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$
(resp.
$\mathcal{E}^{\left(1\right)}=L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$).
By Lemma 3.4,[14] $\mathcal{Q}^{\left(1\right)}$ is an open convex cone in
$L_{\mathbb{F}}^{2}\left(\Omega;C\left(\left[0,T\right];H\right)\right)$.

Define the Hamiltonian

$\mathbb{H}\left(t,x,\alpha,p,q,\omega\right)=\left\langle
p,a\left(t,x,\alpha,\omega\right)\right\rangle_{H}+\left\langle
q,b\left(t,x,\alpha,\omega\right)\right\rangle_{\mathcal{L}_{2}^{0}},$

where $\left(t,x,\alpha,p,q,\omega\right)\in\left[0,T\right]\times H\times
H_{1}\times H\times\mathcal{L}_{2}^{0}\times\Omega$.

Now, we state a first-order necessary optimality condition as presented in
Theorem 3.1.[14]

######  Proposition 4.3.

Let Assumptions 4.3-4.4 and 3.7-3.8 hold. Let
$\left(\bar{x}\left(\cdot\right),\bar{\alpha}\left(\cdot\right),\bar{x}_{0}\right)$
be an optimal triple of Problem (31). If
$E\left|g_{x}^{0}\left(\bar{x}\left(t\right)\right)\right|_{H}\neq 0$ for any
$t\in\mathcal{I}_{0}\left(\bar{x}\right)$, then there exists
$\lambda_{0}\in\left\\{0,1\right\\}$, $\lambda_{j}\geq 0$ for
$j\in\mathcal{I}\left(\bar{x}\right)$ and
$\psi\in\left(\mathcal{Q}^{\left(1\right)}\right)^{-}$ with
$\psi\left(0\right)=0$ verifying

$\lambda_{0}+\sum_{j\in\mathcal{I}\left(\bar{x}\right)}\lambda_{j}+\left|\psi\right|_{L_{\mathbb{F}}^{2}\left(\Omega;BV\left(0,T;H\right)\right)}\neq
0,$

such that the corresponding transposition solution
$\left(y\left(\cdot\right),Y\left(\cdot\right)\right)$ of the first order
adjoint equation (33) with
$y\left(T\right)=-\lambda_{0}h_{x}\left(\bar{x}\left(T\right)\right)-\sum_{j\in\mathcal{I}\left(\bar{x}\right)}\lambda_{j}g_{x}^{j}\left(\bar{x}\left(T\right)\right)$
satisfies the variational inequality

$E\left\langle
y\left(0\right),\nu\right\rangle_{H}+E\int_{0}^{T}\left\langle\mathbb{H}_{\alpha}\left[t\right],v\left(t\right)\right\rangle_{H_{1}}dt\leq
0,\textrm{
}\forall\nu\in\mathcal{T}_{\mathcal{K}_{a}}\left(\bar{x}_{0}\right),\textrm{
}\forall v\left(\cdot\right)\in\mathcal{T}_{\Phi}\left(\bar{\alpha}\right),$

where
$\mathbb{H}_{\alpha}\left[t\right]=\mathbb{H}_{\alpha}\left(t,\bar{x}\left(t\right),\bar{\alpha}\left(t\right),y\left(t\right),Y\left(t\right),\omega\right)$.
Furthermore, if
$\mathcal{G}^{\left(1\right)}\cap\mathcal{Q}^{\left(1\right)}\cap\mathcal{E}^{\left(1\right)}\neq\emptyset$,
the above holds with $\lambda_{0}=1$.

##  5 Proof of Results

######  Proof 5.1 (Theorem 3.1).

(i) and (ii) follow from the theory presented in subsection 4.1.

Let

$\Phi\left[u,v\right]=\int_{\mathcal{O}}L\left(u,v\right)dx,$

where $L=\frac{1}{2}\left(d_{1}\left|\partial u/\partial
x\right|^{2}+d_{2}\left|\partial v/\partial
x\right|^{2}\right)+F\left(u,v\right)$. Then, we can write the deterministic
part of (2), i.e., $\sigma_{1}=\sigma_{2}=0$, as

$\displaystyle\frac{\partial u}{\partial t}=$
$\displaystyle-\frac{\delta\Phi\left(u,v\right)}{\delta u},$
$\displaystyle\frac{\partial v}{\partial t}=$
$\displaystyle-\frac{\delta\Phi\left(u,v\right)}{\delta v},$

where $\delta$ stands for the variational derivative; this shows the steepest
descent (gradient flow) feature of the deterministic system. Note that
$\varphi_{0}\left(x\right)\equiv\hat{\varphi}^{\left(k\right)}$,
$k\in\left\\{1,\ldots,4\right\\}$ is an equilibrium point of the deterministic
epigenetic system. Moreover,
$\det\left(\mathbf{H}\left(F\right)\left(\hat{\varphi}^{k}\right)\right)=4a_{k}^{2}>0$
, $k=1,\ldots 4$; it follows that $\varphi_{0}$ is asymptotically stable. On
the other hand, Corollary 14.8 (b) in Ref. smoller2012shock ensures that the
phenotype (or cellular type) regions associated to each
$\hat{\varphi}^{\left(k\right)}$ are invariant.

Since our field is potential, i.e., $\tilde{f}=-\partial
F\left(u,v\right)/\partial u$, $\tilde{g}=-\partial F\left(u,v\right)/\partial
v$, then the field $B\left(\varphi\right)$ is also potential (see Section 3 of
Ref. freidlin1988random):

$B\left(\varphi\right)=-\nabla U\left(\varphi\right).$

Note that $U\left(\varphi\right)$ is finite on $H^{1}$. Following the ideas of
Lemma 3 in Ref. freidlin1988random we can show that if $U\left(\varphi\right)$
is extended onto $C\left(\mathbb{T};\mathbb{R}^{2}\right)$ with
$U\left(\varphi\right)=+\infty$ for $\varphi\in
C\left(\mathbb{T};\mathbb{R}^{2}\right)\backslash H^{1}$, then
$U\left(\varphi\right)$ is lower semicontinuous in
$C\left(\mathbb{T};\mathbb{R}^{2}\right)$. Moreover, for any
$a\in\left(0,\infty\right)$, the set $\Phi_{a}=\left\\{\varphi\in
C\left(\mathbb{T};\mathbb{R}^{2}\right)\textrm{ }:\textrm{
}U\left(\varphi\right)\leq a\right\\}$ is compact in
$C\left(\mathbb{T};\mathbb{R}^{2}\right)$ (Arzela-Ascoli theorem). Hence,
$U\left(\varphi\right)$ is a regular functional and (iii) follows from
Proposition 4.1.

######  Proof 5.2 (Theorem 3.2).

Both $\tilde{f}$ and $\tilde{g}$ are smooth with bounded derivatives of all
orders since we are considering the mollified version of $F$. The
corresponding (autonomous) Nemytskii operators on $H^{1}$ are Lipschitz and
continuously (Fréchet) differentiable; see Theorem 1.4 Ref.
efendiev2023linear. The result is also valid if we consider $H^{k}$, $k>1/2$;
see Ref. ambrosetti1995primer,appell1981implicit,appell1990nonlinear for
details. The $L^{2}$ case is more complicated, and the associated Nemytskii
operators are only Gateaux differentiable; see Theorem 2.7 Ref.
ambrosetti1995primer. Then, (i) and (ii) follow if we consider in addition
Assumptions 3.1-3.4; see Subsection 4.2.

######  Proof 5.3 (Theorem 3.3).

We proceed similarly to the proof of Theorem 3.2 for the regularity of
$\tilde{f}$ and $\tilde{g}$. Then, (i) and (ii) follow if we consider in
addition Assumptions 3.5-3.8.

##  6 Numerical Simulations

We now present some numerical simulations showing the dynamics of our
stochastic system. In particular, the exit from a basin of attraction and the
evolution afterward. Our results are based on the ideas presented in Chapter
10 of Ref. lord2014introduction.

We use the finite difference method for our stochastic epigenetic system with
homogeneous Neumann boundary conditions on $\mathcal{O}=\left(0,1\right)$. Set
$\sigma_{1}=\sigma_{2}=\sigma>0$ and let $W^{1}\left(t\right)$ and
$W^{2}\left(t\right)$ be two independent $Q-$Wiener processes on
$L^{2}\left(\mathcal{O}\right)$ with kernel
$q\left(x,y\right)=\exp\left(-\left|x-y\right|/l\right)$ for a correlation
length $l>0$; the white noise case is more difficult to handle, and we will
avoid it for our illustrative purposes.

Consider the grid points $x_{j}=jh$ for $h=1/J$ and $j=0,\ldots,J$. Let
$\boldsymbol{u}_{J}\left(t\right)$ and $\boldsymbol{v}_{J}\left(t\right)$ be
the finite difference approximations to
$\left[u\left(t,x_{1}\right),\ldots,u\left(t,x_{J-1}\right)\right]^{\mathrm{T}}$
and
$\left[v\left(t,x_{1}\right),\ldots,v\left(t,x_{J-1}\right)\right]^{\mathrm{T}}$,
respectively, resulting from the centered difference approximation $A^{N}$ of
the (negative of the) Laplacian (see eq. (3.51) in Ref. lord2014introduction).
That is, $\boldsymbol{u}_{J}\left(t\right)$ and
$\boldsymbol{v}_{J}\left(t\right)$ are the solution of

$\displaystyle d\boldsymbol{u}_{J}=$
$\displaystyle\left[-d_{1}A^{N}\boldsymbol{u}_{J}+\tilde{f}\left(\boldsymbol{u}_{J},\boldsymbol{v}_{J}\right)\right]dt+\sigma\boldsymbol{W}_{J}^{1}\left(t\right),$
$\displaystyle d\boldsymbol{v}_{J}=$
$\displaystyle\left[-d_{2}A^{N}\boldsymbol{v}_{J}+\tilde{g}\left(\boldsymbol{u}_{J},\boldsymbol{v}_{J}\right)\right]dt+\sigma\boldsymbol{W}_{J}^{2}\left(t\right),$

with
$\boldsymbol{u}_{J}\left(0\right)=\left[u_{0}\left(x_{1}\right),\ldots,u_{0}\left(x_{J-1}\right)\right]^{\mathrm{T}}$,
$\boldsymbol{v}_{J}\left(0\right)=\left[v_{0}\left(x_{1}\right),\ldots,v_{0}\left(x_{J-1}\right)\right]^{\mathrm{T}}$,
and
$\boldsymbol{W}_{J}^{k}=\left[W^{k}\left(t,x_{1}\right),\ldots,W^{k}\left(t,x_{J-1}\right)\right]^{\mathrm{T}}$,
$k=1,2.$ To discretize in time, we apply the semi-implicit Euler-Maruyama
method (see eq. (8.121) in Ref. lord2014introduction) with time step
$\triangle t>0$, which gives the approximations $\boldsymbol{u}_{J,n}$ to
$\boldsymbol{u}_{J}\left(t_{n}\right)$ and $\boldsymbol{v}_{J,n}$ to
$\boldsymbol{v}_{J}\left(t_{n}\right)$ at $t_{n}=n\triangle t$ defined by

$\displaystyle\boldsymbol{u}_{J,n+1}=$ $\displaystyle\left(I+\triangle
t\textrm{
}d_{1}A^{N}\right)^{-1}\left[\boldsymbol{u}_{J,n}+\tilde{f}\left(\boldsymbol{u}_{J,n},\boldsymbol{v}_{J,n}\right)\triangle
t+\sigma\triangle\boldsymbol{W}_{n}^{1}\right],$
$\displaystyle\boldsymbol{v}_{J,n+1}=$ $\displaystyle\left(I+\triangle
t\textrm{
}d_{2}A^{N}\right)^{-1}\left[\boldsymbol{v}_{J,n}+\tilde{g}\left(\boldsymbol{u}_{J,n},\boldsymbol{v}_{J,n}\right)\triangle
t+\sigma\triangle\boldsymbol{W}_{n}^{2}\right],$

with $\boldsymbol{u}_{J,0}=\boldsymbol{u}_{J}\left(0\right)$,
$\boldsymbol{v}_{J,0}=\boldsymbol{v}_{J}\left(0\right)$, and
$\triangle\boldsymbol{W}_{n}^{k}=\boldsymbol{W}_{J}^{k}\left(t_{n+1}\right)-\boldsymbol{W}_{J}^{k}\left(t_{n}\right)$,
$k=1,2$. Furthermore, $\boldsymbol{W}_{J}^{k}\left(t\right)\sim
N\left(\mathbf{0},tC\right)$, where $C$ is the matrix with entries
$q\left(x_{i},x_{j}\right)$ for $i,j=1,\ldots,J-1$. We use the circulant
embedding method to generate the increments $\triangle\boldsymbol{W}_{n}^{k}$
(see Algorithms 6.9 and 10.7 in Ref. lord2014introduction). For the stability
and convergence analysis of this method, see Chapter 10 in Ref.
lord2014introduction.

We use Python for our implementation and the function gaussian_filter of
scipy.ndimage to smooth out the epigenetic potential (1). Fig. (1) shows the
colors corresponding to the four organs in the flower arabidopsis-thaliana.
Fig (2) shows the epigenetic landscape and the system’s evolution (sample
path) for different values of $t$ (see also the animation in the online
version); the white dot at $t=0.00$ corresponds to the initial condition. We
compute the $L^{2}-$average (Avg) of the sample path at each time:

$Avg\textrm{
}u\left(t\right)=\sqrt{\frac{1}{\left|\mathcal{O}\right|}\int_{\mathcal{O}}\left|u\left(t,x\right)\right|^{2}dx},\textrm{
}Avg\textrm{
}v\left(t\right)=\sqrt{\frac{1}{\left|\mathcal{O}\right|}\int_{\mathcal{O}}\left|v\left(t,x\right)\right|^{2}dx},$
(38)

where $\left|\mathcal{O}\right|$ is the length of $\mathcal{O}$. We use the
trapezoidal method for the integrals and the initial conditions are constant
functions.

Fig (2) shows that, due to the effect of the noise, the system eventually
exits the local minima, traversing the epigenetic landscape in the spatial
order that corresponds to the correct architecture of the flower, that is,
following the observed geometrical features of the meristem
($sepals\rightarrow petals\rightarrow stamens\rightarrow carpels$). We remark
that both the depths and the locations of the minima of the basins of
attraction play a crucial role in describing the correct epigenetic dynamics
of a biological system like arabidopsis-thaliana. To see this, Fig 3 displays
the system’s evolution for a different arrangement of the epigenetic
potential, which is not the one expected from observed mutations of
arabidopsis-thaliana. Furthermore, these numerical simulations show the
significance of our rigorous results and how they can help to study other
biological systems (see Theorem 3.1).

Figs (4)-(5) show the distribution of Avg $u\left(t\right)$ and Avg
$v\left(t\right)$ (see (38)) of a sample path up to time $t=T$ with the same
initial condition as in Fig (2) and for different values of $\sigma$ (see also
the animations in the online version). Notice that, as the value of $\sigma$
is halved, so is the standard deviation of the distribution. Moreover, the
mean of the distributions is equal to the (Avg of the) minimum of the basin of
attraction (sepals in this case). Both observations agree with our results in
Theorem 3.1.

Figure 1: Colors corresponding to the four organs in the flower arabidopsis-thaliana. |
---|---
|
|
Figure 2: Epigenetic landscape and system’s evolution (sample path) for different values of $t$ with $\sigma=0.0140$ and $d_{1}=d_{2}=1$. The white dot at $t=0.00$ corresponds to the initial condition and Avg is given by (38) |
---|---
|
|
Figure 3: Epigenetic landscape (different configuration; cf. Fig (2)) and
system’s evolution (sample path) for different values of $t$ with
$\sigma=0.0140$ and $d_{1}=d_{2}=1$. The white dot at $t=0.00$ corresponds to
the initial condition and Avg is given by (38)
---
Figure 4: Distribution of Avg $u\left(t\right)$ and Avg $v\left(t\right)$ (see
(38)) of a sample path up to time $t=T$ with the same initial condition as in
Fig (2) and for $\sigma=0.0120$ (top), $\sigma=0.0060$ (bottom).
---
Figure 5: Distribution of Avg $u\left(t\right)$ and Avg $v\left(t\right)$ (see
(38)) of a sample path up to time $t=T$ with the same initial condition as in
Fig (2) and for $\sigma=0.0030$ (top), $\sigma=0.0015$ (bottom).

## References

  * [1] E. Alvarez-Buylla, A. Chaos, M. Aldana, M. Benítez, Y. Cortes-Poza, C. Espinosa-Soto, D.A. Hartasánchez, R.B. Lotto, D. Malkin, G.J. Escalera Santos, et al.  Floral morphogenesis: stochastic explorations of a gene network epigenetic landscape.  Plos one, 3(11):e3626, 2008.
  * [2] A. Ambrosetti and G. Prodi.  A primer of nonlinear analysis.  Number 34. Cambridge University Press, 1995.
  * [3] J. Appell.  Implicit functions, nonlinear integral equations, and the measure of noncompactness of the superposition operator.  Journal of Mathematical Analysis and Applications, 83(1):251–263, 1981.
  * [4] J. Appell and P. Zabrejko.  Nonlinear superposition operators.  Cambridge university press, 1990.
  * [5] K.J. Åström.  Introduction to stochastic control theory.  Courier Corporation, 2012.
  * [6] J.P. Aubin and H. Frankowska.  Set-valued analysis.  Springer Science & Business Media, 2009.
  * [7] Y. Cortés-Poza and P. Padilla-Longoria.  A variational approach to morphogenesis: Recovering spatial phenotypic features from epigenetic landscapes.  Bulletin of Mathematical Biology, 84(3):33, 2022.
  * [8] G. Da Prato and J. Zabczyk.  Ergodicity for infinite dimensional systems, volume 229.  Cambridge university press, 1996.
  * [9] G. Da Prato and J. Zabczyk.  Stochastic equations in infinite dimensions.  Cambridge university press, 2014.
  * [10] M. Efendiev.  Linear and Nonlinear Non-Fredholm Operators: Theory and Applications.  Springer Nature, 2023.
  * [11] M. Esteller.  Epigenetics in cancer.  New England Journal of Medicine, 358(11):1148–1159, 2008.
  * [12] K.A. Feldmann and S.A. Goff.  The first plant genome sequence-arabidopsis thaliana.  In Advances in Botanical Research, volume 69, pages 91–117. Elsevier, 2014.
  * [13] W.H. Fleming and R.W. Rishel.  Deterministic and stochastic optimal control, volume 1.  Springer Science & Business Media, 2012.
  * [14] H. Frankowska and Q. Lü.  First and second order necessary optimality conditions for controlled stochastic evolution equations with control and state constraints.  Journal of Differential Equations, 268(6):2949–3015, 2020.
  * [15] M. Freidlin.  Random perturbations of reaction-diffusion equations: the quasideterministic approximation.  Transactions of the American Mathematical Society, 305(2):665–697, 1988.
  * [16] M.I. Freidlin and A.D. Wentzell.  Random Perturbations of Dynamical Systems.  Springer, 1998.
  * [17] S.J. Gould.  The structure of evolutionary theory.  Harvard university press, 2002.
  * [18] M. Hairer.  Ergodic theory for stochastic PDEs.  preprint, 2008.
  * [19] Emory Winship Cancer Institute.  Cancer epigenetics.  https://cancerquest.org/cancer-biology/cancer-epigenetics, October 2023.
  * [20] L.S. Kristensen, H.M. Nielsen, and L.L. Hansen.  Epigenetics and cancer treatment.  European journal of pharmacology, 625(1-3):131–142, 2009.
  * [21] G.J. Lord, C.E. Powell, and T. Shardlow.  An introduction to computational stochastic PDEs, volume 50.  Cambridge University Press, 2014.
  * [22] Q. Lü and X. Zhang.  Mathematical control theory for stochastic partial differential equations.  Springer, 2021.
  * [23] S. Sharma, T.K. Kelly, and P.A. Jones.  Epigenetics in cancer.  Carcinogenesis, 31(1):27–36, 2010.
  * [24] J. Smoller.  Shock waves and reaction-diffusion equations, volume 258.  Springer Science & Business Media, 2012.
  * [25] S.R.S. Varadhan.  Large deviations, volume 27.  American Mathematical Soc., 2016.