User Tools

Site Tools


examples:differential_equations

Ordinary differential equation models

ODE model: Cell cycle

Time plots of ODE model of Xenopus embryonic cell cycle

Introduction

This model is a simple three-species ODE model of the Xenopus embryonic cell cycle (Ferrell et al., 2011). It exhibits sustained limit cycle oscillations.

Model description

One CellType is created that has three variables of Properties representing the concentrations of APC, Plk1, and CDK1. These variables are coupled in a System of DiffEqns.

In the System, a number of Constants are defined whose symbols are used in the DiffEqn. Note that the equations are entered in simple plain text.

The System uses the runga-kutta (4th order) solver for the differential equations and speficies a particular time step (here $ht = 10^{-2}$) which is interpreted in global time steps.

The global time is defined in the Time element and runs from StartTime to StopTime ($0 - 25$). In this non-spatial model, Space defines a Lattice of size $(x,y,z)=(1,0,0)$.

Results are written to a file using the Analysis plugin Logger. The Logger also visualizes the time plot to screen (in “interactive” mode) or to PNG files (in “local” mode).

Things to try

  • Change the dynamics by altering System/time-scaling.

Model

CellCycle.xml

<?xml version='1.0' encoding='UTF-8'?>
<MorpheusModel version="3">
    <Description>
        <Title>Example-CellCycle</Title>
        <Details>James Ferrell, Tony Yu-Chen Tsai and Qiong Yang (2011) Modeling the Cell Cycle: Why Do Certain Circuits Oscillate?, Cell 144, p874-885. http://dx.doi.org/10.1016/j.cell.2011.03.006</Details>
    </Description>
    <Global/>
    <Space>
        <Lattice class="square">
            <Size symbol="size" value="1,1,0"/>
            <Neighborhood>
                <Order>1</Order>
            </Neighborhood>
        </Lattice>
        <SpaceSymbol symbol="space"/>
    </Space>
    <Time>
        <StartTime value="0"/>
        <StopTime value="25"/>
        <TimeSymbol symbol="time"/>
    </Time>
    <CellTypes>
        <CellType class="biological" name="cells">
            <Property symbol="APC" value="0"/>
            <Property symbol="Plk1" value="0"/>
            <Property symbol="CDK1" value="0"/>
            <System solver="runge-kutta" time-scaling="1" time-step="1e-2">
                <Constant symbol="n" value="8"/>
                <Constant symbol="K" value="0.5"/>
                <Constant symbol="α1" value="0.1"/>
                <Constant symbol="α2" value="3.0"/>
                <Constant symbol="α3" value="3.0"/>
                <Constant symbol="β1" value="3.0"/>
                <Constant symbol="β2" value="1.0"/>
                <Constant symbol="β3" value="1.0"/>
                <DiffEqn symbol-ref="CDK1">
                    <Expression>α1 - β1 * CDK1 * (APC^n) / (K^n + APC^n)</Expression>
                </DiffEqn>
                <DiffEqn symbol-ref="Plk1">
                    <Expression>α2*(1-Plk1) * ((CDK1^n) / (K^n + CDK1^n)) - β2*Plk1</Expression>
                </DiffEqn>
                <DiffEqn symbol-ref="APC">
                    <Expression>α3*(1- APC) * ((Plk1^n) / (K^n + Plk1^n)) - β3*APC</Expression>
                </DiffEqn>
            </System>
        </CellType>
    </CellTypes>
    <CellPopulations>
        <Population size="0" type="cells">
            <InitCellLattice/>
        </Population>
    </CellPopulations>
    <Analysis>
        <Logger time-step="1e-2">
            <Restriction>
                <Celltype celltype="cells"/>
            </Restriction>
            <Input>
                <Symbol symbol-ref="APC"/>
                <Symbol symbol-ref="CDK1"/>
                <Symbol symbol-ref="Plk1"/>
            </Input>
            <Output>
                <TextOutput file-format="csv"/>
            </Output>
            <Plots>
                <Plot time-step="-1">
                    <Style style="lines" line-width="4.0"/>
                    <Terminal terminal="png"/>
                    <X-axis>
                        <Symbol symbol-ref="time"/>
                    </X-axis>
                    <Y-axis>
                        <Symbol symbol-ref="APC"/>
                        <Symbol symbol-ref="CDK1"/>
                        <Symbol symbol-ref="Plk1"/>
                    </Y-axis>
                </Plot>
                <Plot time-step="-1">
                    <Style style="lines" line-width="4.0"/>
                    <Terminal terminal="png"/>
                    <X-axis>
                        <Symbol symbol-ref="CDK1"/>
                    </X-axis>
                    <Y-axis>
                        <Symbol symbol-ref="APC"/>
                    </Y-axis>
                    <Color-bar>
                        <Symbol symbol-ref="Plk1"/>
                    </Color-bar>
                </Plot>
            </Plots>
        </Logger>
    </Analysis>
