# How to solve forward Kolmogorov birth-death equations for unspecified number of populations

I'm attempting to use a forward Kolmogorov differential equation to model a birth-death process. This is fairly trivial when there's only one population, but I'm working with an unspecified and time-varying number of populations, so simply developing a system of equations doesn't seem viable.

I'm trying to solve for the dynamics equivalent to this paper by Desai and Fisher (genetics.org/content/genetics/176/3/1759.full.pdf), but using a slightly different approach. More importantly, I'm trying to extend their work, but can't seem to get the basic formatting down for solving their original system.

Supposing a set of populations each undergoing a separate birth-death process (with mutations feeding in from less fit populations to more fit ones) with fitness denoted as f, I have some set of population couts n(f,t) and probabilities of the a population being at that number at time p(n(f,t)).

My first and second attempts at solving this in Mathematica are below:

DSolve[Prob[n[f,t]]'==\[Mu]*n[f-1,t]*Prob[n[f-1,t]] + (f-Mean[f] - \[Mu])*n[f,t]*Prob[n[f,t]], Prob[n[f,t]],t]


And:

DSolve[Prob[n,f,t]'==\[Mu]*n*Prob[n,f-1,t] + (f-Mean[f] - \[Mu])*n*Prob[n,f,t], Prob[n,f,t],t]


In words, the equations depict what happens as individuals from a population with fitness f-1 feed into the population of fitness f by acquiring mutations at rate Mu. Meanwhile, individuals from population of fitness f are growing or dying off at rate proportional to their fitness relative to the mean fitness, and are also being drained into the population of fitness f+1 at rate Mu.

I'm not sure how to get DSolve to solve for the general dynamics of all populations. Any assistance or advice would be greatly appreciated!

• I'm interested in this problem but won't have time to read the paper for a few days. From your question, I'm a bit confused by the equations. If you have multiple populations, won't each need its own n? So, solve for Prob[n1, n2, n3, ...], which is out of control for more than a few n's. You might end up needing NDSolve, or do eigenanalysis of the transition matrix or just simulate the stochastic process. Anyhow, I hope to look at this more in a week or so. – Chris K Oct 2 at 6:07
• In the nonlinear FEM tests for NDSolve there is an example with analytical solution of a Fisher equation. See this link FEMDocumentation/tutorial/NonlinearFiniteElementVerificationTests# 703651640. Maybe that is useful. – user21 Oct 2 at 7:52