I am having hard time to convert following set of differential equation to state space equation. I am a biologist and my math skills fall short as I don't know where to start. Any suggestion or feedback is highly appreciated. Thanks in advance.

Updates: Based on comments I have updated the question (which is not a good practice, my sincere apologies to everyone),
I am modelling biological system with a cascade of signal transduction steps with feedback loops. Each level in the cascade has corresponding equation $\dot{y_i}$ or $dy_i/dt$ given by \begin{equation*} \frac{dy_i}{dt}= g_i(v_i, y_i)+ c_i + e_i \end{equation*} For example \begin{equation*} \frac{dy_0}{dt}= g_0(a_{00}y_0+a_{02}y_2, y_0)+ c_0 \end{equation*} \begin{equation*} \frac{dy_1}{dt}= g_1(a_{10}y_0+a_{11}y_1+a_{12}y_2, y_1)\\ \end{equation*} \begin{equation*} \frac{dy_2}{dt}= g_2(a_{20}y_0+a_{21}y_1+a_{22}y_2+a_{23}y_3+a_{24}y_4+a_{25}y_5, y_2)+e_2 \end{equation*} \begin{equation*} \frac{dy_3}{dt}= g_3(a_{32}y_2+a_{33}y_3, y_3) \end{equation*} \begin{equation*} \frac{dy_4}{dt}= g_4(a_{42}y_2+a_{44}y_4, y_4) \end{equation*} \begin{equation*} \frac{dy_5}{dt}= g_5(a_{25}x_2-a_{25}x_5) \end{equation*}

and so on. $y_i$ represents state variable and only $y_5$ can be observed ("output") others are hidden. $e_2$ is input variable. $a_{jk}$ represents parameters. Function $g_i$ is given by multiplication of $h_i$ and $r_i$ \begin{equation*} g_i({v_i}, y_i)= h_i ({v_i})\cdot r_i({v_i}, y_i)\\ \end{equation*} where function $h_i$ and $r_i$ are given by
\begin{equation*} h_i({v_i}) = \begin{cases} \frac{{v_i}}{1+\frac{{v_i}}{S_i}(1-exp(-{v_i}/S_i))} & \mbox{when } {v_i}> 0, \cr {v_i} & \mbox{when } {v_i}\leq 0, \end{cases} \end{equation*} \begin{equation*} r_i({v_i}, y) = \begin{cases} 1-exp(\frac{{v_i}^2S_i}{{v_i}(\varepsilon -y)^2}) & \mbox{when } y<\varepsilon\ \&\ {v_i}< 0, \cr 1 & \mbox{otherwise.} \end{cases} \end{equation*}


Function $h_i$ puts physiologically relevant soft upper limit.

Function $r_i$ ensures non-negative ligand concentration.

$\dot v$ is corresponding linear rate.

$S_i$ is max size of the pool.

$\varepsilon$ is a small positive constant.

Update-2: Ok, I want to know if I am doing right thing here. I need state space equation for my system of differential equation with nonlinear function. Technically I can have two different versions

A non-linear version \begin{equation*} \frac{dy_i}{dt}= g_i(v_i, y_i)+ c_i + e_i \end{equation*}

A linear version \begin{equation*} \frac{dy_i}{dt}= v_i+ c_i + e_i \end{equation*}

For linear version I can write the state equation as $y'= Ay+c+e$ (see the matrix equation as Image)

Taking this to next step can I write like following?