</MorpheusModel>


In Morpheus GUI: Examples → ODE → CellCycle.xml.

Reference

Ferrell JE Jr, Tsai TY, Yang Q. Modeling the cell cycle: why do certain circuits oscillate? Cell, 18:144(6), 2011.



SBML import: MAPK signaling

Time plots of oscillations in MAPK signaling (Kholodenko, 2000)

Introduction

This model has been imported and automatically converted from SBML format by Morpheus. It shows oscillations in the MAPK signaling cascade (Kholodenko, 2000).

Model description

Upon importing an SBML file, a Morpheus model is automatically created. A System of DiffEqns is generated, based on the function and reactions defined in the SBML file and defined as part of a CellType. Additionally, a Logger is generated to record and visualize the output.

Simulation details, such as StartTime and StopTime, as well as the time-step of System, need to be specified manually.

Things to try

Model

Original SBML model:

Failed to retrieve: http://imc.zih.tu-dresden.de/morpheus/examples/ODE/BIOMD0000000010.xml

Generated Morpheus model:

MAPK_SBML.xml

<?xml version='1.0' encoding='UTF-8'?>
<MorpheusModel version="3">
    <Description>
        <Title>BIOMD0000000010</Title>
        <Details>Imported SBML model: http://www.ebi.ac.uk/biomodels-main/BIOMD0000000010

