A deterministic method for estimating free energy genetic network landscapes with applications to cell commitment and reprogramming paths

Depicting developmental processes as movements in free energy genetic landscapes is an illustrative tool. However, exploring such landscapes to obtain quantitative or even qualitative predictions is hampered by the lack of free energy functions corresponding to the biochemical Michaelis–Menten or Hill rate equations for the dynamics. Being armed with energy landscapes defined by a network and its interactions would open up the possibility of swiftly identifying cell states and computing optimal paths, including those of cell reprogramming, thereby avoiding exhaustive trial-and-error simulations with rate equations for different parameter sets. It turns out that sigmoidal rate equations do have approximate free energy associations. With this replacement of rate equations, we develop a deterministic method for estimating the free energy surfaces of systems of interacting genes at different noise levels or temperatures. Once such free energy landscape estimates have been established, we adapt a shortest path algorithm to determine optimal routes in the landscapes. We explore the method on three circuits for haematopoiesis and embryonic stem cell development for commitment and reprogramming scenarios and illustrate how the method can be used to determine sequential steps for onsets of external factors, essential for efficient reprogramming.


The synthetic single-gene motif
The free energy representing this system ( Figure 2A) is given by F (X; T ) = − 1 2 ω X,X X 2 + 1 2 (l X X) where X is the expression level of gene X. The corresponding dynamical equation is given by Here k is the degradation rate. Parameters values are found below in Table S1. Table S1: Parameter values for the synthetic single-gene motif ( Figure 2A).
Parameter Process Value ωX,X self-activation of X 1.50 lX external signal on X 1.50 k degradation 1.00 T low low temperature 0.10 T high high temperature 0.50

Sensitivity Analysis
We conducted a simple sensitivity analysis for the synthetic single-gene switch. First, we calculated the free energy for the parameters values in Table S1 with the T high value (blue lines in Figure S1). The magenta lines in Figure S1 show how variations of the self interaction strength ω X,X affect the free energy landscape leading to modifications of the two attractors heights. When the self interaction strength is increased the high concentration attractors is most affected while when it is decrease the low concentration attractor height increases. The opposite effect was obtained when we varied the external signal strength l X (orange lines in figure below). We also varied the decay rate k and observed major modification of the free energy as shown in Figure S1 with green lines. The sensitivity analysis shows that decay rates, self interactions and external signal strength have a major impact on the free energy landscape and slight variations of these parameters are able to influence cell fate. Energy ω X,X =1.4 ω X,X =1.5 ω X,X =1.6 k=0.3 k=1 k=1.6 l X =1.4 l X =1.5 l X =1.6 Figure S1: Sensitivity analysis for the synthetic single-gene switch. The free energy when the self interaction strength ωX,X is varied (magenta), decay rate k is varied (green) and the external signal strength lX is varied (orange). The blue curves show the free energy for the parameters values in Table S1. All graphs are for T high .

The synthetic two-gene mutual repressor motif
The free energy representing this system ( Figure 2C) is given by The corresponding dynamical equations are given by assuming the same degradation rate k for both genes. X and Y are the expression levels of genes X and Y respectively. Parameter values are found below in Table S2. Table S2: Parameter values for the synthetic two-gene mutual repressor motif ( Figure 2C). The hematopoietic gene regulatory network topology shown in Figure 3A was first introduced in [1]. We provide our corresponding free energy model with its associated sigmoidal rate equations. All links in the network are described by single parameters. Activations and repressions correspond to positive and negative parameter values respectively. The free energy representing the hematopoietic motif in Figure 3A is then given by: where G2, G1, and Gb are the expression levels of Gata2, Gfi1, and Gfi1b respectively. The corresponding dynamical equations are given by assuming the same degradation rate k for Gata2, Gfi1, and Gfi1b. Parameter values are found below in Table S3.

The Gata2-Gata1-Pu.1 regulatory motif
The hematopoietic motif shown in Figure 4A was first introduced in [2]. We provide our corresponding free energy model with its associated sigmoidal rate equations. All links in the network are described by single parameters. Activations and repressions correspond to positive and negative parameter values respectively. The free energy representing the hematopoietic motif in Figure 4A is then given by: where G2, G1, and P 1 are the expression levels of Gata2, Gata1, and Pu.1 respectively and EP O the Erythropoietin external signal. The corresponding dynamical equations are given by assuming the same degradation rate k for Gata2, Gata1, and Pu.1. Parameter values are found below in Table S4.

The Oct4/Sox2-Nanog-Fgf4-G regulatory motif
The topology of the gene regulatory network in Figure 5A is retrieved from [3] to which we added the experimentally proven Nanog self-repression [4]. We provide our corresponding free energy model with its associated sigmoidal rate equations. All links in the network are described by single parameters. Activations and repressions correspond to positive and negative parameter values respectively. The free energy representing the embryonic stem cell motif in Figure 5A is then given by: where N , O, F , and G are the expression levels of Nanog, Oct4-Sox2, Fgf4-Gsk3 and G respectively and L and I are the external signals Lif-Bmp4 and 2i-3i. The corresponding dynamical equations are given by assuming the same degradation rate k for Nanog, Oct4-Sox2, Fgf4-Gsk3, and G. Parameter values are found below in Table S5.