# LAPLACE TRANSFORMS FOR REACTION NETWORKS

# Preliminaries


```python
!pip -q install controlSBML
import controlSBML as ctl
import control
from controlSBML.util import makeSimulationTimes
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
import sympy as sy
import tellurium as te
import sympy
print("controlSBML version: " + ctl.__version__)
```

    controlSBML version: 0.2.8



```python
s, S1, S2, S3, k1, k2, k3, k4 = sympy.symbols("s, S1, S2, S3, k1, k2, k3, k4 ")
```

# Sequential Pathway

$S_1 \xrightarrow[\text{}]{k_1 S_1} S_2$,
$S_2 \xrightarrow[\text{}]{k_2 S_2} S_3$

Find the transfer function for $\frac{S_3}{S_1}$

## Steps
1. Specify the input and output for the transfer function.
1. Write the system equations.
1. Convert to Laplace transforms.
1. Solve for the output LT in terms of the input LT and form the ratio.

### Specify the input and output for the transfer function

This is specified in the problem. The input is $S_1$ and the output is $S_3$.

### Write the system equations

The reactions are:

    S1 -> S2; k1*S1
    S2 -> S3; k2*S2

From inspections, the system equations are:

\begin{eqnarray}
\dot{S}_1 & = & - k_1 S_1 \\
\dot{S}_2 & = & k_1 S_1 - k_2 S_2
\end{eqnarray}

### Convert the system equations to Laplace transforms

This is a mechanical process based on the fact that $L[\dot{g}(t)] = s G(s)$,
where $L[g(t)] = G(s)$.

\begin{eqnarray}
s S_1(s) & = & - k_1 S_1(s) \\
sS_2 (s) & = & k_1 S_1(s) - k_2 S_2(s)
\end{eqnarray}

### Solve for LT of the output in terms of the input and form the ratio

## Simulate to check

## Analysis

### Poles

### Response to initial conditions

### DC Gain (Step Response)

# Branching

$S_1 \xrightarrow[\text{}]{k_1 S_1} S_2$,
$S_2 \xrightarrow[\text{}]{k_2 S_2} S_3$,
$S_2 \xrightarrow[\text{}]{k_3 S_2} S_4$


Find the transfer function for $\frac{S_3}{S_1}$

# Reversible Reaction

$S_1 \xrightarrow[\text{}]{k_1 S_1} S_2$,
$S_2 \xrightarrow[\text{}]{k_2 S_2} S_1$


Find the transfer function for $\frac{S_2}{S_1}$

### Write the system equations

\begin{eqnarray}
\dot{S_1} & = & -k_1 S_1 + k_2 S_2 \\
\dot{S_2} & = & k_1 S_1 - k_2 S_2 \\
\end{eqnarray}

### Find the transfer function


```python
expr = -k1*S1 + k2*S2 - s*S1
tf = sympy.solve(expr, S2)[0] / S1
tf
```




$\displaystyle \frac{k_{1} + s}{k_{2}}$



What is $\frac{k_1}{k_2}$ in terms of first year chemistry?

# Nonlinear

$S_1 \xrightarrow[\text{}]{k_1 S_1 } S_2$,
$S_2 \xrightarrow[\text{}]{k_2 S_1 S_2} S_3$

Find the transfer function for $\frac{S_3}{S_1}$

# Sequential Pathway With a Reversible Reaction

$S_1 \xrightarrow[\text{}]{k_1 S_1} S_2$,
$S_2 \xrightarrow[\text{}]{k_2 S_2} S_3$,
$S_3 \xrightarrow[\text{}]{k_3 S_3} S_2$,
$S_3 \xrightarrow[\text{}]{k_4 S_3} \emptyset$

Find the transfer function for $\frac{S_3}{S_1}$

### Write the state equations

\begin{eqnarray}
\dot{S_1} & = & -k_1 S_1 \\
\dot{S_2} & = & k_1 S_1 - k_2 S_2 + k_3 S_3 \\
\dot{S_3} & = & k_2 S_2 - k_3 S_3 \\
\end{eqnarray}

### Convert the system equations to Laplace Transforms

\begin{eqnarray}
s S_1 (s) & = & -k_1 S_1 (s) \\
s S_2 (s)  & = & k_1 S_1(s) - k_2 S_2(s) + k_3 S_3(s) \\
s S_3(s)  & = & k_2 S_2 (s) - k_3 S_3(s) -k_4 S_3(s) \\
\end{eqnarray}

### Express the output in terms of the input

``sympy`` has its limits. It works best if there is one expression with the two variables of interest.


