Consider a jump process $X(t)$ that has the following state-space transtion rates $ \begin{aligned} x \to x+1:& \text{ at rate } \lambda \\ x \to x-1:& \text{ at rate } \mu x \end{aligned} $ Our goal is to be given a set of population observations and use pathwise-gradient estimation to identify a region in parameter space that is most likely to have produced such observations. To fix ideas, suppose that the following observations $(t,x)$ have been made for a single path: $ (\mathsf{t},\mathsf{x}) = \{(t_i, x_i)\}_{i=1}^3 = \{(1,5), (2.5, 8), (5,10)\}. $ One vital aspect of the scenario we are considering is whether or not we are assuming that there is observation error. If there is observation error, then the activities among the increments are correlated through the latent variables. On this particular page, let's assume there is no observation error, and we will pursue an approach that uses independence of the intervals. So in that sense, this can be seen as three observations of the form $(\delta \,:\, x_1 \to x_2))$, which is to say: the elapsed time, the starting value, the ending value and the elapsed time. The particular data written above would be $ (\mathsf{t},\mathsf{x}) = \{(1: 0 \to 5), (1.5 : 5 \to 8), (2.5 : 8 \to 10)\}. $ Our approach to conducting inference will be to simulate the intervals separately, conditioned on the values at the start of each interval. For the sake of the thought-problem, we will proceed from there using exclusively the pathwise-estimation of gradients. #### Likelihood function for a single interval We assume that a path is expressed in terms of the conditionally independent intervals. >[!NOTE] Representations of pathwise likelihood >Let $(t,x)$ be an initial condition of the system and let $(t',x')$ be the final. Then the likelihood function is >$ >L\big(\theta \, ; \, (t,x) \to (t',x')\big) = \mathbb{P}_{\theta}\big(X(t') = t' \, | \, X(t) = x\big). >$ >Let $\Omega$ be the space of all point processes on the nonnegative real line. We can thinking of $X(\cdot \, ; \, \omega)$ as a map from the space $\Omega$ to elements of the space $\mathcal{X}$ of all RCLL (right-continuous with left-hand limits) piecewise integer-valued functions on $[0,\infty)$. >Let $\mathbb{P}_\theta$ be the probability measure on $\Omega$ induced by the continuous-time birth-death Markov chain with parameter vector $\theta = (\lambda, \mu)$. Then we can express >$ >L\big(\theta \, ; \, (t,x) \to (t',x')\big) = \int_\Omega 1_{X(t) = x}(\omega) 1_{X(t') = x'}(\omega) \mathbb{P}_\theta(d \omega) >$ Now, there are two things worth noting. The first is that there exists a measure (called the Jannossy measure) that plays the role for point processes that Lebesgue measure plays for $\mathbb{R}^d$-valued random variables. And second, there is a density function that can be written in very familiar terms. That density gives the likelihood function we are looking for. ##### Jump representation of a path. Any birth-death process can be written in terms of a sequence of jump/duration pairs. To this end, let $\{J_i\}_{i = 1}^\infty$ be a sequence of values in $\{-1,1\}$, and let $\{\tau_i\}_{i =1}^\infty$ be an increasing sequence of values in $(0,\infty)$. We can these pairs to define the sample space $ \Omega = \{(J_i, \tau_i)\}_{i = 1}^\infty $ For each $i$, define the sequence of event values $X_i \coloneqq \sum_{i'\leq i} J_i$. When we wish to emphasize the function mapping elements of $\omega$ to real values, we write $J_i(\omega)$ and $\tau_i(\omega)$. Each element of $\Omega$ is a marked point process and for every $\omega \in \Omega$, we can define $ X(t \, ; \, \omega) = \sum_{i = 1}^\infty J_i(\omega) 1_{[0,\tau_i(\omega)]}(t). $ ##### Density function <--> Statistical Likelihood THIS SUBSECTION NEEDS TO BE EDITED (but it is right in spirit) The density function can be derived from a limit of probabilities of the following kind. For $\epsilon > 0$, define the event (analogous to (or exactly?) cylinder sets in Brownian motion theory) $ (\mathsf{x},\mathsf{t} \,)_{\epsilon,[a,b]} \coloneqq \{\omega \in \Omega \, : \, \bigcap_{i \, : \, t_i \in [a,b]} \big(T_i(\omega) \in (t_i-\epsilon,t_i+\epsilon) \cap X_i(\omega) = x_i\big)\} $ CTMC theory yields that $ \mathbb{P}_\theta\big(X \in (\mathsf{x},\mathsf{t})_{\epsilon, [0,T]}\big) = \prod_{\{i \, : \, t_i \leq T\}} (\lambda 1_{J_i = 1} + \mu 1_{J_i = -1}) \int_{t_i - \epsilon- T_{i-1}}^{t_i + \epsilon - T_{i-1}}\exp(-(\lambda + X(t_{i-1})\mu)s) \mathrm{d} s. $ And the density function, which also serves as the likelihood is $ L\big((\mathsf{x},\mathsf{t} \, ; \, \theta)\big) = \lim_{\epsilon \to 0} \frac{1}{\epsilon} \mathbb{P}_\theta\big(X \in (\mathsf{x},\mathsf{t})_{\epsilon, [0,T]}\big) $ END OF BAD SUBSECTION ##### Likelihood of a given path >[!NOTE] Event-by-event likelihood >Define $L(\theta\,;\, (t,x)\to(t + \tau,x))$ to be the likelihood of the next system value being $(t+\tau,x')$ given $X(t) = x$ and the parameter vector $\theta$. Then >$ >L\big(\theta\,;\, (t,x) \to (t+\tau, x+1))\big) = >(\lambda + \mu x) e^{-(\lambda + \mu x)\tau} \frac{\lambda}{\lambda + \mu x} = \lambda e^{-(\lambda + \mu x)\delta}. >$ >Similarly >$ >\begin{aligned} >L\big(\theta\,;\,(t,x) \to (t+\tau,x-1))\big) &= \mu e^{-(\lambda + \mu x)\tau}, \text{ and} \\ >L\big(\theta\,;\,(t,x) \to (t+\tau,x)) &= e^{-(\lambda + \mu x)\tau}. >\end{aligned} >$ >[!NOTE] Likelihood for a path segment >Define $L(\theta \, ; \, (t,x)\to(t',x') \,|\, (\mathsf{t},\mathsf{x}))$ to be the likelihood being the probability density function for transitioning from $(x,t)$ to $(x't')$ through a sequence of jump events summarized by the vectors $\mathsf{t}$ and $\mathsf{x}$. Then >$ >L(\theta \, ; \, (\mathsf{t},\mathsf{x})) = \lambda^{n^\mathrm{birth}} e^{-\lambda (t'-t)} \mu^{n^{\mathrm{death}}} e^{-\mu \int_t^{t'} X(s) \mathrm{d}s} >$ >where $n^\mathrm{birth}$ and $n^\mathrm{death}$ are the numbers of births and deaths in the interval $(t,t')$ respectively. ##### Example For intuition, suppose that a path is given as follows | Step Duration | Event | $t$ | $X(t)$ | Log-Likelihood | | ------------- | ----- | ---- | ------ | --------------------------------------- | | 0.1 | B | 0.1 | 1 | $\ln(\lambda) - (\lambda)(0.1)$ | | 0.2 | B | 0.3 | 2 | $\ln(\lambda) - (\lambda + \mu)(0.2)$ | | 0.15 | D | 0.45 | 1 | $\ln(\mu) - (\lambda + 2\mu)(0.15)$ | | 0.1 | B | 0.55 | 2 | $\ln(\lambda) - (\lambda + \mu)(0.1)$ | | 0.2 | B | 0.75 | 3 | $\ln(\lambda) - (\lambda + 2\mu)(0.2)$ | | 0.25 | ~ | 1.00 | 3 | $-(\lambda + 3 \mu)(0.25)$ | Next we run a simulation of the interval $(1,2.5)$ with $X(1)$ = 5. | Step Duration | Event | $t$ | $X(t)$ | Log-likelihood | | ------------- | ----- | ---- | ------ | --------------------------------------- | | 0.05 | D | 1.05 | 4 | $\ln(\mu) - (\lambda + 5\mu)(0.05)$ | | 0.25 | B | 1.30 | 5 | $\ln(\lambda) - (\lambda + 4\mu)(0.25)$ | | 0.3 | B | 1.35 | 6 | $\ln(\lambda) - (\lambda + 5\mu)(0.3)$ | | 0.1 | D | 1.45 | 5 | $\ln(\mu) - (\lambda + 6\mu)(0.1)$ | | 0.15 | B | 1.60 | 6 | $\ln(\lambda) - (\lambda + 5\mu)(0.15)$ | | 0.1 | B | 1.70 | 7 | $\ln(\lambda) - (\lambda + 6\mu)(0.1)$ | | 0.2 | B | 1.90 | 8 | $\ln(\lambda) - (\lambda + 7\mu)(0.2)$ | | 0.15 | D | 2.05 | 7 | $\ln(\mu) - (\lambda + 7\mu)(0.15)$ | | 0.1 | B | 2.15 | 8 | $\ln(\lambda) - (\lambda + 7\mu)(0.1)$ | | 0.2 | B | 2.35 | 9 | $\ln(\lambda) - (\lambda + 8\mu)(0.1)$ | | 0.15 | ~ | 2.50 | 9 | $-(\lambda + 9\mu)(0.05)$ | So, for this particular path, $X(1) = 3$ and $X(2.5) = 6$, which means that it undershot the target at both time points. For the sake of argument, let's suppose that this path was generated with the parameter set $\theta_1 = (\lambda_1 = 4, \mu = 1)$, so to hit the desired targets, either the birth rate is too low, or the death rate is too high, or some combination thereof. >[!NOTE] Summary: Log-Likelihood of a birth-death path and sensitivities >Suppose that all transitions of a birth-death process with the above mentioned rates are observed over a time interval $[a,b]$. Let $n^\mathrm{birth}$ denote the number of birth events and $n^\mathrm{death}$ denote the number of death events in this interval. Then the log-likelihood of a path $\mathsf{X} = \{X(t)\}_{t \in [a,b]}$ can be written >$ >\ell(\mathsf{X} \, ; \, \theta) = n^\mathrm{birth} \ln(\lambda) - \lambda T + n^\mathrm{death} \ln(\mu) - \mu \int_0^T \!\! X(t) \mathrm{d} t. >$ >The partial derivatives with respect to $\lambda$ and $\mu$ are then >$ >\begin{aligned} >\partial_\lambda \ell(\mathsf{X} \, ; \, \theta) &= \frac{n^\mathrm{birth}}{\lambda} - T; \\ >\partial_\mu \ell(\mathsf{X} \, ; \, \theta) &= \frac{n^\mathrm{death}}{\mu} - \int_0^T X(t) \mathrm{d} t. \\ >\end{aligned} >$ #### Loss and gradient descent Suppose that we are given $n$ population averages that will be used as targets at times $(t_i, x_i)$ Define $ \Psi(\theta) \coloneqq \sum_{i=1}^n \frac{1}{2} \Big((x_i - x_{i-1})- \mathbb{E}_\theta\big(X(t_i) - X(t_{i-1}) \, | \, X(t_{i-1}) = x_{i-1}\big)\Big)^2 $ which is equivalent to $ \Psi(\theta) \coloneqq \sum_{i=1}^n \frac{1}{2} \Big(\Delta x_i - \mathbb{E}_\theta^{x_{i-1}}\big(X(\Delta t_i)\big)\Big)^2 $ where $\Delta x_i \coloneqq x_i - x_{i-1}$ and $\Delta t_i \coloneqq t_i - t_{i-1}$ and the superscript of the expectation indicates the initial condition. It follows that $ \nabla_\theta \Psi(\theta) = -\sum_{i=1}^n \Big(\Delta x_i - \mathbb{E}_\theta^{x_{i-1}}\big(X(\Delta t_i)\big)\Big) \nabla_\theta \mathbb{E}_\theta^{x_{i-1}}\big(X(\Delta t_i)\big). $ Then, using the notation adopted at the top of this note (NEEDS REFINING) $ \mathbb{E}_{\theta}^x(X(t)) = \int_\Omega X(t \, ; \, \omega) L\big(\theta \, ; \, \mathsf{X}(\omega),x)\big) \mathrm{d} \Lambda(\omega) $ where $\Lambda(\omega)$ is the Janossy measure (fill in details from point process theory literature). Because the likelihood function has the nice form that it does, derivatives pass through the integral and we have $ \begin{aligned} \nabla_\theta \mathbb{E}_\theta(X(t)) &= \int_\Omega X(t \, ; \, \omega) \nabla_\theta L\big(\theta \, ; \, \mathsf{X}(\omega),x)\big) L\big(\theta \, ; \, \mathsf{X}(\omega),x)\big) \mathrm{d} \Lambda(\omega) \\ &= \mathbb{E}_\theta^x\big(X(t) q_\theta(\mathsf{X})\big) \end{aligned} $ where $ q_\theta(X(t)) = \nabla_\theta \, \ell(\theta \, ; \, \mathsf{X}). $ which was derived in the previous section. So the gradient descent algorithm (with learning rate $\alpha$) is $ \hat \theta_{k + 1} = \hat \theta_k + \alpha \sum_{i=1}^n \big((\Delta x_i - \mathbb{E}_\theta^{x_{i-1}}(X(\Delta t_i)) \,\mathbb{E}_\theta^{x_{i-1}}\big(X(\Delta t_i) q_\theta(\mathsf{X})\big) $ #### Monte Carlo approximation for expected value and gradient *Sorry, running out of steam.* Finally, we can write $ \mathbb{E}_\theta^{x}\big(X(t)\big) \approx \frac{1}{M} \sum_{m=1}^M X_m(t) $ and $ \nabla_\theta \mathbb{E}_\theta^{x}(X(t)) = \mathbb{E}_\theta^{x}(X(t) q_\theta(\mathsf{X})) \approx \frac{1}{M} \sum_{m=1}^M X_m(t) \binom{\frac{n_{m}^\mathrm{births}}{\lambda} - t}{\frac{n_{m}^\mathrm{deaths}}{\mu} - \int_0^t \! X_m(s) \mathrm{d} s}. $ where the sample paths have been generated with the appropriate initial condition.