Return to Question

The following system of coupled ODEs arises in the study of DNA sequence evolution:

\begin{eqnarray*} \frac{da}{dt} & = & \frac{\mu (1-y) b u}{S - y(S-b-v)} - (\lambda +\mu ) a \\ \frac{db}{dt} & = & -\frac{\mu (b+v) b}{S - y(S-b-v)} + \lambda (1-b) \\ \frac{du}{dt} & = & -\frac{\mu (b+v) u}{S - y(S-b-v)} + \lambda a \\ \frac{dv}{dt} & = & \frac{\mu (b+v) (S-v)}{S - y(S-b-v)} \\ S & = & \exp\left(\frac{\lambda t}{1-x}\right)-1 \end{eqnarray*} for $t>0$, with $a(0)=1$, $b(0)=u(0)=v(0)=0$, $a'(0)=-\lambda-\mu$, $b'(0)=u'(0)=\lambda$, and $v'(0)=0$.

Specifically, they are a moment-matched approximation to the General Geometric Indel model, the simplest continuous-time Markov chain on strings that allows insertions and deletions (a.k.a. "indels") of length > 1; details are here.

The parameters are $\lambda,\mu \in \Re^+$, the indel rates, and $x,y \in [0,1)$, the corresponding probability parameters for the geometric distributions over indel lengths. The mean insertion length is $1/(1-x)$ and the mean deletion length is $1/(1-y)$.

The special case $x=y=0$ (where indels all have unit length) is known as the TKF91 model, and has the following exact solution:

\begin{eqnarray*} a & = & \frac{\mu - \lambda}{\mu L - \lambda M} LM \\ b & = & \frac{\lambda}{\mu L - \lambda M} (L - M) \\ u & = & \frac{\mu - \lambda}{\mu L - \lambda M} M(1-L) \\ v & = & \frac{\mu - \lambda}{\mu L - \lambda M} - 1 \end{eqnarray*} with $L=\exp(-\lambda t)$ and $M=\exp(-\mu t)$.

Can anyone suggest a way to solve this? Things I have tried/observed:

1. Clearly $(b,v)$ form an independent subsystem of coupled ODEs, so it's presumably worth solving these first.
2. I haven't yet been able to simplify further to just a single ODE.
3. I've tried plugging into Mathematica's DSolve as-is. No joy.
4. Following the TKF91 solutions, I tried using a trial solution that is a rational function of polynomials that are first-order in $L=\exp(-\lambda t/(1-x))$ and $M=\exp(-\mu t/(1-y))$ and share a common denominator, i.e. $b=\frac{b_0 + b_L L + b_M M + b_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ and $v=\frac{v_0 + v_L L + v_M M + v_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ (Mathematica gets stuck; if I constrain $d_0=d_{LM}=0$, as in TKF91, it says there is no solution.)
5. I can of course approximate it numerically, e.g. with Runge Kutta methods, but I'd really like an analytical solution.

The following system of coupled ODEs arises in the study of DNA sequence evolution:

\begin{eqnarray*} \frac{da}{dt} & = & \frac{\mu (1-y) b u}{S - y(S-b-v)} - (\lambda +\mu ) a \\ \frac{db}{dt} & = & -\frac{\mu (b+v) b}{S - y(S-b-v)} + \lambda (1-b) \\ \frac{du}{dt} & = & -\frac{\mu (b+v) u}{S - y(S-b-v)} + \lambda a \\ \frac{dv}{dt} & = & \frac{\mu (b+v) (S-v)}{S - y(S-b-v)} \\ S & = & \exp\left(\frac{\lambda t}{1-x}\right)-1 \end{eqnarray*} for $t>0$, with $a(0)=1$, $b(0)=u(0)=v(0)=0$, $a'(0)=-\lambda-\mu$, $b'(0)=u'(0)=\lambda$, and $v'(0)=0$.

Specifically, they are a moment-matched approximation to the General Geometric Indel model, the simplest continuous-time Markov chain on strings that allows insertions and deletions (a.k.a. "indels") of length > 1; details are here.

The parameters are $\lambda,\mu \in \Re^+$, the indel rates, and $x,y \in [0,1)$, the corresponding probability parameters for the geometric distributions over indel lengths. The mean insertion length is $1/(1-x)$ and the mean deletion length is $1/(1-y)$.

The special case $x=y=0$ (where indels all have unit length) is known as the TKF91 model, and has the following exact solution:

\begin{eqnarray*} a & = & \frac{\mu - \lambda}{\mu L - \lambda M} LM \\ b & = & \frac{\lambda}{\mu L - \lambda M} (L - M) \\ u & = & \frac{\mu - \lambda}{\mu L - \lambda M} M(1-L) \\ v & = & \frac{\mu - \lambda}{\mu L - \lambda M} - 1 \end{eqnarray*} with $L=\exp(-\lambda t)$ and $M=\exp(-\mu t)$.

Can anyone suggest a way to solve this? Things I have tried/observed:

1. Clearly $(b,v)$ form an independent subsystem of coupled ODEs, so it's presumably worth solving these first.
2. I haven't yet been able to simplify further to just a single ODE.
3. I've tried plugging into Mathematica's DSolve as-is. No joy.
4. Following the TKF91 solutions, I tried using a trial solution that is a rational function of polynomials that are first-order in $L=\exp(-\lambda t/(1-x))$ and $M=\exp(-\mu t/(1-y))$ and share a common denominator, i.e. $b=\frac{b_0 + b_L L + b_M M + b_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ and $v=\frac{v_0 + v_L L + v_M M + v_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ (Mathematica gets stuck; if I constrain $d_0=d_{LM}=0$, as in TKF91, it says there is no solution.)
5. I can of course approximate it numerically, e.g. with Runge Kutta methods, but I'd really like an analytical solution.

No rush, I've been working on this for 20 years ;)

The following system of coupled ODEs arises in the study of DNA sequence evolution:

\begin{eqnarray*} \frac{da}{dt} & = & \frac{\mu (1-y) b u}{S - y(S-b-v)} - (\lambda +\mu ) a \\ \frac{db}{dt} & = & -\frac{\mu (b+v) b}{S - y(S-b-v)} + \lambda (1-b) \\ \frac{du}{dt} & = & -\frac{\mu (b+v) u}{S - y(S-b-v)} + \lambda a \\ \frac{dv}{dt} & = & \frac{\mu (b+v) (S-v)}{S - y(S-b-v)} \\ S & = & \exp\left(\frac{\lambda t}{1-x}\right)-1 \end{eqnarray*} for $t>0$, with $a(0)=1$, $b(0)=u(0)=v(0)=0$, $a'(0)=-\lambda-\mu$, $b'(0)=u'(0)=\lambda$, and $v'(0)=0$.

Specifically, they are a moment-matched approximation to the General Geometric Indel model, the simplest continuous-time Markov chain on strings that allows insertions and deletions (a.k.a. "indels") of length > 1; details are here.

The parameters are $\lambda,\mu \in \Re^+$, the indel rates, and $x,y \in [0,1)$, the corresponding probability parameters for the geometric distributions over indel lengths. The mean insertion length is $1/(1-x)$ and the mean deletion length is $1/(1-y)$.

The special case $x=y=0$ (where indels all have unit length) is known as the TKF91 model, and has the following exact solution:

\begin{eqnarray*} a & = & \frac{\mu - \lambda}{\mu L - \lambda M} LM \\ b & = & \frac{\lambda}{\mu L - \lambda M} (L - M) \\ u & = & \frac{\mu - \lambda}{\mu L - \lambda M} M(1-L) \\ v & = & \frac{\mu - \lambda}{\mu L - \lambda M} - 1 \end{eqnarray*} with $L=\exp(-\lambda t)$ and $M=\exp(-\mu t)$.

Can anyone suggest a way to solve this? Things I have tried/observed:

1. Clearly $(b,v)$ form an independent subsystem of coupled ODEs, so it's presumably worth solving these first.
2. I haven't yet been able to simplify further to just a single ODE.
3. I've tried plugging into Mathematica's DSolve as-is. No joy.
4. Following the TKF91 solutions, I tried using a trial solution that is a rational function of polynomials that are first-order in $L=\exp(-\lambda t/(1-x))$ and $M=\exp(-\mu t/(1-y))$ and share a common denominator, i.e. $b=\frac{b_0 + b_L L + b_M M + b_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ and $v=\frac{v_0 + v_L L + v_M M + v_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ (Mathematica gets stuck; if I constrain $d_0=d_{LM}=0$, as in TKF91, it says there is no solution.)
5. I can of course simulate it numerically, but I'd really like an analytical solution.

No rush, I've been working on this for 20 years ;)

The following system of coupled ODEs arises in the study of DNA sequence evolution:

\begin{eqnarray*} \frac{da}{dt} & = & \frac{\mu (1-y) b u}{S - y(S-b-v)} - (\lambda +\mu ) a \\ \frac{db}{dt} & = & -\frac{\mu (b+v) b}{S - y(S-b-v)} + \lambda (1-b) \\ \frac{du}{dt} & = & -\frac{\mu (b+v) u}{S - y(S-b-v)} + \lambda a \\ \frac{dv}{dt} & = & \frac{\mu (b+v) (S-v)}{S - y(S-b-v)} \\ S & = & \exp\left(\frac{\lambda t}{1-x}\right)-1 \end{eqnarray*} for $t>0$, with $a(0)=1$, $b(0)=u(0)=v(0)=0$, $a'(0)=-\lambda-\mu$, $b'(0)=u'(0)=\lambda$, and $v'(0)=0$.

Specifically, they are a moment-matched approximation to the General Geometric Indel model, the simplest continuous-time Markov chain on strings that allows insertions and deletions (a.k.a. "indels") of length > 1; details are here.

The parameters are $\lambda,\mu \in \Re^+$, the indel rates, and $x,y \in [0,1)$, the corresponding probability parameters for the geometric distributions over indel lengths. The mean insertion length is $1/(1-x)$ and the mean deletion length is $1/(1-y)$.

The special case $x=y=0$ (where indels all have unit length) is known as the TKF91 model, and has the following exact solution:

\begin{eqnarray*} a & = & \frac{\mu - \lambda}{\mu L - \lambda M} LM \\ b & = & \frac{\lambda}{\mu L - \lambda M} (L - M) \\ u & = & \frac{\mu - \lambda}{\mu L - \lambda M} M(1-L) \\ v & = & \frac{\mu - \lambda}{\mu L - \lambda M} - 1 \end{eqnarray*} with $L=\exp(-\lambda t)$ and $M=\exp(-\mu t)$.

Can anyone suggest a way to solve this? Things I have tried/observed:

1. Clearly $(b,v)$ form an independent subsystem of coupled ODEs, so it's presumably worth solving these first.
2. I haven't yet been able to simplify further to just a single ODE.
3. I've tried plugging into Mathematica's DSolve as-is. No joy.
4. Following the TKF91 solutions, I tried using a trial solution that is a rational function of polynomials that are first-order in $L=\exp(-\lambda t/(1-x))$ and $M=\exp(-\mu t/(1-y))$ and share a common denominator, i.e. $b=\frac{b_0 + b_L L + b_M M + b_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ and $v=\frac{v_0 + v_L L + v_M M + v_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ (Mathematica gets stuck; if I constrain $d_0=d_{LM}=0$, as in TKF91, it says there is no solution.)
5. I can of course approximate it numerically, e.g. with Runge Kutta methods, but I'd really like an analytical solution.

The following system of coupled ODEs arises in the study of DNA sequence evolution:

\begin{eqnarray*} \frac{da}{dt} & = & \frac{\mu (1-y) b u}{S - y(S-b-v)} - (\lambda +\mu ) a \\ \frac{db}{dt} & = & -\frac{\mu (b+v) b}{S - y(S-b-v)} + \lambda (1-b) \\ \frac{du}{dt} & = & -\frac{\mu (b+v) u}{S - y(S-b-v)} + \lambda a \\ \frac{dv}{dt} & = & \frac{\mu (b+v) (S-v)}{S - y(S-b-v)} \\ S & = & \exp\left(\frac{\lambda t}{1-x}\right)-1 \end{eqnarray*} for $t>0$, with $a(0)=1$, $b(0)=u(0)=v(0)=0$, $a'(0)=-\lambda-\mu$, $b'(0)=u'(0)=\lambda$, and $v'(0)=0$.

Specifically, they are a moment-matched approximation to the General Geometric Indel model, the simplest continuous-time Markov chain on strings that allows insertions and deletions (a.k.a. "indels") of length > 1; details are here.

The parameters are $\lambda,\mu \in \Re^+$, the indel rates, and $x,y \in [0,1)$, the corresponding probability parameters for the geometric distributions over indel lengths. The mean insertion length is $1/(1-x)$ and the mean deletion length is $1/(1-y)$.

The special case $x=y=0$ (where indels all have unit length) is known as the TKF91 model, and has the following exact solution:

\begin{eqnarray*} a & = & \frac{\mu - \lambda}{\mu L - \lambda M} LM \\ b & = & \frac{\lambda}{\mu L - \lambda M} (L - M) \\ u & = & \frac{\mu - \lambda}{\mu L - \lambda M} M(1-L) \\ v & = & \frac{\mu - \lambda}{\mu L - \lambda M} - 1 \end{eqnarray*} with $L=\exp(-\lambda t)$ and $M=\exp(-\mu t)$.

Can anyone suggest a way to solve this? Things I have tried/observed:

1. Clearly $(b,v)$ form an independent subsystem of coupled ODEs, so it's presumably worth solving these first.
2. I haven't yet been able to simplify further to just a single ODE.
3. I've tried plugging into Mathematica's DSolve as-is. No joy.
4. Following the TKF91 solutions, I tried using a trial solution that is a rational function of polynomials that are first-order in $L=\exp(-\lambda t/(1-x))$ and $M=\exp(-\mu t/(1-y))$ and share a common denominator, i.e. $b=\frac{b_0 + b_L L + b_M M + b_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ and $v=\frac{v_0 + v_L L + v_M M + v_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ (Mathematica gets stuck; if I constrain $d_0=d_{LM}=0$, as in TKF91, it says there is no solution.)

No rush, I've been working on this for 20 years ;)

The following system of coupled ODEs arises in the study of DNA sequence evolution:

\begin{eqnarray*} \frac{da}{dt} & = & \frac{\mu (1-y) b u}{S - y(S-b-v)} - (\lambda +\mu ) a \\ \frac{db}{dt} & = & -\frac{\mu (b+v) b}{S - y(S-b-v)} + \lambda (1-b) \\ \frac{du}{dt} & = & -\frac{\mu (b+v) u}{S - y(S-b-v)} + \lambda a \\ \frac{dv}{dt} & = & \frac{\mu (b+v) (S-v)}{S - y(S-b-v)} \\ S & = & \exp\left(\frac{\lambda t}{1-x}\right)-1 \end{eqnarray*} for $t>0$, with $a(0)=1$, $b(0)=u(0)=v(0)=0$, $a'(0)=-\lambda-\mu$, $b'(0)=u'(0)=\lambda$, and $v'(0)=0$.

Specifically, they are a moment-matched approximation to the General Geometric Indel model, the simplest continuous-time Markov chain on strings that allows insertions and deletions (a.k.a. "indels") of length > 1; details are here.

The parameters are $\lambda,\mu \in \Re^+$, the indel rates, and $x,y \in [0,1)$, the corresponding probability parameters for the geometric distributions over indel lengths. The mean insertion length is $1/(1-x)$ and the mean deletion length is $1/(1-y)$.

The special case $x=y=0$ (where indels all have unit length) is known as the TKF91 model, and has the following exact solution:

\begin{eqnarray*} a & = & \frac{\mu - \lambda}{\mu L - \lambda M} LM \\ b & = & \frac{\lambda}{\mu L - \lambda M} (L - M) \\ u & = & \frac{\mu - \lambda}{\mu L - \lambda M} M(1-L) \\ v & = & \frac{\mu - \lambda}{\mu L - \lambda M} - 1 \end{eqnarray*} with $L=\exp(-\lambda t)$ and $M=\exp(-\mu t)$.

Can anyone suggest a way to solve this? Things I have tried/observed:

1. Clearly $(b,v)$ form an independent subsystem of coupled ODEs, so it's presumably worth solving these first.
2. I haven't yet been able to simplify further to just a single ODE.
3. I've tried plugging into Mathematica's DSolve as-is. No joy.
4. Following the TKF91 solutions, I tried using a trial solution that is a rational function of polynomials that are first-order in $L=\exp(-\lambda t/(1-x))$ and $M=\exp(-\mu t/(1-y))$ and share a common denominator, i.e. $b=\frac{b_0 + b_L L + b_M M + b_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ and $v=\frac{v_0 + v_L L + v_M M + v_{LM} LM}{d_0 + d_L L + d_M M + d_{LM} LM}$ (Mathematica gets stuck; if I constrain $d_0=d_{LM}=0$, as in TKF91, it says there is no solution.)
5. I can of course simulate it numerically, but I'd really like an analytical solution.

No rush, I've been working on this for 20 years ;)