Kholodenko BN, Negative feedback and ultrasensitivity can bring about oscillations in the mitogen-activated protein kinase cascades.</Details>
    </Description>
    <Global/>
    <Space>
        <Lattice class="linear">
            <Size symbol="size" value="1 0 0"/>
            <Neighborhood>
                <Order>1</Order>
            </Neighborhood>
        </Lattice>
        <SpaceSymbol symbol="space"/>
    </Space>
    <Time>
        <StartTime value="0"/>
        <StopTime value="1.0"/>
        <TimeSymbol symbol="time" name="time"/>
    </Time>
    <CellTypes>
        <CellType class="biological" name="sbml_ct">
            <System solver="runge-kutta" time-scaling="4000" time-step="0.01">
                <Constant symbol="V1" value="2.5"/>
                <Constant symbol="Ki" value="9"/>
                <Constant symbol="n" value="1"/>
                <Constant symbol="K1" value="10"/>
                <DiffEqn symbol-ref="MKKK" name="gained from reactions">
                    <Expression> - uVol * V1 * MKKK / ((1 + pow(MAPK_PP / Ki, n)) * (K1 + MKKK)) + uVol * V2 * MKKK_P / (KK2 + MKKK_P)</Expression>
                </DiffEqn>
                <DiffEqn symbol-ref="MKKK_P" name="gained from reactions">
                    <Expression> uVol * V1 * MKKK / ((1 + pow(MAPK_PP / Ki, n)) * (K1 + MKKK)) - uVol * V2 * MKKK_P / (KK2 + MKKK_P)</Expression>
                </DiffEqn>
                <Constant symbol="V2" value="0.25"/>
                <Constant symbol="KK2" value="8"/>
                <Constant symbol="k3" value="0.025"/>
                <Constant symbol="KK3" value="15"/>
                <DiffEqn symbol-ref="MKK" name="gained from reactions">
                    <Expression> - uVol * k3 * MKKK_P * MKK / (KK3 + MKK) + uVol * V6 * MKK_P / (KK6 + MKK_P)</Expression>
                </DiffEqn>
                <DiffEqn symbol-ref="MKK_P" name="gained from reactions">
                    <Expression> uVol * k3 * MKKK_P * MKK / (KK3 + MKK) - uVol * k4 * MKKK_P * MKK_P / (KK4 + MKK_P) + uVol * V5 * MKK_PP / (KK5 + MKK_PP) - uVol * V6 * MKK_P / (KK6 + MKK_P)</Expression>
                </DiffEqn>
                <Constant symbol="k4" value="0.025"/>
                <Constant symbol="KK4" value="15"/>
                <DiffEqn symbol-ref="MKK_PP" name="gained from reactions">
                    <Expression> uVol * k4 * MKKK_P * MKK_P / (KK4 + MKK_P) - uVol * V5 * MKK_PP / (KK5 + MKK_PP)</Expression>
                </DiffEqn>
                <Constant symbol="V5" value="0.75"/>
                <Constant symbol="KK5" value="15"/>
                <Constant symbol="V6" value="0.75"/>
                <Constant symbol="KK6" value="15"/>
                <Constant symbol="k7" value="0.025"/>
                <Constant symbol="KK7" value="15"/>
                <DiffEqn symbol-ref="MAPK" name="gained from reactions">
                    <Expression> - uVol * k7 * MKK_PP * MAPK / (KK7 + MAPK) + uVol * V10 * MAPK_P / (KK10 + MAPK_P)</Expression>
                </DiffEqn>
                <DiffEqn symbol-ref="MAPK_P" name="gained from reactions">
                    <Expression> uVol * k7 * MKK_PP * MAPK / (KK7 + MAPK) - uVol * k8 * MKK_PP * MAPK_P / (KK8 + MAPK_P) + uVol * V9 * MAPK_PP / (KK9 + MAPK_PP) - uVol * V10 * MAPK_P / (KK10 + MAPK_P)</Expression>
                </DiffEqn>
                <Constant symbol="k8" value="0.025"/>
                <Constant symbol="KK8" value="15"/>
                <DiffEqn symbol-ref="MAPK_PP" name="gained from reactions">
                    <Expression> uVol * k8 * MKK_PP * MAPK_P / (KK8 + MAPK_P) - uVol * V9 * MAPK_PP / (KK9 + MAPK_PP)</Expression>
                </DiffEqn>
                <Constant symbol="V9" value="0.5"/>
                <Constant symbol="KK9" value="15"/>
                <Constant symbol="V10" value="0.5"/>
                <Constant symbol="KK10" value="15"/>
            </System>
            <Constant symbol="uVol" value="1" name="compartment size"/>
            <Property symbol="MKKK" value="90" name="Mos"/>
            <Property symbol="MKKK_P" value="10" name="Mos-P"/>
            <Property symbol="MKK" value="280" name="Mek1"/>
            <Property symbol="MKK_P" value="10" name="Mek1-P"/>
            <Property symbol="MKK_PP" value="10" name="Mek1-PP"/>
            <Property symbol="MAPK" value="280" name="Erk2"/>
            <Property symbol="MAPK_P" value="10" name="Erk2-P"/>
            <Property symbol="MAPK_PP" value="10" name="Erk2-PP"/>
        </CellType>
    </CellTypes>
    <CellPopulations>
        <Population size="1" type="sbml_ct"/>
    </CellPopulations>
    <Analysis>
        <Logger time-step="0.001">
            <Restriction>
                <Celltype celltype="sbml_ct"/>
            </Restriction>
            <Input>
                <Symbol symbol-ref="MAPK"/>
                <Symbol symbol-ref="MAPK_P"/>
                <Symbol symbol-ref="MAPK_PP"/>
                <Symbol symbol-ref="MKK"/>
                <Symbol symbol-ref="MKK_P"/>
                <Symbol symbol-ref="MKK_PP"/>
                <Symbol symbol-ref="MKKK"/>
                <Symbol symbol-ref="MKKK_P"/>
            </Input>
            <Output>
                <TextOutput/>
            </Output>
            <Plots>
                <Plot time-step="-1">
                    <Style grid="true" style="lines" line-width="3.0"/>
                    <Terminal terminal="png"/>
                    <X-axis>
                        <Symbol symbol-ref="time"/>
                    </X-axis>
                    <Y-axis>
                        <Symbol symbol-ref="MAPK"/>
                        <Symbol symbol-ref="MAPK_P"/>
                        <Symbol symbol-ref="MAPK_PP"/>
                        <Symbol symbol-ref="MKK"/>
                        <Symbol symbol-ref="MKK_P"/>
                        <Symbol symbol-ref="MKK_PP"/>
                        <Symbol symbol-ref="MKKK"/>
                        <Symbol symbol-ref="MKKK_P"/>
                    </Y-axis>
                </Plot>
            </Plots>
        </Logger>
    </Analysis>
