Estimating Population Growth with Limited Information

Question

I have been thinking about these problems for a while and think I might have found a way to partly answer them. These problems deal with estimating the birth and death rate of a system in different scenarios.

Scenario 1:

• Suppose at time=0 my population is $k_1$. At some final time =t, my population is $k_2$ ( $k_1, k_2$ are positive integers)
• At each time point from (0,T), I assume that there is a probability of $p_1$ that the current population can increase by $j$ units and a probability of $p_2$ that the current population decrease by $l$ units (assume $p_1, p_2, j, l$ are all constant)
• I only observe the system at time=0 and time=t (i.e. I know the exact values of $k_1, k_2$.
• Assuming that the population can never go below 0, what are the most likely values of $p_1, p_2, j, l$ ?

Scenario 2:

• Suppose at time=0 my population is $k_1$. At some final time =t, my population is $k_2$ ( $k_1, k_2$ are positive integers)
• At each time point from (0,T), I assume that there is a probability of $p_1$ that the current population can increase by $j$ units and a probability of $p_2$ that the current population decrease by $l$ units (assume $p_1, p_2, j, l$ are all constant)
• I observe the population of time =0 and time=T ... as well as at some times $t_i$ (e.g. suppose I have intervals 0,1,2,3,4,5,6,7,8,9,10. I observe the population at times=0, 5,6,9,10 ... at time=0 population is $k_1$, at time=5 population is $a$, time= 6 population is $b$, time = 9 population is $c$, time=10 population is $k_2$)
• Assuming that the population can never go below 0, what are the most likely values of $p_1, p_2, j, l$ ?

My Answers:

For Scenario 1, I used a latent variable approach (e.g. EM algorithm/Gaussian Mixture) wrote:

$$ L(p_1, p_2, j, l) = \sum_{j=0}^{J} \sum_{l=0}^{L} \sum_{n_1=0}^{T} \sum_{n_2=0}^{T} \left[ p_1^{n_1} \cdot (1 - p_1)^{T - n_1} \cdot p_2^{n_2} \cdot (1 - p_2)^{T - n_2} \right] \cdot I(n_1 \cdot j - n_2 \cdot l = k_2 - k_1) $$

• $p_1^{n_1}$ represents the probability of $n_1$ increases in the population. $p_2^{n_2}$ represents the probability of $n_2$ decreases in the population.
• $(1 - p_1 - p_2)^{T - n_1 - n_2}$ represents the probability of the remaining time intervals during which the population neither increases nor decreases. This happens with probability $1 - p_1 - p_2$, and there are $T - n_1 - n_2$ such intervals.
• $I(n_1 \cdot j - n_2 \cdot l = k_2 - k_1)$ is an indicator function that takes values 1 if the condition is true else 0. This is to ensure that the total number of decreases and increases respect the initial and final population
• However, I am not sure if this Likelihood Function prevents the population from going below 0 and some intermediate time point. I am also not sure if a combinatorial/multinomial term is needed in the likelihood

I think its more difficult to write the likelihood function for Scenario 2, even though Scenario 2 has more information compared to Scenario 1.

I think this question can be broken down into two parts: A likelihood function for the times where we have full information, and a likelihood function for the times where we have missing information:

$$ L(p_1, p_2, j, l) = L_{obs}(p_1, p_2, j, l) \cdot L_{unobs}(p_1, p_2, j, l) $$

From this point on, I think we can use a similar approach as in Scenario 1:

$$ L_{obs}(p_1, p_2, j, l) = \prod_{i=0}^{N_{obs}-1} \sum_{n_1=0}^{t_{o_{i+1}} - t_{o_i}} \sum_{n_2=0}^{t_{o_{i+1}} - t_{o_i}} \sum_{j=0}^{J} \sum_{l=0}^{L} \left[ p_1^{n_1} \cdot (1 - p_1)^{t_{o_{i+1}} - t_{o_i - n_1}} \cdot p_2^{n_2} \cdot (1 - p_2)^{t_{o_{i+1}} - t_{o_i - n_2}} \right] \cdot I(n_1 \cdot j - n_2 \cdot l = k_{t_{o_{i+1}}} - k_{t_{o_i}}) $$

$$ L_{unobs}(p_1, p_2, j, l) = \prod_{i=0}^{N_{unobs}-1} \sum_{n_1=0}^{t_{u_{i+1}} - t_{u_i}} \sum_{n_2=0}^{t_{u_{i+1}} - t_{u_i}} \sum_{j=0}^{J} \sum_{l=0}^{L} \left[ p_1^{n_1} \cdot (1 - p_1)^{t_{u_{i+1}} - t_{u_i - n_1}} \cdot p_2^{n_2} \cdot (1 - p_2)^{t_{u_{i+1}} - t_{u_i - n_2}} \right] \cdot I(n_1 \cdot j - n_2 \cdot l = k_{t_{u_{i+1}}} - k_{t_{u_i}}) $$

I have been thinking about how to correctly write the likelihood functions for both Scenario 1 and Scenario 2 for quite some time and find myself getting lost/confused. Can someone please help me write these correctly?

Answer 1

In both scenario's I'm assuming (based on the likelihood functions you've provided) that $p_1$ and $p_2$ represent the probabilities of independent events (at each time) which increase the population by $j$ or decrease it by $l$ respectively. Therefore, it is possible for there to be times where the population increases by $j-l$ (or decreases by $l-j$ if $l > j$). I will also assume your first observation is at time $0$ and your last one is at time $T$. For convenience, let $X_t$ denote the population at time $t$.

Scenario 1: The likelihood is actually quite different than you think. However, you are right. You would need to adjust your likelihood function to avoid the case of a negative population and you do need an extra combinatorial term. I could give you the likelihood function, but I think it may be helpful for you to see how you can play around with your data and your model to construct the correct likelihood function piece by piece. This will not be the most efficient answer, but I hope it's clear!

When you construct a likelihood function, you want to consider two things. What is your data? What are your parameters? In this case, your data is $X_0=k_1$ and $X_T = k_2$ and your parameters are $p_1,p_2,j$ and $l$. So your likelihood function is of the form: $$L_{0,T}(p_1,p_2,j,l;k_1,k_2) = P_{p_1,p_2,j,l}(X_0=k_1,X_T=k_2),$$ where I use $L_{0,T}$ to indicate the likelihood function given the population at time $0$ and time $T$ (this will be useful for scenario 2) and $P_{p_1,p_2,j,l,k_1}$ is the probability given the relevant parameters, assuming that $X_0 = k_1$.

$n_1$ and $n_2$ are useful latent variables, but you need to be careful in how you use them. Fix any $n_1$ and $n_2$ and let $x_i$ be any sequence of integers so that $x_0,x_1,\dots,x_T$ has $n_1$ increases of size $j$ and $n_2$ decreases of size $l$ (where an increase of size $j-l$ counts as both an increase and decrease). For simplicity, I'll call such a sequence an $(n_1,n_2)$-sequence. Then, $$P_{p_1,p_2,j,l,k_1}(X=x) = \begin{cases} \prod_{i=1}^2 p_i^{n_i}(1-p_i)^{T-n_i}&\text{ if } x_i \geq 0\text{ for all }i,\\ 0 &\text{ otherwise.} \end{cases}$$ Then, $$P_{p_1,p_2,j,l,k_1}(X=x,X_0=k_1,X_T=k_2) = \begin{cases} \prod_{i=1}^2 p_i^{n_i}(1-p_i)^{T-n_i}&\text{ if } x_i \geq 0\text{ for all }i\text{ and } x_0=k_1,x_T=k_2,\\ 0 &\text{ otherwise,} \end{cases} = \begin{cases} \prod_{i=1}^2 p_i^{n_i}(1-p_i)^{T-n_i}&\text{ if } x_i \geq 0\text{ for all }i\text{ and } jn_1-ln_2 = k_2-k_1,\\ 0 &\text{ otherwise.} \end{cases}$$ Now consider the following combinatorial function $$f_{j,l}(k_1,k_2,n_1,n_2) = \#\text{ of }(n_1,n_2)\text{-sequences }x\text{ such that }x_i \geq 0\text{ for all }i.$$ I don't know any easy way to compute $f_{j,l}$, but assuming your numbers are small you can probably write some code to directly compute this function. Estimating this function when your numbers are large is more difficult. But, using this function, \begin{align*} P_{p_1,p_2,j,l,k_1}(X_0=k_1,X_T=k_2,n_1,n_2) &= \sum_{(n_1,n_2)\text{-sequences }x}P_{p_1,p_2,j,l,k_1}(X_0=k_1,X_T=k_2,X=x)\\ &=\sum_{(n_1,n_2)\text{-sequences }x}\begin{cases} \prod_{i=1}^2 p_i^{n_i}(1-p_i)^{T-n_i} &\text{ if } x_i \geq 0\text{ for all }i\text{ and } jn_1-ln_2=k_2-k_1,\\ 0 &\text{ otherwise,} \end{cases}\\ &=\begin{cases} f_{j,l}(k_1,k_2,n_1,n_2)\prod_{i=1}^2 p_i^{n_i}(1-p_i)^{T-n_i} &\text{ if } jn_1-ln_2=k_2-k_1,\\ 0 &\text{ otherwise.} \end{cases} \end{align*} Lastly, \begin{align*} L_{0,T}(p_1,p_2,j,l;k_1,k_2) &= P_{p_1,p_2,j,l,k_1}(X_0=k_1,X_T= k_2)\\ &= \sum_{n_1=0}^T\sum_{n_2=0}^T P_{p_1,p_2,j,l,k_1}(X_0=k_1,X_T= k_2,n_1,n_2)\\ &= \sum_{n_1=0}^T\sum_{n_2=0}^T\begin{cases} f_{j,l}(k_1,k_2,n_1,n_2)\prod_{i=1}^2 p_i^{n_i}(1-p_i)^{T-n_i} &\text{ if } jn_1-ln_2=k_2-k_1,\\ 0 &\text{ otherwise.} \end{cases}\\ &=\sum_{n_1=0}^T\sum_{n_2=0}^T f_{j,l}(k_1,k_2,n_1,n_2)\mathbb{I}_{\{jn_1-ln_2=k_2-k_1\}}\prod_{i=1}^2 p_i^{n_i}(1-p_i)^{T-n_i}, \end{align*} where $\mathbb{I}_{A}$ is the indicator function of the event $A$.\

Scenario 2: Surprisingly, given scenario 1 this is much easier. Just for convenience I'm going to change the notation a bit. Let's say you observe the population at times $0 =t_0 < t_1 < \cdots < t_m = T$. At time $t_i$, the population is observed to be $k_i$. Then you may simply notice that (because $X$ is a time-homogeneous Markov chain),

\begin{align*} P_{p_1,p_2,j,l,k_0}(X_{t_i}=k_i\text{ for } i = 0,\dots,m) &= P_{p_1,p_2,j,l,k_0}(X_0=k_0)\prod_{i=1}^{m} P_{p_1,p_2,j,l,k_0}(X_{t_i}=k_i,X_{t_{i-1}}=k_{i-1}|X_{t_{i-1}}=k_{i-1})\\ &=\prod_{i=1}^{m} P_{p_1,p_2,j,l,k_{i-1}}(X_{0}=k_{i-1},X_{t_i-t_{i-1}}=k_i)\\ &=\prod_{i=1}^{m} L_{0,t_i-t_{i-1}}(p_1,p_2,j,l;k_{i-1},k_i), \end{align*} where $L_{0,t_i-t_{i-1}}$ is the likelihood function from scenario 1 with $T = t_i-t_{i-1}$. So, \begin{align*} L_{t_0,\dots,t_m}(p_1,p_2,j,l;k_0,\dots,k_m) &= P_{p_1,p_2,j,l,k_0}(X_{t_i}=k_i\text{ for } i = 0,\dots,m)\\ &=\prod_{i=1}^{m} L_{0,t_i-t_{i-1}}(p_1,p_2,j,l;k_{i-1},k_i). \end{align*}