```python
S2_expr = k1*S1 - k2*S2 + k3*S3 - s*S2
S2_subs = sympy.solve(S2_expr, S2)[0]
S2_subs
```




$\displaystyle \frac{S_{1} k_{1} + S_{3} k_{3}}{k_{2} + s}$




```python
S3_expr = k2*S2 - k3*S3 - k4*S3 - s*S3
new_expr = S3_expr.subs({S2: S2_subs})
result = sympy.solve(new_expr, S3)[0]
result
```




$\displaystyle \frac{S_{1} k_{1} k_{2}}{k_{2} k_{4} + k_{2} s + k_{3} s + k_{4} s + s^{2}}$




```python
tf = result/S1
tf
```




$\displaystyle \frac{k_{1} k_{2}}{k_{2} k_{4} + k_{2} s + k_{3} s + k_{4} s + s^{2}}$




```python
numerator, denominator = sympy.fraction(tf)
numerator, denominator
```




    (k1*k2, k2*k4 + k2*s + k3*s + k4*s + s**2)




```python
# We can find the poles symbolically by solving for s
sim_dct = {k1: 1, k2: 2, k3: 2.5, k4: 1}
sympy.solve(denominator.subs(sim_dct), s)
```




    [-5.10849528301415, -0.391504716985849]




```python
dct = dict(sim_dct)
dct.update({s: 0})
tf.subs(dct)
```




$\displaystyle 1$



What constants can we change to get a DC gain of 4?

## Simulation


```python
MODEL = """
J1: S1 -> S2; k1*S1
J2: S2 -> S3; k2*S2
J3: S3 -> S2; k3*S3
J4: S3 -> ; k4*S3
k1 = 1
k2 = 2
k3 = 2.5
k4 = 1
S1 = 1
S2 = 0
S3 = 0
"""
RR = te.loada(MODEL)
RR.simulate(0, 20)
RR.plot(title="Impulse Response")
```


```python
STEP_MODEL = """
J1: $S1 -> S2; k1*$S1
J2: S2 -> S3; k2*S2
J3: S3 -> S2; k3*S3
J4: S3 -> ; k4*S3
k1 = 1
k2 = 2
k3 = 2.5
k4 = 1
$S1 = 1
S2 = 0
S3 = 0
"""
rr = te.loada(STEP_MODEL)
rr.simulate(0, 20)
rr.plot(title="Step Response")
```


```python
ctlsb = ctl.ControlSBML(MODEL, input_names=["S1"], output_names=["S3"])
tf = ctlsb.makeTransferFunction()
```


```python
tf.dcgain()
```




    0.9999999999999996




```python
tf.pole()
```




    array([-5.10849528, -1.        , -0.39150472])



We got one more pole than solving symbolically beause the symbolic solution simplified the system by eliminating $S_2$.


```python
state_space = ctlsb.makeStateSpace()
times = [0.01*v for v in range(5001)]
result = control.forced_response(tf, times, np.repeat(10, len(times)))
```


```python
tf
```




$$\frac{2}{s^3 + 6.5 s^2 + 7.5 s + 2}$$




```python
plt.plot(result.t, result.outputs)
```


```python
result.y
```




    array([[0.00000000e+00, 3.27974089e-06, 2.58182195e-05, ...,
            9.99999994e+00, 9.99999994e+00, 9.99999994e+00]])




```python
type(tf)
```




    control.xferfcn.TransferFunction



9bbd583ed753990ce61ac7a05e1a5c7cbdd4f807
ipynb
Jupyter Notebook
Lecture_10_11-Practice-Constructing-Laplace-Transforms-For-Reaction-Networks/Practice-Constructing-Laplace-Transforms-For-Reaction-Networks.ipynb
joseph-hellerstein/advanced-controls-lectures
dc43f6c3517616da3b0ea7c93192d911414ee202
MIT
Lecture_10_11-Practice-Constructing-Laplace-Transforms-For-Reaction-Networks/Practice-Constructing-Laplace-Transforms-For-Reaction-Networks.ipynb
joseph-hellerstein/advanced-controls-lectures
dc43f6c3517616da3b0ea7c93192d911414ee202
MIT
Lecture_10_11-Practice-Constructing-Laplace-Transforms-For-Reaction-Networks/Practice-Constructing-Laplace-Transforms-For-Reaction-Networks.ipynb
joseph-hellerstein/advanced-controls-lectures
dc43f6c3517616da3b0ea7c93192d911414ee202
MIT
Qwen/Qwen-72B
1. YES
2. YES