$y'= g(Ay)+c+e$

  • $\begingroup$ I've taken the liberty to modify some of the TeX to improve readability. $\endgroup$ – José Figueroa-O'Farrill Sep 9 '10 at 3:06
  • $\begingroup$ I think that the notation leaves a lot to be desired. You seem to have $g_0$ and $g_1$ and this is what presumably you are calling generically $g_i$, which is then given in terms of $h_i$ and $r_i$. (I'm guessing that the $\cdot$ is just multiplication?) However, there is no $i$ in the definition of $h_i$ nor in that of $r_i$. ??? $\endgroup$ – José Figueroa-O'Farrill Sep 9 '10 at 3:09
  • $\begingroup$ Also, what is $\frac{dx_i}{dt}$? Is this just $y_i$? $\endgroup$ – José Figueroa-O'Farrill Sep 9 '10 at 3:11
  • $\begingroup$ I'm having a difficult time following your notation and your equations. Perhaps you could add a few sentences at the beginning explaining the biological system which you are modeling along with a set of definitions of the variables. This would help me understand your formula and what it is that you are trying to accomplish. $\endgroup$ – sleepless in beantown Sep 9 '10 at 3:24
  • $\begingroup$ Abht, is there anyone with you who can help you in person? $\endgroup$ – Will Jagy Sep 9 '10 at 3:25

Ok, let me give this a shot:

  • Because your system is nonlinear, I'm assuming you want the nonlinear state-space form. You can easily get the linear form by doing a Taylor series expansion on it around some equilibrium point.
  • The fact that most of the states are not measurable is not a big problem. You can estimate them using your output variables (subject to observability conditions), using an state observer such as a Kalman filter, Moving Horizon Estimator (MHE) or a Luenberger observer. Also note for a nonlinear system, only local observability can be checked.
  • Because you have conditional statements, I don't believe you'll be able to write the above as a single state-space system. You have 4 conditions, but they can be reduced to 3 disjunctions, so you'll need 3 state-space systems and a conditional switching equation that "activates" the correct state-space system depending on the values of $v$ and $y$. This is known as a hybrid (or switched) system.
  • As to the treatment of a hybrid system, perhaps you could clarify what the purpose is of getting your model into state-space form. Is it for simulation reasons? Do you need it in order to do analysis (i.e. controllability, observability)? Or do you need to do optimization? If it is the last case, you can write the logic as a disjunctive program, which will allow you generate a very efficient mixed-integer programming (MIP) problem.

Anyway, this is one way of writing your state-space system:

$$ \begin{align} \frac{dy_{i}}{dt} &= g_{i}^{m}(v_{i},y_{i}) + p_{i} c_{i} + q_{i} e_{i},\quad i=0,\ldots,N-1\\ v_{i} &= \sum_{j=0}^{N-1} a_{ij}y_{j},\quad i=0,\ldots,N-1 \end{align} $$ where $p_{i},q_{i} \in \{0,1\}$ = coefficients, $N$ = number of states, and $m \in \{1,2,3\}$ = modes of the system. In addition, you will need a switching function $T(m,v_{i},y_{i}) = 0$ to select the appropriate modes based on the current states of your system. This can be done programmatically through IF-THEN-ELSE clauses (or via integer variables in an optimization problem).

  • For $m = 1$ (where $v_{i} < 0, y_{i} < \varepsilon$):

$$g_{i}^{1}(v_{i},y_{i}) = v_{i} \left[1-\exp\left(\frac{{v_i}^2S_i}{{v_i}(\varepsilon -y_{i})^2}\right)\right]$$

  • For $m = 2$ (where $v_{i} \leq 0, y_{i} \geq \varepsilon$):

$$g_{i}^{2}(v_{i},y_{i}) = v_{i}$$

  • For $m = 3$ (where $v_{i} > 0, y_{i} \in \mathbb{R}$):

$$g_{i}^{3}(v_{i},y_{i}) = \frac{{v_i}}{1+\frac{{v_i}}{S_i}(1-\exp(-{v_i}/S_i))}$$

  • $\begingroup$ Thanks very much Gilead, I am going through your answer. Seems like problem can be solved but require some advanced tricks (my mathematics knowledge is in fancy stage). My purpose is dual estimation, where I can estimate states as well as parameters using a nonlinear variant of Kalman filter. I have discrete time series data for $y_5$. $\endgroup$ – Abhishek Tiwari Sep 9 '10 at 22:39
  • $\begingroup$ I haven't really done any work on the parameter and state estimation of hybrid systems myself, but it seems to me your Kalman filter has to modified significantly (and its basic assumptions checked). I'm almost thinking it's easier to formulate the estimation problem as an MHE problem, where your states and parameters are obtained through optimization (see www-rocq.inria.fr/who/Giancarlo.Ferrari-Trecate/files/…). Perhaps you could get help from someone? $\endgroup$ – Gilead Sep 10 '10 at 1:13
  • $\begingroup$ Thanks Gilead, at this point I am not sure what is best way to go but I think I can accept you answer as it helped me a lot. $\endgroup$ – Abhishek Tiwari Sep 10 '10 at 3:24
  • $\begingroup$ @gilead, informative response, thanks. The hybrid system with the disjunctions is what I could not wrap my head around. @Abht, thanks for clarifying your question. It was of interest to me to see the progression and to see gilead's method of attack. I'll have to look up more about MHE. $\endgroup$ – sleepless in beantown Sep 10 '10 at 3:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.