2
$\begingroup$

Context for this question: I'm in mathematical biology, and while I've got a decent grasp of the concepts, and am a decent hand with simulation and numerical approaches, my analytical abilities and general skill at manipulating equations has always been rough. Probably the root of deciding I was "bad at math" back in the day, though that's neither here nor there.

What I'm trying to do is find the equilibriums for a set of differential equations, similar to the following:

$$\begin{align} \frac{{dS}}{{dt\strut}} =& \mu - (\beta I + \mu )S\\ \frac{{dE}}{{dt\strut}} =& \beta SI - (\mu + \sigma )E\\ \frac{{dI}}{{dt\strut}} =& \sigma E - (\mu + \gamma )I\\ \frac{{dR}}{{dt\strut}} =& \gamma I - \mu R \end{align}$$

I can do this one, though admittedly it takes me quite some time, but more complex systems involve a fair number of hours going down dead ends that really shouldn't be, basic errors etc.

Any chance any of the myriad math software packages out there have a way to simplify the computation/manipulation steps of finding the equilibria for a system like this? Any general tips if there aren't? You've got Matlab, Mathematica, and really any open source package you please to work with.

Improve this question
$\endgroup$
0
3
$\begingroup$

Since at the equilibrium time derivatives are zero, the issue reduces to solving a system of algebraic equations. Mathematica has many commands for equation solving, among them Solve and Reduce.

Here is a snapshot of solving your particular system:

enter image description here

Improve this answer
$\endgroup$
3
  • $\begingroup$ This is exactly what I was talking about, as even solving a system of algebraic equations results in errors because of how genuinely terrible I am at solving anything. $\endgroup$ – Fomite Nov 4 '11 at 1:22
  • $\begingroup$ I understand most of the code above, but could you possibly walk me through the arguments in the Simplify part? Particularly the beta*sigma =/= (gamma+mu) (sigma+mu) bit. $\endgroup$ – Fomite Nov 4 '11 at 1:44
  • $\begingroup$ @EpiGrad The simplify was there to remove additional requirements identified by Solve, which were $\beta \not= \frac{(\gamma+\mu)(\mu+\sigma)}{\sigma}$. $\endgroup$ – Sasha Nov 4 '11 at 1:45
3
$\begingroup$

In Maple, it's pretty straightforward.

eqs:= {mu - (beta*i + mu)*S = 0, beta*S*i - (mu+sigma)*E = 0, sigma*E - (mu+gamma)*i = 0, gamma*i - mu*R = 0}; solve(eqs,{S,E,i,R});

{E = 0, R = 0, S = 1, i = 0}, {E = -mu*(-beta*sigma+mu^2+mu*gamma+sigma*mu+sigma*gamma)/sigma/(mu+sigma)/beta, R = -gamma*(-beta*sigma+mu^2+mu*gamma+sigma*mu+sigma*gamma)/beta/(mu^2+mu*gamma+sigma*mu+sigma*gamma), S = (mu^2+mu*gamma+sigma*mu+sigma*gamma)/beta/sigma, i = -mu*(-beta*sigma+mu^2+mu*gamma+sigma*mu+sigma*gamma)/beta/(mu^2+mu*gamma+sigma*mu+sigma*gamma)}

Improve this answer
$\endgroup$
3
$\begingroup$

To find Equilibrium or Steady State of Differential Equations using MATLAB would be,

    clc, clear all;
    syms S E I R mu beta sigma gamma 'real'
    e1 = eval('mu-(beta*I+mu)*S');
    e2 = eval('beta*S*I-(mu+sigma)*E');
    e3 = eval('sigma*E-(mu+gamma)*I');
    e4 = eval('gamma*I-mu*R');
    sol=solve(e1,S,e2,E,e3,I,e4,R);
    disp([sol.S sol.E sol.I sol.R])

Output;

$$ S=(\gamma\mu + \gamma\sigma + \mu\sigma + \mu^2)/(\beta\sigma),\ \\ E=-(\gamma\mu^2 + \mu^2\sigma +\mu^3 - \beta\mu\sigma + \gamma\mu\sigma)/(\beta\sigma^2 + \beta\mu\sigma),\ \\ I=-(\gamma\mu^2 + \mu^2\sigma + \mu^3 - \beta\mu\sigma + \gamma\mu\sigma)/(\beta\mu^2 + \beta\gamma\mu + \beta\gamma\sigma + \beta\mu\sigma),\ \\ R=-(\gamma^2\mu + \sigma\gamma^2 + \gamma*\mu^2 + \sigma\gamma\mu - \beta\sigma\gamma)/(\beta\mu^2 + \beta\gamma\mu + \beta\gamma\sigma + \beta\mu\sigma) $$

$\endgroup$

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.