</MorpheusModel>


In Morpheus GUI: Examples → ODE → MAPK_SBML.xml.

Reference

Delay differential equations: Cell cycle

Time plots of ODE model of Xenopus embryonic cell cycle, modeled with delay differential equations

Introduction

This model is a two-species version of the Xenopus embryonic cell cycle shown above that uses delay differential equations (Ferrell et al., 2011). It exhibits sustained limit cycle oscillations.

Model description

This model uses two Properties (CDK1 and APC) and two DelayProperties (CDK1_d and APC_d) with delay $\tau$. The latter are properties that return the value that has been assigned at time $t-\tau$.

The updated values of CDK1 and APC are assigned to (the back of) CDK1_d and APC_d using Equations. When these properties used in the DiffEqn, they return the value assigned in the past.

The two variables are logged and both a time plot and a phase plot are drawn.

Things to try

  • Explore the effect of delays by altering the DelayProperty/delay.

Model

CellCycleDelay.xml

<?xml version='1.0' encoding='UTF-8'?>
<MorpheusModel version="3">
    <Description>
        <Title>Example-CellCycleDelay</Title>
        <Details>Example of delay differential equations.

Implements equation 23 and 24 and reproduces figure 7 from:

James Ferrell, Tony Yu-Chen Tsai and Qiong Yang (2011) Modeling the Cell Cycle: Why Do Certain Circuits Oscillate?, Cell 144, p874-885. http://dx.doi.org/10.1016/j.cell.2011.03.006</Details>
    </Description>
    <Global/>
    <Space>
        <Lattice class="linear">
            <Size symbol="size" value="1 0 0"/>
            <Neighborhood>
                <Order>1</Order>
            </Neighborhood>
        </Lattice>
        <SpaceSymbol symbol="space"/>
    </Space>
    <Time>
        <StartTime value="0"/>
        <StopTime value="25"/>
        <SaveInterval value="0"/>
        <TimeSymbol symbol="time"/>
    </Time>
    <CellTypes>
        <CellType class="biological" name="cells">
            <Property symbol="APC" value="0"/>
            <Property symbol="CDK1" value="0"/>
            <DelayProperty symbol="APC_d" value="0" delay="0.5"/>
            <DelayProperty symbol="CDK1_d" value="0" delay="0.5"/>
            <Equation symbol-ref="APC_d">
                <Expression>APC</Expression>
            </Equation>
            <Equation symbol-ref="CDK1_d">
                <Expression>CDK1</Expression>
            </Equation>
            <System solver="runge-kutta" time-scaling="1" time-step="1e-2">
                <Constant symbol="n" value="8"/>
                <Constant symbol="K" value="0.5"/>
                <Constant symbol="α1" value="0.1"/>
                <Constant symbol="α2" value="3.0"/>
                <Constant symbol="β1" value="3.0"/>
                <Constant symbol="β2" value="1.0"/>
                <DiffEqn symbol-ref="CDK1">
                    <Expression>α1 - β1 * CDK1 * (APC_d^n) / (K^n + APC_d^n)</Expression>
                </DiffEqn>
                <DiffEqn symbol-ref="APC">
                    <Expression>α2*(1- APC) * ((CDK1_d^n) / (K^n + CDK1_d^n)) - β2*APC</Expression>
                </DiffEqn>
            </System>
        </CellType>
    </CellTypes>
    <CellPopulations>
        <Population size="1" type="cells"/>
    </CellPopulations>
    <Analysis>
        <Logger time-step="1e-2">
            <Restriction>
                <Celltype celltype="cells"/>
            </Restriction>
            <Input>
                <Symbol symbol-ref="APC"/>
                <Symbol symbol-ref="APC_d"/>
                <Symbol symbol-ref="CDK1"/>
                <Symbol symbol-ref="CDK1_d"/>
            </Input>
            <Output>
                <TextOutput/>
            </Output>
            <Plots>
                <Plot time-step="-1">
                    <Style style="lines" line-width="3.0"/>
                    <Terminal terminal="png"/>
                    <X-axis>
                        <Symbol symbol-ref="time"/>
                    </X-axis>
                    <Y-axis>
                        <Symbol symbol-ref="CDK1"/>
                        <Symbol symbol-ref="CDK1_d"/>
                        <Symbol symbol-ref="APC"/>
                        <Symbol symbol-ref="APC_d"/>
                    </Y-axis>
                </Plot>
            </Plots>
        </Logger>
    </Analysis>
</MorpheusModel>


In Morpheus GUI:
Examples → ODE → CellCycleDelay.xml.

Reference

Ferrell JE Jr, Tsai TY, Yang Q. Modeling the cell cycle: why do certain circuits oscillate? Cell, 18:144(6), 2011.



Coupled ODE lattice: Lateral signaling

Patterning as a result of lateral inhibition and lateral stabilization.

Introduction

This example model cell fate decisions during early patterning of the pancreas (de Back et al., 2012). The simple gene regulatory network of each cell is coupled to adjacent cells by lateral (juxtacrine) signaling.

Model description

The model defines a lattice of cells with a simplified hexagonal epithelial packing. This is specified in Space using a hexagonal lattice structure of size $(x,y,z)=(20,20,0)$ with periodic boundary conditions. The lattice is filled by seeding it with a Population of $400$ cells.

Each cell has two basic Properties X and Y representing the expression levels of Ngn3 and Ptf1a that are coupled in a System of DiffEqns.

The NeighborsReporter plugin is used to couple the cells to their directly adjacent neighbors. This plugin checks the values of X in neighboring cells and outputs its mean value in Property Xn.

This model uses a number of Analysis plugins:

  • Gnuplotter visualizes the values of Y with a ColorMap that maps values to colors. It outputs to screen (interactive mode) or to PNG (local mode).
  • Logger records the values of X and Y expression to file and, at the end of simulation, shows a time plot.
  • The first HistogramLogger records and plots the distribution of X and Y expression cells over time.
  • The second HistogramLogger records and, after simulation, plots the distribution of $\tau$, the time to cell fate decision (see reference).

Model

LateralSignaling.xml

<?xml version='1.0' encoding='UTF-8'?>
<MorpheusModel version="3">
    <Description>
        <Title>Example-LateralSignaling</Title>
        <Details>Reference:

Walter de Back, Joseph X. Zhou, Lutz Brusch, On the Role of Lateral Stabilization during Early Patterning in the Pancreas, Roy. Soc. Interface 10(79): 20120766, 2012.

http://dx.doi.org/10.1098/rsif.2012.0766
</Details>
    </Description>
    <Global>
        <Constant symbol="X" value="0"/>
        <Constant symbol="Y" value="0"/>
    </Global>
    <Space>
        <Lattice class="hexagonal">
            <Size symbol="size" value="20 20 0"/>
            <BoundaryConditions>
                <Condition boundary="x" type="periodic"/>
                <Condition boundary="y" type="periodic"/>
            </BoundaryConditions>
            <Neighborhood>
                <Order>1</Order>
            </Neighborhood>
        </Lattice>
        <SpaceSymbol symbol="space"/>
    </Space>
    <Time>
        <StartTime value="0"/>
        <StopTime value="30"/>
        <TimeSymbol symbol="t"/>
        <!--    <Disabled>
        <RandomSeed value="2"/>
    </Disabled>
-->
    </Time>
    <CellTypes>
        <CellType class="biological" name="cells">
            <Property symbol="X" value="0.0" name="Ngn3"/>
            <Property symbol="Xn" value="0.0" name="Ngn3-Neighbors"/>
            <Property symbol="Y" value="0" name="Ptf1a"/>
            <Property symbol="Yn" value="0" name="Ptf1a-neighbors"/>
            <System solver="heun" time-step="0.02">
                <Constant symbol="a" value="1"/>
                <Constant symbol="b" value="21"/>
                <Constant symbol="c" value="1"/>
                <DiffEqn symbol-ref="X">
                    <Expression>((th / (th + a*Xn^n)) - X) + rand_norm(0.0,noise)</Expression>
                </DiffEqn>
                <DiffEqn symbol-ref="Y">
                    <Expression>(((th + b*(Y * Yn)^n) / (th + c*X^n + b*(Y * Yn)^n))  - Y ) + rand_norm(0.0,noise)</Expression>
                </DiffEqn>
                <Constant symbol="n" value="4"/>
                <Constant symbol="th" value="1e-4"/>
                <Constant symbol="noise" value="1e-4"/>
            </System>
            <NeighborhoodReporter>
                <Input scaling="cell" value="X"/>
                <Output symbol-ref="Xn" mapping="average"/>
            </NeighborhoodReporter>
            <NeighborhoodReporter>
                <Input scaling="cell" value="Y"/>
                <Output symbol-ref="Yn" mapping="average"/>
            </NeighborhoodReporter>
            <Event trigger="on change">
                <Condition>tau == -1 and (X-Xn) > 0.05</Condition>
                <Rule symbol-ref="tau">
                    <Expression>t</Expression>
                </Rule>
            </Event>
            <Property symbol="tau" value="-1" name="time to cell fate decision"/>
        </CellType>
    </CellTypes>
    <CellPopulations>
        <Population size="0" type="cells">
            <InitCellLattice/>
        </Population>
    </CellPopulations>
    <Analysis>
        <Gnuplotter time-step="5">
            <Terminal size="400 400 0" persist="true" name="png"/>
            <Plot>
                <Cells value="X" min="0.0" max="1">
                    <ColorMap>
                        <Color value="1.0" color="blue"/>
                        <Color value="0.5" color="light-blue"/>
                        <Color value="0.0" color="white"/>
                    </ColorMap>
                </Cells>
            </Plot>
            <Plot>
                <Cells value="Y" min="0.0" max="1">
                    <ColorMap>
                        <Color value="1.0" color="red"/>
                        <Color value="0.5" color="light-red"/>
                        <Color value="0.0" color="white"/>
                    </ColorMap>
                </Cells>
            </Plot>
        </Gnuplotter>
        <Logger time-step="0.1">
            <Input>
                <Symbol symbol-ref="X"/>
                <Symbol symbol-ref="Y"/>
            </Input>
            <Output>
                <TextOutput file-separation="cell"/>
            </Output>
            <Plots>
                <Plot time-step="-1">
                    <Style style="lines" line-width="2"/>
                    <Terminal terminal="png"/>
                    <X-axis>
                        <Symbol symbol-ref="t"/>
                    </X-axis>
                    <Y-axis>
                        <Symbol symbol-ref="X"/>
                    </Y-axis>
                    <Color-bar>
                        <Symbol symbol-ref="Y"/>
                    </Color-bar>
                    <Range>
                        <Data increment="3"/>
                    </Range>
                </Plot>
            </Plots>
        </Logger>
        <HistogramLogger minimum="-0.1" normalized="true" maximum="1.1" time-step="5" number-of-bins="20">
            <Plot minimum="0" maximum="1.0" terminal="png"/>
            <Column symbol-ref="X" celltype="cells"/>
            <Column symbol-ref="Y" celltype="cells"/>
        </HistogramLogger>
        <HistogramLogger minimum="0.0" normalized="true" maximum="30" time-step="-1" number-of-bins="30">
            <Plot minimum="0" maximum="1.0" terminal="png"/>
            <Column symbol-ref="tau" celltype="cells"/>
        </HistogramLogger>
    </Analysis>
</MorpheusModel>


In Morpheus-GUI:
Examples → ODE → LateralSignaling.xml.

Things to try

  • Change the lattice structure from hexagonal to square. See Space/Lattice.
  • Change the strength of lateral stabilization b and observe the pattern. See CellTypes/CellType/System.
  • Change the noise amplitude and observe time to cell fate decision ($\tau$).

Reference

W de Back, J X Zhou, L Brusch, On the Role of Lateral Stabilization during Early Patterning in the Pancreas, Journal of the Royal Society Interface, 10:79, 2013.

examples/differential_equations.txt · Last modified: 16:35 02.12.2013 by Walter