Fast Bayesian parameter estimation for stochastic logistic growth models
Jonathan Heydari
Conor Lawless
David A. Lydall
Darren J. Wilkinson
UK
The transition density of a stochastic, logistic population growth model with multiplicative intrinsic noise is analytically intractable. Inferring model parameter values by fitting such stochastic differential equation (SDE) models to data therefore requires relatively slow numerical simulation. Where such simulation is prohibitively slow, an alternative is to use model approximations which do have an analytically tractable transition density, enabling fast inference. We introduce two such approximations, with either multiplicative or additive intrinsic noise, each derived from the linear noise approximation of the logistic growth SDE. After Bayesian inference we find that our fast LNA models, using Kalman filter recursion for computation of marginal likelihoods, give similar posterior distributions to slow arbitrarily exact models. We also demonstrate that simulations from our LNA models better describe the characteristics of the stochastic logistic growth models than a related approach. Finally, we demonstrate that our LNA model with additive intrinsic noise and measurement error best describes an example set of longitudinal observations of microbial population size taken from a typical, genome-wide screening experiment. \\ \noindent Keywords: Kalman Filter; Linear Noise Approximation; Logistic; Population Growth; Stochastic Modelling;
paper.tex
\documentclass[12pt]{article}
\usepackage{numprint}
\usepackage{amsmath}
\usepackage{graphicx}
\usepackage{tabularx}
\usepackage{subfigure}
\usepackage{float}
\usepackage{longtable}
\usepackage{color}
\usepackage{times,natbib}
\usepackage{amsfonts}
\usepackage{verbatim}
\usepackage{fixltx2e}
\usepackage[labelfont=bf]{caption}
\usepackage[hyperindex,breaklinks]{hyperref}%\topmargin=-0.5in
\usepackage{breakurl}

\topmargin=-0.8in %\topmargin=-1.15in
\oddsidemargin=0in
\evensidemargin=0in
\textwidth=6.5in
\textheight=9.5in

%\newcommand{\hl}[1]{\colorbox{yellow}{#1}}
%\newcommand{\hl}[1]{\textcolor{red}{#1}}
\newcommand{\hl}[1]{\textcolor{black}{#1}}

\title{Fast Bayesian parameter estimation for stochastic logistic growth models}
\author{Jonathan Heydari \and Conor Lawless \and David A. Lydall \and Darren J. Wilkinson\\[1ex] Newcastle University, UK}
%\date{\today}

\begin{document}
 \maketitle
%1\address{Newcastle University, UK}
%2 \email{j.j.heydari@ncl.ac.uk}
%\doublespacing
\input{sections/abstract} %summary 
\input{sections/introduction} 
\input{sections/roman} 
\input{sections/LNAM} 
\input{sections/LNAA} 
\input{sections/application} 
\input{sections/conclusion} 
\appendix
\input{sections/appendix} % An appendix.
\section*{Acknowledgments}
This research was supported by grants from both the Biotechnology and Biological Sciences Research Council and the Medical Research Council.\vspace*{-8pt}
\bibliography{references} % References.
\bibliographystyle{Chicago}
\end{document}

abstract.tex
\begin{abstract}
The transition density of a stochastic, logistic population growth model with multiplicative intrinsic noise is analytically intractable.
Inferring model parameter values by fitting such stochastic differential equation (SDE) models to data therefore requires relatively slow numerical simulation.
Where such simulation is prohibitively slow, an alternative is to use model approximations which do have an analytically tractable transition density, enabling fast inference.
We introduce two such approximations, with either multiplicative or additive intrinsic noise, each derived from the linear noise approximation of the logistic growth SDE. 
After Bayesian inference we find that our fast LNA models, using Kalman filter recursion for computation of marginal likelihoods, give similar posterior distributions to slow arbitrarily exact models.
We also demonstrate that simulations from our LNA models better describe the characteristics of the stochastic logistic growth models than a related approach.
Finally, we demonstrate that our LNA model with additive intrinsic noise and measurement error best describes an example set of longitudinal observations of microbial population size taken from a typical, genome-wide screening experiment.
\\
\noindent \textbf{\emph{Keywords}: Kalman Filter; Linear Noise Approximation; Logistic; Population Growth; Stochastic Modelling;}

\end{abstract}

\begin{comment}
The transition density of a stochastic, logistic population growth model is analytically intractable.
Inferring model parameter values by fitting such SDE models to data therefore requires relatively slow numerical simulation.
Where such simulation is prohibitively slow, an alternative is to use model approximations which do have an analytically tractable transition density, enabling inference by Gibbs sampling, for example.
We introduce two such approximations, each derived from the linear noise approximation (LNA) of the logistic population growth SDE with either multiplicative or additive intrinsic noise.  We show that a related approach is a specific case (a zeroth order approximation) of one of our approximate models.
We demonstrate that simulations from our LNA models are better than their zeroth order equivalent and compare posterior distributions after Bayesian inference with all three approximate models, using Kalman filter recursion for computation of marginal likelihoods.
Finally, we demonstrate that our LNA model with additive intrinsic noise and measurement error best describes an example set of longitudinal observations of microbial population size taken from a typical, genome-wide screening experiment.
\end{comment}
appendix.tex
\appendix
\setcounter{figure}{0}
\renewcommand\thefigure{\thesection.\arabic{figure}}    
\setcounter{table}{0}
\renewcommand\thetable{\thesection.\arabic{table}}    

\section{\label{app:LNAM_sol} LNAM Solution}
\setcounter{figure}{0}
\setcounter{table}{0}
First we look to solve $dZ_t$, given in equation (\ref{eq:LNAM_dz}).
We define $f(t)=-\frac{baPe^{aT}}{bP(e^{aT}-1)+a}$ to obtain the following,
\begin{equation*}
dZ_t=f(t)Z_tdt+\sigma dW_t.
\end{equation*}
Define a new process $U_t=e^{-\int^t_{t_0}f(s)ds}Z_t$ and solve the integral,
\begin{equation*}
\int^t_{t_0}f(s)ds=\int^t_{t_0}-\frac{baPe^{aS}}{bP(e^{aS}-1)+a}ds=\log\left(\frac{a}{bP(e^{aT}-1)+a}\right),
\end{equation*}
where, $S=s-{t_0}$ and $T=t-{t_0}$.
Apply the chain rule to $U_t$,
\begin{equation*}
U_t=e^{-\int^t_{t_0}f(s)ds}dZ_t-f(t)e^{-\int^t_{t_0}f(s)ds}Z_tdt.
\end{equation*}
Now substitute in $dZ_t=f(t)Z_tdt+\sigma dW_t$ and simplify to give
\begin{equation*}
U_t= e^{-\int^t_{t_0}f(s)ds}\sigma dW_t.
\end{equation*}
Apply the following notation
$\phi(t)=e^{\int^t_{t_0}f(s)ds}=\frac{a}{bP(e^{aT}-1)+a}$ and $\psi(t)=\sigma$ to give
\begin{equation*}
U_t=\phi(t)^{-1}\psi(t) dW_t.
\end{equation*}
$U_t$, has the following solution,
\begin{equation*}
U_t=U_0+\int^t_{t_0} \phi(s)^{-1} \psi(s)dW_s.
\end{equation*}
As $U_t=\phi(t)^{-1}Z_t$, $Z_t$ then has the following solution,
\begin{equation*}
Z_t=\phi(t)\left[Z_0+\int^t_{t_0}\phi(s)^{-1}\psi(s) dW_s\right].
\end{equation*}
Finally, the distribution at time t is $Z_t|Z_0\sim N(M_t,E_t)$, where \\
$M_t=\phi(t)Z_0$ and 
$E_t=\phi(t)\int^t_{t_0}\left[{\phi(s)}^{-1}\psi(s)\right]\left[{\phi(s)}^{-1}\psi(s)\right]^{T}ds\:\phi(t)^T$.
\\
Further, $M_t=\frac{a}{bP(e^{aT}-1)+a}Z_0 $ and 
$E_t=\sigma^2\left[\frac{a}{bP(e^{aT}-1)+a}\right]^2
\int^t_{t_0}\left[
\frac{a}{bP(e^{aS}-1)+a}
\right]^{-2}$.\\
As 
$\int^t_{t_0}\left[
\frac{a}{bP(e^{aS}-1)+a}
\right]^{-2}ds=\frac{
b^2P^2e^{2aT}
+4bP(a-bP)e^{aT}
+2at(a-bP)^2
}{2a^3}
-
\frac{
b^2P^2
+4bP(a-bP)
+2at_0(a-bP)^2}{2a^3}
$,
\begin{align*}
E_t=&\sigma^2\left[\frac{a}{bP(e^{aT}-1)+a}\right]^2
\left[
\frac{
b^2P^2(e^{2aT}-1)
+4bP(a-bP)(e^{aT}-1)
+2aT(a-bP)^2
}{2a^3}
\right]\\
=&\sigma^2\left[\frac{
b^2P^2(e^{2aT}-1)
+4bP(a-bP)(e^{aT}-1)
+2aT(a-bP)^2
}
{2a\left(bP(e^{aT}-1)+a\right)^2}\right].
\end{align*}
Taking our solutions for $V_t$ (\ref{eq:LNAM_det_sol}) and $Z_t$, we can now write our solution for the LNA to the log of the logistic growth process (\ref{eq:SDE2}).\\
As $Y_t=V_t+Z_t$, 
\begin{equation*}
Y_t|Y_0\sim \mathcal{N}\left(\log\left[\frac{aPe^{aT}}{bP(e^{aT}-1)+a}\right]+M_t,E_t\right).
\end{equation*}
Note: $\frac{aPe^{aT}}{bP(e^{aT}-1)+a}$ has the same functional form as the solution to the deterministic part of the logistic growth process (\ref{eq_det_sde}) and is equivalent when $\sigma=0$ (such that $a=r-\frac{\sigma^2}{2}=r$).\\
\\
Further, as $Y_t$ is normally distributed, we know $X_t=e^{Y_t}$ will be log normally distributed and
\begin{equation*}
X_t|X_0\sim \log\:\mathcal{N}(\log\left(\frac{aPe^{aT}}{bP(e^{aT}-1)+a}\right)+M_t,E_t).
\end{equation*}
%$x_{t_i}|x_{t_{i-1}}\sim\operatorname{N}\left(x_{t_{i-1}}+\log\left(\frac{aPe^{aT}}{bP(e^{aT}-1)+a}\right)-\log\left(\frac{aPe^{aS}}{bP(e^{aS}-1)+a}\right)+\mu_t-\mu_s,\Xi_tt-\Xi_ss)\right)$\\
%$x_{t_i}|x_{t_{i-1}}\sim \mathcal{N}\left(x_{t_{i-1}}+a\left(t-s\right)+\log\left(\frac{bP(e^{aT}-1)+a}{bP(e^{aS}-1)+a}\right)+\mu_t-\mu_s,\Xi_tt-\Xi_ss\right)$\\
Alternatively set $Q=\left(\frac{\frac{a}{b}}{P}-1\right)e^{at_{0}}$,
\begin{equation*}
X_t|X_{0}\sim \log\:\mathcal{N}(\log\left(\frac{\frac{a}{b}}{1+Qe^{-at}}\right)+M_t,E_t).
\end{equation*}

\noindent From our solution to the log process we can obtain the following transition density
\begin{align*}
\begin{split}
(Y_{t_i}|Y_{t_{i-1}}&=y_{t_{i-1}})\sim\operatorname{N}\left(\mu_{t_i},\Xi_{t_i}\right),\\
\text{where } y_{t_{i-1}}&=v_{t_{i-1}}+z_{t_{i-1}}, Q=\left(\frac{\frac{a}{b}}{P}-1\right)e^{at_{0}},\\
\mu_{t_i}&=v_{t_{i-1}}+\log\left(\frac{1+Qe^{-at_{i-1}}}{1+Qe^{-at_i}}\right)+e^{-a(t_i-t_{i-1})}\frac{1+Qe^{-at_{i-1}}}{1+Qe^{-at_i}}z_{t_{i-1}} \text{ and}\\
\Xi_{t_i}&=\sigma^2\left[\frac{4Q(e^{at_i}-e^{at_{i-1}})+e^{2at_i}-e^{2at_{i-1}}+2aQ^2(t_i-t_{i-1})}{2a(Q+e^{at_i})^2}\right].
\end{split}
\end{align*}

\clearpage

\begin{comment}
\section{\label{app:LNAM_tra} LNAM $dY_t$ Transition Density}
\setcounter{figure}{0}
\setcounter{table}{0}
First we look to derive the transition density for $dZ_t$, given in (\ref{eq:LNAM_dz}).
We define $f(t)=-\frac{baPe^{aT}}{bP(e^{aT}-1)+a}$ to obtain the following,
\begin{equation*}
dZ_t=f(t)Z_tdt+\sigma dW_t.
\end{equation*}
Define $\phi(t)=e^{\int^{t_{i}}_{t_{i-1}}f(s)ds}=\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}$ and $\psi(t)=\sigma$.
The transition density for $Z_{t_i}$ at time $t_i$ is $Z_{t_i}|Z_{t_{i-1}}=z_{t_{i-1}}\sim N(\mu_{t_i},\Xi_{t_i})$,\\
where $\mu_{t_i}=\phi(t)z_{t_{i-1}}$ and 
$\Xi_{t_i}=\phi(t)\int^{t_i}_{t_{i-1}}\left[{\phi(s)}^{-1}\psi(s)\right]\left[{\phi(s)}^{-1}\psi(s)\right]^{T}ds\:\phi(t)^T$.
\\
$\mu_{t_i}=\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}z_{t_{i-1}} $ and 
$\Xi_{t_i}=\sigma^2\left[\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right]^2
\int^{t_{i}}_{t_{i-1}}\left[
\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aS}-1)+a}
\right]^{-2}ds$, where $T_i=t_i-t_0$.\\
\\
As 
$\int^{t_{i}}_{t_{i-1}}\left[
\frac{1}{bP(e^{aS}-1)+a}
\right]^{-2}ds$
\begin{align*}
=&\frac{
b^2P^2e^{2aT_{i}}
+4bP(a-bP)e^{aT_{i}}
+2at_i(a-bP)^2
}{2a}\\
&-
\frac{
b^2P^2e^{2aT_{i-1}}
+4bP(a-bP)e^{aT_{i-1}}
+2at_{i-1}(a-bP)^2
}{2a},
\end{align*}
where $S=s-t_0$ and $T_{i}=t_i-t_0$,
\begin{align*}
\Xi_{t_i}=&\sigma^2\left[\frac{1}{bP(e^{aT_{i}}-1)+a}\right]^2\\
&\times
\left[
\frac{
b^2P^2(e^{2aT_{i}}-e^{2aT_{i-1}})
+4bP(a-bP)(e^{aT_{i}}-e^{aT_{i-1}})
+2a(t_{i}-t_{i-1})(a-bP)^2
}{2a}
\right]
\\
=&\sigma^2\left[\frac{
b^2P^2(e^{2aT_{i}}-e^{2aT_{i-1}})
+4bP(a-bP)(e^{aT_{i}}-e^{aT_{i-1}})
+2a(t_{i}-t_{i-1})(a-bP)^2
}
{2a\left(bP(e^{aT_i}-1)+a\right)^2}\right].
\end{align*}
Taking our solutions for $V_t$ (\ref{eq:LNAM_det_sol}) and $Z_t$, we can now write our transition density for the LNA to the log of the logistic growth process (\ref{eq:SDE2}). 
As $Y_t=Z_t+V_t$,
\begin{align}
Y_{t_i}|Y_{t_{i-1}}&\sim N(v_{t_{i-1}}+V_{t}-V_{t_{i-1}}+\mu_{t_i},\Xi_{t_i})\\
Y_{t_i}|Y_{t_{i-1}}&\sim N(v_{t_{i-1}}+\log\left(e^{a(T_i-T_{i-1})}\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)+\mu_{t_i},\Xi_{t_i}).
\label{eq_transition2}
\end{align}
Or alternatively set $Q=(\frac{\frac{a}{b}}{P}-1)e^{at_0}$,\\
\begin{equation*}
Y_{t_i}|Y_{t_{i-1}}=(v_{t_{i-1}}+z_{t_{i-1}})\sim\operatorname{N}\left(v_{t_{i-1}}+\log\left(\frac{1+Qe^{-at_{i-1}}}{1+Qe^{-at_i}}\right)+\mu_{t_i},\Xi_{t_i}\right)
\end{equation*}
where $\quad a=r-\frac{\zeta^2}{2},\quad b=\frac{r}{K},\quad$ \\ 
$\mu_{t_i}=e^{-a(t_i-t_{i-1})}\frac{1+Qe^{-at_{i-1}}}{1+Qe^{-at_i}}z_{t_{i-1}}$ and $\Xi_{t_i}=\sigma^2\left[\frac{4Q(e^{at_i}-e^{at_{i-1}})+e^{2at_i}-e^{2at_{i-1}}+2aQ^2(t_i-t_{i-1})}{2a(Q+e^{at_i})^2}\right].\\$
\clearpage
\end{comment}

\begin{comment}
\section{\label{app:zero_ord}Zero-order Noise Approximation logistic growth diffusion process}
\setcounter{figure}{0}
\setcounter{table}{0}
After obtaining (\ref{eq:NLGstart}) in Section~\ref{sec:LNAM}, we can derive the ``Roman-Roman" logistic growth diffusion process as follows.
First our expression for $\frac{dV_t}{dt}$, given in (\ref{eq:NLGstart}), is approximated by setting $\sigma^2=0$ such that $a=r-\frac{\sigma^2}{2}=r$,
\begin{equation*}
 \frac{dV_t}{dt}=\frac{a(a-bP)}{bP(e^{aT}-1)+a}=\frac{r(K-P)}{K+P(e^{aT}-1)}.
\end{equation*}
We now solve $dZ_t$, where $dZ_t=dY_t-dV_t$. Expressions for $dY_t$, given in (\ref{eq:SDE2}), and $dV_t$ are known so,
\begin{equation*}
dZ_t=
\left(r-\frac{1}{2}\sigma^2-\frac{r}{K}e^{Y_t}-\frac{a(a-bP)}{bP(e^{aT}-1)+a}\right)dt+\sigma dW_t.
\end{equation*}
We can rearrange to give the following,
\begin{equation*}
dZ_t=\left(\frac{r}{K}\left[\frac{rPe^{rT}}{K+P(e^{aT}-1)}-e^{Y_t}\right]-\frac{1}{2}\sigma^2\right)dt+\sigma dW_t.
\end{equation*}
We now substitute in $Y_t=V_t+Z_t$, 
\begin{equation*}
dZ_t=\left(\frac{r}{K}\left[\frac{rPe^{rT}}{K+P(e^{aT}-1)}-e^{V_t+Z_t}\right]-\frac{1}{2}\sigma^2\right)dt+\sigma dW_t.
\end{equation*}
We now apply a zero order LNA by setting $e^{Z_t}=1$ to obtain,
\begin{equation*}
dZ_t=\left(\frac{r}{K}\left[\frac{rPe^{rT}}{K+P(e^{aT}-1)}-e^{V_t}\right]-\frac{1}{2}\sigma^2\right)dt+\sigma dW_t.
\end{equation*}
Recognizing $e^{V_t}=\frac{rPe^{rT}}{bP(e^{rT}-1)+r}$ we obtain the following,
\begin{equation*}
dZ_t=\left(\frac{r}{K}\left[e^{V_t}-e^{V_t}\right]-\frac{1}{2}\sigma^2\right)dt+\sigma dW_t.
\end{equation*}
We can then simplify to give the following,
\begin{equation*}
dZ_t=-\frac{1}{2}\sigma^2 dt+\sigma dW_t.
\end{equation*}
We know $dY_t=dV_t+dZ_t$ and expressions for both $dZ_t$ and $dV_t$ are known so,
\begin{equation*}
dY_t=\left(\frac{r(K-P)}{K+P(e^{aT}-1)}-\frac{1}{2}\sigma^2\right)dt+\sigma dW_t
\end{equation*}
or alternatively by setting $Q=\left(\frac{K}{P}-1\right)e^{at_{0}}$,
\begin{equation*}
dY_t=\left(\frac{Qr}{e^{rt}+Q}-\frac{1}{2}\sigma^2\right)dt+\sigma dW_t.
\end{equation*}
Which is exactly the RRTR logistic diffusion process presented by \cite{roman}.
\clearpage
\end{comment}

\section{\label{app:LNAA_sol} LNAA Solution}
\setcounter{figure}{0}
\setcounter{table}{0}
First we look to solve $dZ_t$, given in (\ref{eq:LNAA_dz}).
We define $f(t)=a-2bV_t$ and $g(t)=aV_t-bV_t^2-\frac{a^2Pe^{aT}(a-bP)}{\left(bP(e^{aT}-1)+a\right)^2}$ to obtain the following,
\begin{equation*}
dZ_t=\left(g(t)+f(t)Z_t\right)dt+\sigma V_t dW_t.
\end{equation*}
Define a new process $U_t=e^{-\int^t_{t_0}f(s)ds}Z_t$ and solve the integral,
\begin{equation*}
\int^t_{t_0}f(s)ds=\int^t_{t_0}(a-2bV_s)ds=aT-2\log\left(\frac{bP(e^{aT}-1)+a}{a}\right),
\end{equation*}
as $\int^t_{t_0}V_sds=\frac{1}{b}\log \left(\frac{bP(e^{aT}-1)+a}{a}\right)$, where $S=s-{t_0}$ and $T=t-{t_0}$.
Apply the chain rule to $U_t$,
\begin{equation*}
U_t=e^{-\int^t_{t_0}f(s)ds}dZ_t-f(t)e^{-\int^t_{t_0}f(s)ds}Z_tdt.
\end{equation*}
Now substitute in $dZ_t=\left(g(t)+f(t)Z_t\right)dt+\sigma V_t dW_t$ and simplify to give,
\begin{equation*}
U_t=e^{-\int^t_{t_0}f(s)ds}g(t)dt+e^{-\int^t_{t_0}f(s)ds}\sigma V_t dW_t.
\end{equation*}
Apply the following notation $\phi(t)=e^{\int^t_{t_0}f(s)ds}=e^{aT}\left(\frac{a}{bP(e^{aT}-1)+a}\right)^2$ and $\psi(t)=\sigma V_t$ to give, 
\begin{equation*}
U_t=\phi(t)^{-1}g(t)+\phi(t)^{-1}\psi(t) dW_t.
\end{equation*}
$U_t$ has the following solution,
\begin{equation*}
U_t=U_0+\int^t_{t_0}\phi(s)^{-1}g(s)ds +\int^t_{t_0} \phi(s)^{-1} \psi(s)dW_s.
\end{equation*}
As $U_t=\phi(t)^{-1}Z_t$, $Z_t$ has the following solution,
\begin{equation*}
Z_t=\phi(t)\left[Z_0+\int^t_{t_0}\phi(s)^{-1}g(s)ds+\int^t_{t_0}\phi(s)^{-1}\psi(s) dW_s\right].
\end{equation*}
Finally the distribution at time t is $Z_t|Z_0\sim N(M_t,E_t)$, where \\
$M_t=\phi(t)\left(Z_0+\int^t_{t_0}\phi(s)^{-1}g(s)ds\right)$ and 
$E_t=\phi(t)\int^t_{t_0}\left[{\phi(t)}^{-1}\psi(t)\right]\left[{\phi(t)}^{-1}\psi(t)\right]^{T}ds\:\phi(t)^T$.
\begin{align*}
M_t&=\phi(t)\left(Z_0+ 
\int^t_{t_0}\left[\left(\phi(t)\right)^{-1} \left(aV_s-bV_s^2-\frac{a^2Pe^{aS}(a-bP)}{\left(bP(e^{aS}-1)+a\right)^2}\right)\right] ds
\right)  \\
&=\phi(t)\left(Z_0+ 
\int^t_{t_0}\left[\left(\phi(t)\right)^{-1} \left(0\right)\right] ds
\right) \\
&=e^{aT}\left(\frac{a}{bP(e^{aT}-1)+a}\right)^2Z_0
\end{align*}
 and 
\begin{equation*}
E_t=\left(e^{aT}\left(\frac{a}{bP(e^{aT}-1)+a}\right)^2\right)^2
\int^t_{t_0}\left[e^{aS}\left(\frac{a}{bP(e^{aS}-1)+a}\right)^2\right]^{-2} \sigma^2 V_s^2 ds
\end{equation*}
\begin{equation*}
=\sigma^2\left(e^{aT}\left(\frac{a}{bP(e^{aT}-1)+a}\right)^2\right)^2
\int^t_{t_0}\left[
e^{aS}\left(\frac{a}{bP(e^{aS}-1)+a}\right)^2
\right]^{-2} \left[
\frac{aPe^{aS}}{bP(e^{aS}-1)+a}
\right]^{2} ds.
\end{equation*}
\begin{equation*}
=\sigma^2\left(e^{aT}\left(\frac{a}{bP(e^{aT}-1)+a}\right)^2\right)^2
\int^t_{t_0}\left[
e^{-2aS}\left(\frac{a}{bP(e^{aS}-1)+a}\right)^{-4}
\right] \left[
\frac{aPe^{aS}}{bP(e^{aS}-1)+a}
\right]^{2} ds
\end{equation*}
\begin{equation*}
=\sigma^2\left(e^{aT}\left(\frac{1}{bP(e^{aT}-1)+a}\right)^2\right)^2
\int^t_{t_0}\left[a^2P^2
\left(\frac{1}{bP(e^{aS}-1)+a}\right)^{-2}
\right] ds,
\end{equation*}
as $\int^{t_{i}}_{t_{0}}\left(\frac{1}{bP(e^{aS}-1)+a}\right)^{-2}ds
=\frac{
b^2P^2e^{2aT}
+4bP(a-bP)e^{aT}
+2at(a-bP)^2
}{2a}
-
\frac{
b^2P^2
+4bP(a-bP)
+2at_0(a-bP)^2}{2a}$,
\begin{align*}
E_t=&\frac{1}{2}\sigma^2aP^2e^{2aT}\left(\frac{1}{bP(e^{aT}-1)+a}\right)^4\\
&\times
\left[
b^2P^2(e^{2aT}-1)
+4bP(a-bP)(e^{aT}-1)
+2aT(a-bP)^2
\right].
\end{align*}
\noindent From our solution to the process we can obtain the following transition density
\begin{align*}
\begin{split}
(X_{t_i}|X_{t_{i-1}}=&x_{t_{i-1}})\sim N(\mu_{t_i},\Xi_{t_i}),\\
\text{where }x_{t_{i-1}}=&v_{t_{i-1}}+z_{t_{i-1}},\\
\mu_{t_i}=&v_{t_{i-1}}+\left(\frac{aPe^{aT_i}}{bP(e^{aT_i}-1)+a}\right)-\left(\frac{aPe^{aT_{i-1}}}{bP(e^{aT_{i-1}}-1)+a}\right)\\
&+e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2Z_{t_{i-1}}\text{ and}\\
 \Xi_t=&\frac{1}{2}\sigma^2aP^2e^{2aT_i}\left(\frac{1}{bP(e^{aT_i}-1)+a}\right)^4\\
&\times[
b^2P^2(e^{2aT_i}-e^{2aT_{i-1}})
+4bP(a-bP)(e^{aT_i}-e^{aT_{i-1}})\\
&\;\:\:\:\:+2a(t_i-t_{i-1})(a-bP)^2
].
\end{split}
\end{align*}
\clearpage

\begin{comment}
\section{\label{app:LNAA_tra} LNAA $dY_t$ Transition Density}
\setcounter{figure}{0}
\setcounter{table}{0}
First we look to find the transition density for $dZ_t$, given in (\ref{eq:LNAA_dz}).
We define $f(t)=a-2bV_t$ and $g(t)=aV_t-bV_t^2-\frac{a^2Pe^{aT}(a-bP)}{\left(bP(e^{aT}-1)+a\right)^2}$ to obtain the following,
\begin{equation*}
dZ_t=\left(g(t)+f(t)Z_t\right)dt+\sigma V_t dW_t.
\end{equation*}
Define $\phi(t)=e^{\int^{t_{i}}_{t_{i-1}}f(s)ds}=e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2$ and $\psi(t)=\sigma V_t$. The transition density for $Z_t$ at time $t_i$ is
$Z_{t_i}|Z_{t_{i-1}}\sim N(\mu_{t_i},\Xi_{t_i})$,\\
where $\mu_{t_i}=\phi(t)\left(Z_{t_{i-1}}+ 
\int^t_{t_{i-1}}\left[\left(\phi(t)\right)^{-1} \left(0\right)\right] ds
\right)=\phi(t)Z_{t_{i-1}}$\\ 
and
$\Xi_{t_i}=\phi(t)\int^{t_i}_{t_{i-1}}\left[{\phi(s)}^{-1}\psi(s)\right]\left[{\phi(s)}^{-1}\psi(s)\right]^{T}ds\:\phi(t)^T.$
\\
$\mu_{t_i}=e^{\int^{t_{i}}_{t_{i-1}}f(s)ds}=e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2Z_{t_{i-1}} $ and
\\ 
$\Xi_{t_i}=\left(e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2\right)^2
\int^{t_{i}}_{t_{i-1}}\left[e^{a(s-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aS}-1)+a}\right)^2\right]^{-2} \sigma^2 V_s^2 ds$\\
$=\sigma^2\left(e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2\right)^2
\int^{t_{i}}_{t_{i-1}}\left[e^{a(s-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aS}-1)+a}\right)^2\right]^{-2}\left[
\frac{aPe^{aS}}{bP(e^{aS}-1)+a}
\right]^{2} ds$\\
$=\sigma^2\left(e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2\right)^2
\int^{t_{i}}_{t_{i-1}}\left[e^{-2a(s-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aS}-1)+a}\right)^{-4}\right]
\left[
\frac{aPe^{aS}}{bP(e^{aS}-1)+a}
\right]^{2} ds$\\
$=\sigma^2\left(e^{a(t_i-t_{i-1})}\left(\frac{1}{bP(e^{aT_i}-1)+a}\right)^2\right)^2
\int^{t_{i}}_{t_{i-1}}\left[a^2P^2e^{2aT_{i-1}}(\frac{1}{bP(e^{aS}-1)+a})^{-4}\right]\left(\frac{1}{bP(e^{aS}-1)+a}\right)^{2} ds$\\
$=\sigma^2\left(e^{2a(t_i-t_{i-1})}\left(\frac{1}{bP(e^{aT_i}-1)+a}\right)^4\right)
\int^{t_{i}}_{t_{i-1}}\left[a^2P^2e^{2aT_{i-1}}(\frac{1}{bP(e^{aS}-1)+a})^{-4}\right]\left(\frac{1}{bP(e^{aS}-1)+a}\right)^{2} ds$\\
\\
$=\sigma^2\left(e^{2aT_i}\left(\frac{1}{bP(e^{aT_i}-1)+a}\right)^4\right)
\int^{t_{i}}_{t_{i-1}}\left[a^2P^2\left(\frac{1}{bP(e^{aS}-1)+a}\right)^{-2}\right] ds$,\\
\\
As 
\begin{align*}
\int^{t_i}_{t_{i-1}}\left[
\frac{1}{bP(e^{aS}-1)+a}
\right]^{-2}ds
=&\frac{
b^2P^2e^{2aT_i}
+4bP(a-bP)e^{aT_i}
+2at_i(a-bP)^2
}{2a}\\
&-
\frac{
b^2P^2e^{2aT_{i-1}}
+4bP(a-bP)e^{aT_{i-1}}
+2at_{i-1}(a-bP)^2
}{2a},
\end{align*}
where $S=s-t_0$ and $T_{i}=t_i-t_0$,
\begin{align*}
\Xi_t=&\sigma^2a^2P^2\left(e^{2aT_i}\left(\frac{1}{bP(e^{aT_i}-1)+a}\right)^4\right)\\
&\times\left[\frac{
b^2P^2(e^{2aT_i}-e^{2aT_{i-1}})
+4bP(a-bP)(e^{aT_i}-e^{aT_{i-1}})
+2a(t_i-t_{i-1})(a-bP)^2
}{2a}\right]\\
=&\frac{1}{2}\sigma^2aP^2e^{2aT_i}
\left(\frac{1}{bP(e^{aT_i}-1)+a}\right)^4
\\
&\times
\left[
b^2P^2(e^{2aT_i}-e^{2aT_{i-1}})
+4bP(a-bP)(e^{aT_i}-e^{aT_{i-1}})
+2a(t_i-t_{i-1})(a-bP)^2
\right].
\end{align*}
As $X_t=Z_t+V_t$,
\begin{align*}
X_{t_i}|X_{t_{i-1}}&\sim N(x_{t_{i-1}}+V_{t}-V_{t_{i-1}}+\mu_{t_i},\Xi_{t_i})\\
X_{t_i}|X_{t_{i-1}}&\sim N(x_{t_{i-1}}+\left(\frac{aPe^{aT_i}}{bP(e^{aT_i}-1)+a}\right)-\left(\frac{aPe^{aT_{i-1}}}{bP(e^{aT_{i-1}}-1)+a}\right)+\mu_{t_i},\Xi_{t_i})
\end{align*} 
\clearpage
\end{comment}

\section{\label{app:prior_hyp}Prior Hyper-parameters for Bayesian State Space Models}
\input{tables/SDE_priors}

	\clearpage


\section{\label{app:kalman_fil}Kalman Filter}
\setcounter{figure}{0}
\setcounter{table}{0}
To find $\pi(y_{t_{1:N}})$ for the LNAA with normal measurement error we can use the following Kalman Filter algorithm. First we assume the following:
\begin{align*}
\theta_{{t_{i}}}|y_{1:{{t_{i}}}}&\sim \operatorname{N}(m_{t_{i}},C_{t_{i}}),\\
m_{t_{i}}&=a_{t_{i}}+R_{t_{i}}F(F^{T}R_{{t_{i}}}F+U)^{-1}[y_{t_{i}}-F^{T}a_{t_{i}}],\\
C_{t_{i}}&=R_{t_{i}}-R_{t_{i}}F(F^TR_{t_{i}}F+U)^{-1}F^{T}R_{t_{i}}
\end{align*}
and initialize with $m_0=P$ and $C_0=0$. Now suppose that,
\begin{align*}
\theta_{t_{i}}|y_{1:{t_{i-1}}}&\sim \operatorname{N}(a_{t_{i}},R_{t_{i}}),\\
a_{t_{i}}&=G_{{t_{i}}}m_{{t_{i-1}}}\\  
\text{and }R_{t_{i}}&=G_{{t_{i}}}C_{{t_{i-1}}}G_{t_{i}}^T+W_{t_{i}}.
\end{align*}
The transition density distribution, see (\ref{eq:LNAA_tran}) is as follows:
\begin{align*}
\theta_{{t_{i}}}|\theta_{{t_{i-1}}}&\sim\operatorname{N}(G_{{t_{i}}}\theta_{{t_{i-1}}},W_{t_{i}})\\
\text{or equivalently }(X_{t_i}|X_{t_{i-1}}=x_{t_{i-1}})&\sim\operatorname{N}\left(\mu_{t_i},\Xi_{t_i}\right),\text{ where }x_{t_{i-1}}=v_{t_{i-1}}+z_{t_{i-1}},\\
\theta_{t}
 &=
\begin{pmatrix}
  1 \\
	X_{t_{i}}
 \end{pmatrix}
 =
 \begin{pmatrix}
  1 & 0 \\
  H_{\alpha,t_{i}} & H_{\beta,t_{i}}
 \end{pmatrix}
 \begin{pmatrix}
  1 \\
	X_{t_{i-1}}
 \end{pmatrix}
\\&=
 G_{t_{i}}\theta_{{t_{i-1}}},\\
 G_{t_{i}}&= \begin{pmatrix}
 1 & 0\\
 H_{\alpha,t_{i}} & H_{\beta,t_{i}}
 \end{pmatrix}, \quad
  W_{t_{i}}= \begin{pmatrix}
  0 & 0 \\
  0 & \Xi_{t_i}
	\end{pmatrix} \\
\text{where }H_{\alpha,t_{i}}=H_\alpha({t_{i}},{t_{i-1}})=&V_t-V_{t-1}e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2\\
\text{and }H_{\beta,t_{i}}=&H_\beta({t_{i}},{t_{i-1}})=e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2.
\end{align*}
The measurement error distribution is as follows:
\begin{align*}
y_{t_{i}}|\theta_{{t_{i}}}{\sim}&\operatorname{N}(F^T\theta_{{t_{i}}},U)\\
\text{or equivalently }y_{t_{i}}|\theta_{{t_{i}}}{\sim}&\operatorname{N}(X_{{t_{i}}},\sigma_{\nu}^2),\\
 \text{where }
 F=& \begin{pmatrix}
  0 \\
	1 
 \end{pmatrix}\text{ and }
  U=  \sigma_{\nu}^2.
  \end{align*}  
  Matrix Algebra:
\begin{align*}
a_{t_{i}}=&G_{{t_{i}}}m_{{t_{i-1}}}\\
=&\begin{pmatrix}
 1 & 0\\
H_{\alpha,t_{i}} & H_{\beta,t_{i}}
 \end{pmatrix}
 \begin{pmatrix}
  1\\
	m_{{t_{i-1}}} 
 \end{pmatrix}
=\begin{pmatrix}
  1\\
	H_{\alpha,t_{i}}+H_{\beta,t_{i}}m_{{t_{i-1}}}
 \end{pmatrix}
  \end{align*}  
  \begin{align*}
R_{t_{i}}&=G_{{t_{i}}}C_{{t_{i-1}}}G_{t_{i}}^T+W_{t_{i}}
\\
 &=
  \begin{pmatrix}
  0 & 0\\
  0 & {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2
 \end{pmatrix}
 +
  \begin{pmatrix}
  0 & 0\\
  0 & \Xi_{t_i}
 \end{pmatrix}
   =\begin{pmatrix}
  0 & 0\\
  0 & {H_{\beta,t_{i}}}^2 c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}
\end{align*}
  \begin{align*}
  C_{t_{i-1}}&=
  \begin{pmatrix}
	0 & 0\\
  0 & c_{{t_{i-1}}}^2
 \end{pmatrix}
    \end{align*}
 \begin{align*}
  R_{t_{i}}F(F^{T}R_{{t_{i}}}F+U)^{-1}=&
  \begin{pmatrix}
	0 & 0\\
  0 & {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}
 \begin{pmatrix}
  0 \\
	1
 \end{pmatrix}
\\
&\times
 \left[
 \begin{pmatrix}
 0 & 1\\
 \end{pmatrix}
   \begin{pmatrix}
	0 & 0\\
  0 & {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}
  \begin{pmatrix}
  0 \\
	1
 \end{pmatrix}
 +\sigma_{\nu}^2
 \right]^{-1}
 \\
  =&\left[
    \begin{pmatrix}
{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i} +\sigma_{\nu}^2
 \end{pmatrix}
 \right]^{-1}
 \begin{pmatrix}
	0 \\
  {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}
\end{align*}
\begin{align*}
m_{t_{i}}=&a_{t_{i}}+R_{t_{i}}F(F^{T}R_{{t_{i}}}F+U)^{-1}[y_{t_{i}}-F^{T}a_{t_{i}}]\\
=&\begin{pmatrix}
  1\\
	H_{\alpha,t_{i}}+H_{\beta,t_{i}}m_{{t_{i-1}}}
 \end{pmatrix}\\
 &+
 \left[
    \begin{pmatrix}
{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i} +\sigma_{\nu}^2
 \end{pmatrix}
 \right]^{-1}
 \begin{pmatrix}
	0 \\
  {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}%partA
 \left[
 y_{t_{i}}-
 \begin{pmatrix}
  0 & 1 \\
 \end{pmatrix}
 \begin{pmatrix}
  1\\
	H_{\alpha,t_{i}}+H_{\beta,t_{i}}m_{{t_{i-1}}}
 \end{pmatrix}
 \right]\\
   =&\begin{pmatrix}
	0\\
   H_{\alpha,t_{i}}+H_{\beta,t_{i}}m_{{t_{i-1}}}+\frac{{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}}{ {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}+\sigma_{\nu}^2}\left[y_{t_{i}}-H_{\alpha,t_{i}}-H_{\beta,t_{i}}m_{{t_{i-1}}}\right]
 \end{pmatrix}
 \end{align*}
   \begin{align*}
   C_{t_{i}}=&R_{t_{i}}-R_{t_{i}}F(F^TR_{t_{i}}F+U)^{-1}F^{T}R_{t_{i}}\\
 =&\begin{pmatrix}
	0 & 0\\
  0 & {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}\\
&-
 \left[
    \begin{pmatrix}
{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i} +\sigma_{\nu}^2
 \end{pmatrix}
 \right]^{-1}
 \begin{pmatrix}
	0\\
  {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}%partA
 \left[ 
  \begin{pmatrix}
 0 & 1
 \end{pmatrix}
  \begin{pmatrix}
  0 & 0\\
  0 & {H{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}
 \end{pmatrix}
 \right]
 \\
     =&
       \begin{pmatrix}
  0 & 0 \\
  0 & {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}-\frac{\left({H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}\right)^2}{{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}+\sigma_{\nu}^2}
 \end{pmatrix}
      \end{align*}
\\
With $m_{t_{i}}$ and $C_{t_{i}}$ for $i=1:T$, we can evaluate $a_{t_{i}}$, $R_{t_{i}}$ and $\pi(x_{t_{i}}|y_{t_{1:(i-1)}})$ for $i=1:T$.
We are interested in $\pi(y_{t_{1:i}})=\prod^T_{i=1}\pi(y_{t_{i}}|y_{t_{1:(i-1)}})$, where
$\pi(y_{t_{i}}|y_{t_{1:(i-1)}})=\int_{x}\pi(y_{t_{i}}|x_{t_{i}})\pi(x_{t_{i}}|y_{t_{1:(i-1)}})dx_{t_{i}}$ gives a tractable Gaussian integral.
Finally, 
\begin{align*}
\log\pi(y_{{t_{1:(i-1)}}})&=\sum^T_{i=1}\log\pi(y_{t_{i}}|y_{{t_{1:(i-1)}}})\\
&=\sum^T_{i=1}\left[-\log\left({\sqrt{2\pi(\sigma_{f}^2+\sigma_{g}^2)}}\right){-\frac{(\mu_f-\mu_g)^2}{2(\sigma_{f}^2+\sigma_{g}^2)}}\right],
\end{align*}
\begin{align*}
\text{where }\mu_f-\mu_g=&y_{t_{i}}-a_{t_{i}}=y_{t_{i}}-H_{\alpha,t_{i}}-H_{\beta,t_{i}}m_{t_{i-1}}\\
\text{and }\sigma_{f}^2+\sigma_{g}^2=&\sigma_{\nu}^2+R_{t_{i}}=\sigma_{\nu}^2+{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}.
\end{align*}
%$\mu_f-\mu_g=y_t-\left[H+m_{t-1}+\frac{c_{t-1}^2+\zeta^{2}(t-s)}{ c_{t-1}^2+\zeta^{2}(t-s)+\sigma_{\nu}^2}\left[y_t-(H+m_{t-1})\right]\right]$\\
%$\sigma_{f}^2+\sigma_{g}^2=\sigma_{\nu}^2+\left[c_{t-1}^2+\zeta^{2}(t-s)-\frac{\left(c_{t-1}^2+\zeta^{2}(t-s)\right)^2}{c_{t-1}^2+\zeta^{2}(t-s)+\sigma_{\nu}^2} \right]$\\
  \clearpage
\underline{Procedure}\\
\\
1. Set $i=1$. Initialize $m_0=\log(P)$ and $C_0=0$. \\
\\
2. Evaluate and store the following log likelihood term: 
\begin{align*}
\log\pi(y_{t_{i}}|y_{t_{1:(i-1)}})=&\left[-\log\left({\sqrt{2\pi(\sigma_{f}^2+\sigma_{g}^2)}}\right){-\frac{(\mu_f-\mu_g)^2}{2(\sigma_{f}^2+\sigma_{g}^2)}}\right],\\
\text{where }\mu_f-\mu_g=&y_{t_{i}}-H_{\alpha,t_{i}}-H_{\beta,t_{i}}m_{{t_{i-1}}}
\text{ and }\sigma_{f}^2+\sigma_{g}^2=\sigma_{\nu}^2+{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}.
\end{align*}
3. Create and store both $m_{t_i}$, and $C_{t_{i}}$,
\begin{align*}
 \text{where }m_{t_{i}}=&H_{\alpha,t_{i}}+H_{\beta,t_{i}}m_{t_{i-1}}+\frac{{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}}{ {H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}+\sigma_{\nu}^2}\left[y_{t_{i}}-H_{\alpha,t_{i}}-H_{\beta,t_{i}}m_{{t_{i-1}}}\right]\\
\text{and }c_{{t_{i}}}^2=&{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}-\frac{\left({H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}\right)^2}{{H_{\beta,t_{i}}}^2c_{{t_{i-1}}}^2+\Xi_{t_i}+\sigma_{\nu}^2}.
 \end{align*}
\\
4. Increment $i$, $i$=$(i+1)$ and repeat steps 2-3 till $\log\pi(y_{t_{N}}|y_{t_{1:(N-1)}})$ is evaluated.\\
\\
5. Calculate the sum:
\begin{equation*}
\log\pi(y_{{t_{1:N}}})=\sum^N_{i=1}\log\pi(y_{t_{i}}|y_{{t_{1:(i-1)}}}).
\end{equation*}


\clearpage

\section{\label{app:diag_sim_n}Parameter posterior means}
\setcounter{figure}{0}
\setcounter{table}{0}
\begin{figure}[h!]
  \centering
\includegraphics[width=14cm]{img_sde/diag_sim_n}
\caption{Trace, auto-correlation and density plots for the LNAA model parameter posteriors (sample size = 1000, thinning interval = 4000), see Figure~\ref{sim}D.
Posterior density (black), prior density (dashed blue) and true parameter values (red) are shown in the right hand column.
}
\end{figure}	
	
	


application.tex
\section{\label{sec:SDE_application}Simulation and Bayesian inference for logistic SDE and approximations}
To test which of the three approximate models best represent the SLGM, we first compare simulated forward trajectories from the RRTR, LNAM and LNAA with simulated forward trajectories from the SLGM (Figure~\ref{4nonu}).
We use the Euler-Maruyama method \citep{embook} with very fine discretisation to give arbitrarily exact simulated trajectories from each SDE.

The LNAA and LNAM trajectories are visually indistinguishable from the SLGM (Figure~\ref{4nonu} A,C \& D).
On the other hand, population sizes simulated with the RRTR display large deviations from the mean as the population approaches stationary phase (Figure~\ref{4nonu}A \& B). 
Figure~\ref{4nonu}E further highlights the increases in variation as the population approaches stationary phase for simulated trajectories of the RRTR, in contrast to the SLGM and LNA models.%from the SLGM, LNAM or LNAA. 
\begin{figure}[h!]
  \centering
\includegraphics[width=14cm]{img_sde/4nonuALT}
\caption{\label{4nonu} 
Forward trajectories (No. of simulations=100) for logistic SDE and approximations.
See Table~\ref{app:sde_val_fur} for parameter values. 
A) The stochastic logistic growth model (SLGM). 
B) The \cite{roman} (RRTR) approximation. 
C) The linear noise approximation with multiplicative intrinsic noise (LNAM).
D) The linear noise approximation with additive intrinsic noise (LNAA).
E) Standard deviations of simulated trajectories over time for the SLGM (black), RRTR (red), LNAM (green) and LNAA (blue).
}
\end{figure}			

\subsection{\label{sec:simulation_stu}Bayesian parameter inference with approximate models}
To compare the quality of parameter inference using each of these approximations we simulated synthetic time-course data from the SLGM and combined this with either lognormal or normal measurement error. % [[PERHAPS WE SHOULD GIVE THE APPROXIMATE SDES PLUS MATCHING ERROR MODELS (REQUIRED TO MAINTAIN LINEAR GAUSSIAN STRUCTURE) NAMES?  ALSO, PERHAPS THEY SHOULD BE PRESENTED AT THE END OF THE SECTIONS DESCRIBING THE APPROXIMATIONS (SINCE THEY ARE REQUIRED FOR THE TYPE OF INFERENCE WE CONSIDER)]]
Carrying out Bayesian inference with broad priors (see (\ref{app:LNAM_sta_spa_mod}) and (\ref{app:LNAA_sta_spa_mod})) we compared the parameters recovered using each approximation with those used to generate the synthetic dataset. 
The synthetic time-course datasets consist of 27 time points generated using the Euler-Maruyama method with very fine intervals \citep{euler_maruyama}.

We formulate our inference problem as a dynamic linear state space model.  
To allow fast parameter inference we take advantage of a linear Gaussian structure and construct a Kalman filter recursion for marginal likelihood computation (Appendix~\ref{app:kalman_fil}).
We therefore assume lognormal (multiplicative) error for the RRTR and LNAM, and for the LNAA we assume normal (additive) measurement error.
Dependent variable $y_{t_i}$ and independent variable $\{t_{i},i=1,...,N\}$ are data input to the model (where $t_i$ is the time at point $i$ and $N$ is the number of time points). 
$X_t$ is the state process, describing the population size.
See Table~\ref{table:SDE_priors} for prior hyper-parameter values.

For the RRTR and LNAM,
\begin{align}
\log(y_{t_i}) &\sim \operatorname{N}(X_{t_i},{\nu}^{2} ),\notag\\
(X_{t_i}|X_{t_{i-1}}=x_{t_{i-1}})&\sim\operatorname{N}\left(\mu_{t_i},\Xi_{t_i}\right), \text{ where } x_{t_{i}}=v_{t_{i}}+z_{t_{i}},\label{app:LNAM_sta_spa_mod}
\end{align}
$\mu_{t_i}$ and $\Xi_{t_i}$ are given by (\ref{eq:RRTR_tran}) and (\ref{eq:LNAM_tran}) for the RRTR and LNAM respectively. Priors are as follows:
\begin{align*}
\log X_0 \equiv \log~P &\sim \operatorname{N}({\mu}_P,{\tau_P}^{-1}),& \quad\quad
\log~K &\sim \operatorname{N}({\mu}_K,{\tau_K}^{-1}),& \quad\quad
\log~r &\sim \operatorname{N}({\mu}_r,{\tau_r}^{-1}),\\
\log~\nu^{-2} &\sim \operatorname{N}({\mu}_\nu,{\tau_\nu}^{-1}),&
\log~{\sigma^{-2}} &\sim \operatorname{N}({\mu}_\sigma,{\tau_\sigma}^{-1})I_{[1,\infty]}.&
\end{align*}

For the LNAA, 
\begin{align}
y_{t_i} &\sim \operatorname{N}(X_{t_i},{\nu}^{2} ),\notag\\
(X_{t_i}|X_{t_{i-1}}=x_{t_{i-1}})&\sim\operatorname{N}\left(\mu_{t_i},\Xi_{t_i}\right)
, \text{ where } x_{t_{i}}=v_{t_{i}}+z_{t_{i}},\label{app:LNAA_sta_spa_mod}
\end{align}
$\mu_{t_i}$ and $\Xi_{t_i}$ are given by (\ref{eq:LNAA_tran}). Priors are as in (\ref{app:LNAM_sta_spa_mod}).
Our prior for $\log~{\sigma^{-2}}$ is truncated below 1 to avoid unnecessary exploration of extremely low probability regions, which could be caused when there are problems identifying $\nu$, for example when $\log~\nu^{-2}$ takes large values, and to ensure that intrinsic noise does not dominate the process.
The truncation limit was chosen by visual inspection.
%. The choice of the value 1 is made by visual inspection and ensures intrinsic noise does not dominate the process.

To see how the inference from our approximate models compares with slower ``exact'' models, we consider Euler-Maruyama approximations \citep{euler_maruyama} of (\ref{eq_det_sde}) and of the log transformed process, using fine intervals.
Given these approximations we can construct a state space model for an ``exact'' SLGM with lognormal measurement error (SLGM+L) and similarly for the SLGM with normal measurement error (SLGM+N), priors are as in (\ref{app:LNAM_sta_spa_mod}).

Model fitting is carried out using standard MCMC techniques (the Gibbs sampler) \citep{gamerman}.
Posterior means are used to obtain point estimates and standard deviations for describing variation of inferred parameters. %[[REMOVE]]
%Standard model diagnostics are used to check the output from our Bayesian models. 
The Heidelberger and Welch convergence diagnostic \citep{Heidelberger} and effective sample size diagnostic \citep{coda} are used to determine whether convergence has been achieved for all parameters. 

Computational times for convergence of our MCMC schemes (code is available at \url{https://github.com/jhncl/LNA.git}) can be compared using estimates for the minimum effective sample size per second (ESS\textsubscript{min}/sec). 
The average ESS\textsubscript{min}/sec of our approximate model (coded in C) is $\sim$100 and ``exact'' model $\sim$1 (coded in JAGS \citep{rjags} with 15 imputed states between time points, chosen to maximise ESS\textsubscript{min}/sec).
We find that our C code is typically twice as fast as the simple MCMC scheme used by JAGS, indicating that our inference is ${\sim}50\times$  faster than an ``exact'' approach.
A more efficient ``exact'' approach could speed up further, say by another factor of 5, but our approximate approach will at least be an order of magnitude faster.
We use a burn-in of 600,000 and a thinning of 4,000 to obtain a final posterior sample size of 1,000 for MCMC convergence of all our models.

%In Figure~\ref{self} we can investigate how well the models can recover simulated data with the measurement error assumptions chosen.
%Further, Table~\ref{tab:real_par} demonstrates that all models are able to recover the true parameter values. 
%\begin{figure}[h!]
%  \centering
%\includegraphics[width=14cm]{img_sde/9selfALT}
%\caption{\label{self} 
%Simulated time-courses from each of the approximations to the logistic growth SDE with measurement error (black). Forward simulated time-courses with measurement error, where parameter values are taken from posterior samples (sample size=1000). See (\ref{app:LNAM_sta_spa_mod}) or (\ref{app:LNAA_sta_spa_mod}) for model and Table~\ref{table:SDE_priors} for prior hyper-parameter values. See Table~\ref{tab:real_par} for parameter posterior means and true values.  
%A) Roman-Roman, Torres-Ruiz approximation (RRTR) with lognormal error (red)
%B) Linear noise approximation with multiplicative intrinsic noise (LNAM) and lognormal measurement error (green).
%C) Linear noise approximation with additive intrinsic noise (LNAA) and normal measurement error (blue).
%}
%\end{figure}				
%\input{tables/SDE_values}

%\subsection{\label{sec:application_log}Model application to logistic growth with multiplicative error and lognormal random error}
To compare the approximate models ability to recover parameters from the SLGM with simulated lognormal measurement error, we simulate data and carry out Bayesian inference. Figure~\ref{simlog} shows that all three approximate models can capture the synthetic time-course well, but that the RRTR model is the least representative with the largest amount of drift occurring at the saturation stage, a property not found in the SLGM or the two new LNA models.
Comparing forwards trajectories with measurement error (Figure~\ref{simlog}), the ``exact'' model is visually similar to all our approximate models, but least similar to the RRTR.
Further, Table~\ref{app:sde_val_fur} demonstrates that parameter posterior means are close to the true values and that standard deviations are small for all models and each parameter set.
By comparing posterior means and standard deviations to the true values, Table~\ref{app:sde_val_fur} shows that all our models are able to recover the three different parameter sets considered. 
	\begin{figure}[h!]
  \centering
\includegraphics[width=14cm]{img_sde/simlog}
\caption{\label{simlog}
Forward trajectories with measurement error, simulated from parameter posterior samples (sample size=1000). 
Model fitting is carried out on SLGM forward trajectories with lognormal measurement error (black), for three different sets of parameters (see Table~\ref{app:sde_val_fur}).  See (\ref{app:LNAM_sta_spa_mod}) or (\ref{app:LNAA_sta_spa_mod}) for model and Table~\ref{table:SDE_priors} for prior hyper-parameter values.
See Table~\ref{app:sde_val_fur} for parameter posterior means and true values.
A), E) \& I) SLGM+L (orange).
B), F) \& J) RRTR model with lognormal error (red).
C), G) \& K) LNAM model with lognormal error (green).
D), H) \& L) LNAA model with normal error (blue).
}
\end{figure}

%\subsection{\label{sec:application_nor}Model application to logistic growth with multiplicative error and normal random error}

To compare the approximations to the SLGM with simulated normal measurement error, we simulate data and carry out Bayesian inference. 
Figure~\ref{sim} shows that of our approximate models, only the LNAA model can appropriately represent the simulated time-course as both our models with lognormal measurement error, the RRTR and LNAM do not closely bound the data. 
Comparing forwards trajectories with measurement error (Figure~\ref{sim}), the ``exact'' model is most visually similar to the LNAA, which shares the same measurement error structure.
Further, Table~\ref{app:sde_val_fur} demonstrates that only our models with normal measurement error have posterior means close to the true values and that standard deviations are larger in the models with lognormal measurement error. 
Observing the posterior means for $K$ for each parameter set (Table~\ref{app:sde_val_fur}), we can see that the RRTR has the largest standard deviations and that, of the approximate models, its posterior means are furthest from both the true values and the ``exact'' model posterior means. 
Comparing LNA models to the ``exact'' models with matching measurement error, we can see in Table~\ref{app:sde_val_fur} that they share similar posterior means and only slightly larger standard deviations.
Example posterior diagnostics given in Appendix~\ref{app:diag_sim_n}, demonstrate that posteriors are distributed tightly around true values for our LNAA and data from the SLGM with Normal measurement error.

	\begin{figure}[h!]
  \centering
\includegraphics[width=14cm]{img_sde/sim}
\caption{\label{sim}
Forward trajectories with measurement error, simulated from inferred parameter posterior samples (sample size=1000).  
Model fitting is carried out on SLGM forward trajectories with normal measurement error (black), for three different sets of parameters (see Table~\ref{app:sde_val_fur}). 
See (\ref{app:LNAM_sta_spa_mod}) or (\ref{app:LNAA_sta_spa_mod}) for model and Table~\ref{table:SDE_priors} for prior hyper-parameter values.
See Table~\ref{app:sde_val_fur} for parameter posterior means.
A) SLGM+N (pink).
B) RRTR model with lognormal error (red).
C) LNAM model with lognormal error (green).
D) LNAA model with normal error (blue).
}
\end{figure}
\input{tables/SDE_values}

\subsection{\label{sec:application_obs}Application to observed yeast data}
We now consider which diffusion equation model can best represent observed microbial population growth curves taken from a Quantitative Fitness Analysis (QFA) experiment \citep{QFA1,jove}, see Figure~\ref{real}. 
The data consists of scaled cell density estimates over time for budding yeast \emph{Saccharomyces cerevisiae}. 
Independent replicate cultures are inoculated on plates and photographed over a period of 5 days.
The images captured are then converted into estimates of integrated optical density (IOD, which we assume are proportional to cell population size), by the software package Colonyzer \citep{Colonyzer}.
%The growth measures are then scaled to normalise for camera resolution and give a maximum normalised growth of 1.
The dataset chosen for our model fitting is a representative set of 10 time-courses, each with 27 time points.

As in Figure~\ref{sim}, we see that the LNAA model is the only approximation that can appropriately represent the time-course and that both the RRTR and LNAM fail to bound the data as tightly as the LNAA (Figure~\ref{real}).
Our two ``exact'' models are visually similar to our approximate models with the same measurement error, with the SLGM+N most similar to the LNAA and the SLGM+L to the RRTR and LNAM. This is as  expected due to matching measurement error structures.
Table~\ref{app:sde_val_fur} summarises parameter estimates for the observed yeast data using each model. The variation in the the LNAA model parameter posteriors is much smaller than the RRTR and LNAM, indicating a more appropriate model fit. 
Comparing the LNA models and ``exact'' models with matching measurement error, we can see in Table~\ref{app:sde_val_fur} that they share similar posterior means and standard deviations for all parameters and in particular, they are very similar for both $K$ and $r$, which are important phenotypes for calculating fitness \citep{QFA1}.

	\begin{figure}[h!]
  \centering
\includegraphics[width=14cm]{img_sde/real}
\caption{\label{real}
Forward trajectories with measurement error, simulated from inferred parameter posterior samples (sample size=1000).  
Model fitting is carried out on observed yeast time-course data (black). 
See (\ref{app:LNAM_sta_spa_mod}) or (\ref{app:LNAA_sta_spa_mod}) and Table~\ref{table:SDE_priors} for prior hyper-parameter values.
See Table~\ref{app:sde_val_fur} for parameter posterior means.
A) SLGM+N (pink).
B) SLGM+L (orange).
A) RRTR model with lognormal error (red).
B) LNAM model with lognormal error (green).
C) LNAA model with normal error (blue).
}
\end{figure}		

In Table~\ref{tab:MSE}, to compare quality of parameter inference for 10 observed yeast time-courses with each approximate model. Mean square error (MSE) for 1000 posterior sample forward simulations are calculated for each yeast time course and summed to give a Total MSE for each model.
It is clear that the RRTR is the worst overall representation of the 10 yeast time courses, with the highest total MSE and a much larger total MSE than the ``exact'' SLGM+L. It is interesting to see there is a very similar total MSE for the SLGM+L and LNAM, and similarly for the SLGM+N and LNAA, demonstrating that our approximations perform well.
\input{tables/MSE}




conclusion.tex
\clearpage
\section{\label{sec:conclusions}Conclusion}
We have presented two new diffusion processes for modelling logistic growth data where fast inference is required: the linear noise approximation (LNA) of the stochastic logistic growth model (SLGM) with multiplicative noise and the LNA of the SLGM with additive intrinsic noise. 
Both the LNAM and LNAA are derived from the linear noise approximation of the stochastic logistic growth model (SLGM).
The new diffusion processes approximate the SLGM more closely than an alternative approximation (RRTR) proposed by \cite{roman}.
The RRTR lacks a mean reverting property that is found in the SLGM, LNAM and LNAA, resulting in increasing variance during the stationary phase of population growth (see Figure~\ref{4nonu}).

We compared the ability of each of the three approximate models and the SLGM to recover parameter values from simulated datasets using standard MCMC techniques.  
When modelling stochastic logistic growth with lognormal measurement error we find that our approximate models are able to represent data simulated from the original process and that the RRTR is least representative, with large variation over the stationary phase (see Figure~\ref{simlog}).
When modelling stochastic logistic growth with normal measurement error we find that only our models with normal measurement error can appropriately bound data simulated from the original process (see Figure~\ref{sim}).
We also compared parameter posterior distribution summaries with parameter values used to generate simulated data after inference using both approximate and ``exact'' models (see Table~\ref{app:sde_val_fur}).  We find that, when using the RRTR model, posterior distributions for the carrying capacity parameter $K$ are less precise than for the LNAM and LNAA approximations.  We also note that it is not possible to model additive measurement error while maintaining a linear Gaussian structure (which allows fast inference with the Kalman filter) when carrying out inference with the RRTR.  We conclude that when measurement error is additive, the LNAA model is the most appropriate approximate model.  

To test model performance during inference with real population data, we fitted our approximate models and the ``exact'' SLGM to microbial population growth curves generated by quantitative fitness analysis (QFA) (see Figure~\ref{real}).
We found that the LNAA model was the most appropriate for modelling experimental data.  It seems likely that this is because a normal error structure best describes this particular dataset, placing the LNAM and RRTR models at a disadvantage.  We demonstrate that arbitrarily exact methods and our fast approximations perform similarly during inference for 10 diverse, experimentally observed, microbial population growth curves (see Table~\ref{tab:MSE}) which show that, in practise, our fast approximations are as good as ``exact'' methods.
We conclude that our LNA models are preferable to the RRTR for modelling QFA data. 

It is interesting to note that, although the LNAA is not a better approximation of the original SGLM process than the LNAM, it is still quite reasonable. Figure~\ref{4nonu}A~and~\ref{4nonu}D shows that the SLGM and LNAA processes are visually similar. Figure~\ref{4nonu}E demonstrates that forward trajectories of the LNAA also share similar levels of variation over time with the SLGM and LNAM.  

Fast inference with the LNAA gives us the potential to develop large hierarchical Bayesian models which simultaneously describe thousands of independent time-courses from QFA with a diffusion equation, allowing us to infer the existence of genetic interactions on a genome-wide scale using realistic computational resources. 

Here, we have concentrated on a biological model of population growth.  However, we expect that the approach we have demonstrated: generating linear noise approximations of stochastic processes to allow fast Bayesian inference with Kalman filtering for marginal likelihood computation, will be useful in a wide range of other applications where simulation is prohibitively slow.

introduction.tex
\section{\label{sec:Introduction}Introduction}

Stochastic models simultaneously describe dynamics and noise or heterogeneity in real systems \citep{calibayes}.  For example, stochastic models are increasingly recognised as necessary tools for understanding the behaviour of complex biological systems \citep{wilkinson2012stochastic,wilkinson_nature} and are also used to capture uncertainty in financial market behaviour \citep{stochfinance,stochfinance2}.  Many such models are written as continuous stochastic differential equations (SDEs) which often do not have analytical solutions and are slow to evaluate numerically compared to their deterministic counterparts.  Simulation speed is often a particularly critical issue when inferring model parameter values by comparing simulated output with observed data \citep{woodtrees}.  

For SDE models where no explicit expression for the transition density is available, it is possible to infer parameter values by simulating a latent process using a data augmentation approach \citep{darren2005}. However, this method is computationally intensive and not practical for all applications. When fast inference for SDEs is important, for example for real-time analysis as part of decision support systems or for big data inference problems where we simultaneously fit models to many thousands of datasets (e.g. \cite{heydari}), we need an alternative approach.
Here we demonstrate one such approach: developing an analytically tractable approximation to the original SDE, by making linear noise approximations (LNAs).  We apply this approach to a SDE describing logistic population growth for the first time.
  
The logistic model of population growth, an ordinary differential equation (ODE) describing the self-limiting growth of a population of size $X_t$ at time $t$, was developed by \cite{Verhulst1847}
\begin{align}
\label{eq_det}
\frac{dX_t}{dt}&=rX_t-\frac{r}{K}X_t^2.
\end{align}
The ODE has the following analytic solution:
\begin{align}
\label{eq_det_sol}
X_t&=\frac{K}{1+Qe^{-r (t-t_0)}},
\end{align}
where $Q=\left(\frac{K}{P}-1\right)e^{rt_0}$ and $P=X_{t_0}$.
The model describes a population growing from an initial size $P$ with an intrinsic growth rate $r$, undergoing approximately exponential growth which slows as the availability of some critical resource (e.g. nutrients or space) becomes limiting \citep{theoryoflogisticgro}.
Ultimately, population density saturates at the carrying capacity (maximum achievable population density) $K$, once the critical resource is exhausted.
Where further flexibility is required, generalized forms of the logistic growth process \citep{analysisoflogistic,logisticrevisited} may be used instead.

To account for uncertainty about processes affecting population growth which are not explicitly described by the deterministic logistic model, we can include a term describing intrinsic noise and consider a SDE version of the model \citep{optimal_harvesting,roman}.  Here we extend the ODE in (\ref{eq_det}) by adding a term representing multiplicative intrinsic noise (\ref{eq_det_sde}) to give a model which we refer to as the stochastic logistic growth model (SLGM).  
\begin{align}
\label{eq_det_sde}
dX_t&=\left[rX_t-\frac{r}{K}X_t^2\right]dt+\sigma X_t dW_t, 
\end{align}
where $X_{t_0}=P$ and is independent of $W_t$, $t\geq t_0$, 

Alternative stochastic formulations of the logistic ODE can be generated \citep{stolog}.
The Kolmogorov forward equation has not been solved for (\ref{eq_det_sde}) (or for any similar formulation of a logistic SDE) and so no explicit expression for the transition density is available.

\cite{roman} introduce a diffusion process approximating the SLGM (which we label RRTR) with a transition density that can be derived explicitly.  
The Bayesian approach can be applied in a natural way to carry out parameter inference for state space models with tractable transition densities \citep{dynamicmodels}.  
A state space model describes the probabilistic dependence between an observation process variable $X_t$ and state process $S_t$.
The transition density is used to describe the state process $S_t$ and a measurement error structure is chosen to describe the relationship between $X_t$ and $S_t$. 

The Kalman filter \citep{kalmanoriginal} is typically used to infer the hidden state process of interest $S_t$ and is an optimal estimator, minimising the mean square error of estimated parameters when all noise in the system can be assumed to be Gaussian.
We use the Kalman filter to reduce computational time in a parameter inference algorithm by recursively computing the marginal likelihood \citep{dynamicmodels}.

The RRTR can be fit to data within an acceptable time frame by assuming multiplicative measurement error to give a linear Gaussian structure, allowing us to use a Kalman filter for inference. 


We introduce two new first order linear noise approximations (LNAs) \citep{LNA,komorowski} of (\ref{eq_det_sde}), one with multiplicative and one with additive intrinsic noise, which we label LNAM and LNAA respectively.
The LNA reduces a SDE to a linear SDE with additive noise, which can be solved to give an explicit expression for the transition density. 
We derive transition densities for the two approximate models and construct a Kalman filter by choosing measurement noise to be either multiplicative or additive to construct a linear Gaussian structure. 
%We show that a zeroth order LNA of the SDE with multiplicative intrinsic noise is equivalent to the RRTR.
Exact simulations from the SGLM are compared with each of the three approximate models.
We compare the utility of each of the approximate models during parameter inference by comparing simulations with both synthetic and real datasets.
%In order to model both simulated data from (\ref{eq_det_sde}) and time-course data from a QFA experiment, we compare Bayesian models with the RRTR and our new LNA diffusion equations.

LNAA.tex
\section{\label{sec:LNAA} Linear noise approximation with additive noise}
As in Section~\ref{sec:LNAM}, we start from the SLGM, given in (\ref{eq_det_sde}).
Without first log transforming the process, the LNA will lead to a worse approximation to the diffusion term of the SLGM, but we will see in the coming sections that there are nevertheless advantages.
%, in particular we can assume additive normal measurement error in our state space model and use a Kalman filter to compute the marginal likelihood.   %[[NOT VERY CLEAR]]
We separate the process $X_t$ into a deterministic part $V_t$ and a stochastic part $Z_t$ so that $dX_t=dV_t+dZ_t$ and consequently $X_t=V_t+Z_t$.
We choose $V_t$ to be the solution of the deterministic part
\begin{equation*}
dV_t=\left(rV_t-\frac{r}{K}V_t^2\right)dt.
\end{equation*}
After redefining our previous notation as follows: $a=r$ and $b=\frac{r}{K}$, we solve $dV_t$ to give:
\begin{equation*}
V_t=\frac{aPe^{aT}}{bP(e^{aT}-1)+a}.
\end{equation*}
Differentiating w.r.t.~$t$, we obtain: 
\begin{equation*}
\frac{dV_t}{dt}=\frac{a^2Pe^{aT}(a-bP)}{\left(bP(e^{aT}-1)+a\right)^2}.
\end{equation*}
We now solve $dZ_t$, where $dZ_t=dX_t-dV_t$. Expressions for both $dX_t$ and $dV_t$ are known:
\begin{equation*}
dZ_t=\left(rX_t-\frac{r}{K}X_t^2\right)dt+\sigma X_t dW_t-\frac{a^2Pe^{aT}(a-bP)}{\left(bP(e^{aT}-1)+a\right)^2}dt.
\end{equation*}
Simplifying:
\begin{equation*}
dZ_t=\left(aX_t-bX_t^2-\frac{a^2Pe^{aT}(a-bP)}{\left(bP(e^{aT}-1)+a\right)^2}\right)dt+\sigma X_t dW_t.
\end{equation*}
We then substitute in $X_t=V_t+Z_t$ and rearrange to give
\begin{align*}
dZ_t=&\left(aV_t-bV_t^2-\frac{a^2Pe^{aT}(a-bP)}{\left(bP(e^{aT}-1)+a\right)^2}+(a-2bV_t)Z_t-bZ_t^2\right)dt\\
&+\left( \sigma V_t +\sigma Z_t\right) dW_t.
\end{align*}
We now apply the LNA, by setting second-order terms $-bZ^2dt=0$ and $\sigma Z_t dW_t=0$ to obtain
\begin{equation}\label{eq:LNAA_dz}
dZ_t=\left(aV_t-bV_t^2-\frac{a^2Pe^{aT}(a-bP)}{\left(bP(e^{aT}-1)+a\right)^2}+(a-2bV_t)Z_t\right)dt+\sigma V_t  dW_t.
\end{equation}
This process is a particular case of the Ornstein-Uhlenbeck process, which can be solved. 
The transition density for $X_t$ (derivation in Appendix~\ref{app:LNAA_sol}) is then
\begin{align}
\begin{split}\label{eq:LNAA_tran}
(X_{t_i}|X_{t_{i-1}}&=x_{t_{i-1}})\sim N(\mu_{t_i},\Xi_{t_i}),\\
\text{where } x_{t_{i-1}}&=v_{t_{i-1}}+z_{t_{i-1}},\\
\mu_{t_i}&=v_{t_{i-1}}+\left(\frac{aPe^{aT_i}}{bP(e^{aT_i}-1)+a}\right)-\left(\frac{aPe^{aT_{i-1}}}{bP(e^{aT_{i-1}}-1)+a}\right)\\
&+e^{a(t_i-t_{i-1})}\left(\frac{bP(e^{aT_{i-1}}-1)+a}{bP(e^{aT_i}-1)+a}\right)^2Z_{t_{i-1}}\text{ and}\\
 \Xi_t&=\frac{1}{2}\sigma^2aP^2e^{2aT_i}\left(\frac{1}{bP(e^{aT_i}-1)+a}\right)^4\\
&\times[
b^2P^2(e^{2aT_i}-e^{2aT_{i-1}})
+4bP(a-bP)(e^{aT_i}-e^{aT_{i-1}})\\
&\;\:\:\:\:+2a(t_i-t_{i-1})(a-bP)^2
].
\end{split}
\end{align}
The LNA of the SLGM, with additive intrinsic noise (LNAA) can then be written as
\begin{align*}
dX_t=\left[b{V_t}^2+\left(a-2bV_t\right)X_t\right]dt+\sigma V_t dW_t,
\end{align*}
where $P=X_{t_0}$ and is independent of $W_t$, $t\geq t_0$.
%Combining the LNAM with a normal (additive) measurement error model, we can take advantage of the linear Gaussian structure of the log process and construct a Kalman filter for marginal likelihood computation during parameter inference.
LNAM.tex
\section{\label{sec:LNAM} Linear noise approximation with multiplicative noise}
We now take a different approach to approximating the SLGM (\ref{eq_det_sde}), which will turn out to be closer to the exact solution of the SLGM than the RRTR (\ref{eq_sde}). 
Starting from the original model (\ref{eq_det_sde}),
we apply It\^{o}'s lemma with the transformation $f(t,X_t)\equiv Y_t=\log X_t$
% and the following partial derivatives: $\frac{df}{dt}=0,\:\frac{df}{dx}=\frac{1}{X_t},\:\frac{d^2f}{dx^2}=-\frac{1}{X_t^2}$.
%\begin{align}
%\label{eq:ITO}
%df(t,X_t)=\frac{df}{dt}dt+\mu\frac{df}{dx}dt+\frac{1}{2}\sigma^2\frac{d^2f}{dx^2}dt+\sigma\frac{df}{dx}dW_t
%\end{align}
to obtain the following It\^{o} drift-diffusion process:
%\begin{equation*}
%dY_t=\left(rX_t-\frac{r}{K}X_t^2\right)\frac{1}{X_t} dt+\frac{1}{2}\sigma^2 X_t^2\left(-\frac{1}{X_t^2}\right)dt+\sigma X_t\left(\frac{1}{X_t}\right)dW_t
%\end{equation*}
%and simplify to give
\begin{align}
\label{eq:SDE2}
dY_t=\left(r-\frac{1}{2}\sigma^2-\frac{r}{K}e^{Y_t}\right)dt+\sigma dW_t.
\end{align}
The log transformation from multiplicative to additive noise, gives a constant diffusion term, so that the LNA will give a good approximation to (\ref{eq_det_sde}).  
We now separate the process $Y_t$ into a deterministic part $V_t$ and a stochastic part $Z_t$ so that $Y_t=V_t+Z_t$ and consequently $dY_t=dV_t+dZ_t$.
We choose $V_t$ to be the solution of the deterministic part
\begin{equation*}
dV_t=\left(r-\frac{1}{2}\sigma^2-\frac{r}{K}e^{V_t}\right)dt.
\end{equation*}
After redefining our notation as follows: $a=r-\frac{\sigma^2}{2}$ and $b=\frac{r}{K}$, we solve (\ref{eq:SDE2}) for $V_t$:
\begin{equation}\label{eq:LNAM_det_sol}
V_t=\log\left(\frac{aPe^{aT}}{bP(e^{aT}-1)+a}\right).
\end{equation}
Differentiating w.r.t. t, we obtain: 
\begin{equation*} \label{eq:LNAM_dv_t}
\frac{dV_t}{dt}=\frac{a(a-bP)}{bP(e^{aT}-1)+a}.
\end{equation*}
Writing down an expression for $dZ_t$, where $dZ_t=dY_t-dV_t$.
\begin{equation*} 
dZ_t=\left(r-\frac{1}{2}\sigma^2-\frac{r}{K}e^{Y_t}\right)dt+\sigma dW_t-\frac{a(a-bP)}{bP(e^{aT}-1)+a}dt
\end{equation*}
and simplifying:
\begin{equation*}
dZ_t=\left(\frac{baPe^{aT}}{bP(e^{aT}-1)+a}-be^{Y_t}\right)dt+\sigma dW_t.
\end{equation*}
Recognizing that $e^{V_t}=\frac{aPe^{aT}}{bP(e^{aT}-1)+a}$:
\begin{equation*} 
dZ_t=b\left(e^{V_t}-e^{Y_t}\right)dt+\sigma dW_t.
\end{equation*}
We substitute in $Y_t=V_t+Z_t$ to give
\begin{equation*} 
dZ_t=b\left(e^{V_t}-e^{V_t+Z_t}\right)dt+\sigma dW_t.
\end{equation*} 
We now apply the LNA, by making a first-order approximation of $e^{Z_t}\approx 1+Z_t$ and then simplify to give $dZ_t=-be^{V_t}Z_tdt+\sigma dW_t$.
Finally, substitute $e^{V_t}=\frac{aPe^{aT}}{bP(e^{aT}-1)+a}$  to obtain
\begin{equation}\label{eq:LNAM_dz}
dZ_t=-\frac{baPe^{aT}}{bP(e^{aT}-1)+a}Z_tdt+\sigma dW_t.
\end{equation}
This process is a particular case of the time-varying Ornstein-Uhlenbeck process, which can be solved explicitly. %, see Appendix~\ref{app:LNAM_sol} for the solution and transition density of the log process $Y_t$ respectively.\\
The transition density for $Y_t$ (derivation in Appendix~\ref{app:LNAM_sol}) is then:
\begin{align}
\begin{split}\label{eq:LNAM_tran}
(Y_{t_i}|Y_{t_{i-1}}&=y_{t_{i-1}})\sim\operatorname{N}\left(\mu_{t_i},\Xi_{t_i}\right),\\
\text{redefine } y_{t_{i-1}}&=v_{t_{i-1}}+z_{t_{i-1}}, Q=\left(\frac{\frac{a}{b}}{P}-1\right)e^{at_{0}},\\
\mu_{t_i}&=v_{t_{i-1}}+\log\left(\frac{1+Qe^{-at_{i-1}}}{1+Qe^{-at_i}}\right)+e^{-a(t_i-t_{i-1})}\frac{1+Qe^{-at_{i-1}}}{1+Qe^{-at_i}}z_{t_{i-1}} \text{ and}\\
\Xi_{t_i}&=\sigma^2\left[\frac{4Q(e^{at_i}-e^{at_{i-1}})+e^{2at_i}-e^{2at_{i-1}}+2aQ^2(t_i-t_{i-1})}{2a(Q+e^{at_i})^2}\right].
\end{split}
\end{align}
The LNA of the SLGM with multiplicative intrinsic noise (LNAM) can then be written as
\begin{align*}
d\log X_t=\left[dV_t+be^{V_t}V_t-be^{V_t}\log X_t\right]dt+\sigma dW_t,
\end{align*}
where $P=X_{t_0}$ and is independent of $W_t$, $t\geq t_0$.
\\
%Similar to the RRTR, combining the LNAM with a lognormal (multiplicative) measurement error model, we can take advantage of the linear Gaussian structure of the log process and construct a Kalman filter recursion for marginal likelihood computation during parameter inference.
Note that the RRTR given in (\ref{eq_sde}) can be similarly derived using a zero-order noise approximation ($e^{Z_t}\approx 1$) instead of the LNA.

roman.tex
\section{\label{sec:Roman}The \cite{roman} diffusion process}

\cite{roman} present a logistic growth diffusion process (RRTR) which has a transition density that can be written explicitly, allowing inference of model parameter values from discrete sampling trajectories.
 \\
The RRTR is derived from the following ODE:
\begin{align}
\label{eq_ode}
\frac{dx_t}{dt}&=\frac{Qr}{e^{rt}+Q}x_t.
\end{align}
The solution to (\ref{eq_ode}) is given in (\ref{eq_det_sol}) (it has the same solution as (\ref{eq_det})).

\cite{roman} see (\ref{eq_ode}) as a generalisation of the Malthusian growth model with a deterministic, time-dependent fertility $h(t)=\frac{Qr}{e^{rt}+Q}$, and replace this with $\frac{Qr}{e^{rt}+Q}+\sigma W_t$ to obtain the following approximation to the SLGM:
\begin{align}
\label{eq_sde}
dX_t&=\frac{Qr}{e^{r{t}}+Q}X_td{t}+{\sigma}X_tdW_t,
\end{align}
where  $Q=\left(\frac{K}{P}-1\right)e^{rt_0}$, $P=X_{t_0}$ and is independent of $W_t$, $t\geq t_0$.
%Equation~\ref{eq_sde} is refereed to as a new logistic diffusion process by \cite{roman}. 
The process described in (\ref{eq_sde}) is a particular case of the lognomal process with exogenous factors, therefore an exact transition density is available \citep{gutierrez}.
The transition density for $Y_t$, where $Y_t=\log(X_t)$, can be written:
\begin{align}
\begin{split}\label{eq:RRTR_tran}
(Y_{t_i}|Y_{t_{i-1}} &= y_{t_{i-1}}) \sim\operatorname{N}\left(
\mu_{t_i}
,
\Xi_{t_i}
\right),\\
\text{where } a &= r,\quad b=\frac{r}{K}, \\
\mu_{t_i}=& \log(y_{t_{i-1}})
+\log(\frac{1+be^{-at_i}}{1+be^{-at_{i-1}}})
-\frac{\sigma^2}{2}(t_i-t_{i-1}) \text{ and}\\
\Xi_{t_i} &= \sigma^2(t_i-t_{i-1}).
\end{split}
\end{align}
%Combining the RRTR with a lognormal (multiplicative) measurement error model, allows us to take advantage of the linear Gaussian structure of the log process and construct a Kalman filter for marginal likelihood computation during parameter inference.

MSE.tex
\begin{table}[h!]
\caption{\label{tab:MSE} Total mean squared error (MSE) for 10 observed yeast growth time courses, each with 1000 forward simulated time-courses with measurement error. Parameter values are taken from posterior samples. Standard Deviations give the variation between the sub-total MSEs for each yeast time course fit (n=10). }
\centering
\npdecimalsign{.}
\nprounddigits{3}
\begin{tabular}{c n{1}{3} n{1}{3} n{1}{3} n{1}{3} n{1}{3}}
\hline
\emph{Model} & \multicolumn{1}{c}{{SLGM+N}}  & \multicolumn{1}{c}{{SLGM+L}}  &  \multicolumn{1}{c}{{RRTR}} & \multicolumn{1}{c}{{LNAM}} & \multicolumn{1}{c}{{LNAA}} \\ 
\hline
\noalign{\vskip 0.4mm} 
Total MSE & 29.84698 & 100.165 & 600.6007 & 99.39689 &  30.95872 \\
Standard Deviation & 1.688649 & 8.390729 &55.72009 &9.263094& 2.029686\\
\noalign{\vskip 0.4mm} 
\hline
\end{tabular}
\npnoround%
\end{table}

\begin{comment}
 SLGM
 [1] 1.593582 4.618251 1.046460 1.385996 3.576298 5.240944 4.299443 1.637226
 [9] 1.543953 4.904824
> SLGML
 [1]  5.688341 24.274925  1.330012  2.993272 11.128968 10.813295  7.254462
 [8]  3.805309  7.598946 25.277488
> ROMAN
 [1] 187.87644  40.51193   5.96081  67.65204  25.62293  26.38231  64.22651
 [8] 119.86702  13.01995  49.48071
> LNALOG
 [1]  4.865441 24.883842  1.379867  2.411647 11.078207 10.253571  6.489727
 [8]  3.138613  6.937851 27.958122
> LNANORMAL
 [1] 1.6094509 5.3628841 0.9569713 1.2687224 3.6645595 6.1215305 3.2885825
 [8] 1.5801627 1.3638529 5.7420076
> sum(SLGM)
[1] 29.84698
> sd(SLGM)
[1] 1.688649
> sum(SLGML)
[1] 100.165
> sd(SLGML)
[1] 8.390729
> sum(ROMAN)
[1] 600.6007
> sd(ROMAN)
[1] 55.72009
> sum(LNALOG)
[1] 99.39689
> sd(LNALOG)
[1] 9.263094
> sum(LNANORMAL)
[1] 30.95872
> sd(LNANORMAL)
[1] 2.029686
>
\end{comment}



SDE_priors.tex
\begin{table}[h!]
\caption{Prior hyper-parameters for Bayesian sate space models, Lognormal with mean ($\mu$) and precision ($\tau$) \label{table:SDE_priors}}
\centering     
\begin{tabular}{c c}
    \hline
		\noalign{\vskip 0.4mm} 
    Parameter Name  & \multicolumn{1}{c}{Value} \\ \hline
 ${\mu}_K$ & $\log(0.1)$   \\ 
 $\tau_{K}$ & 2 \\ 
 ${\mu}_r$ & $\log(3)$  \\ 
 $\tau_{r}$ & 5\\ 
 ${\mu}_P$ &  $\log(0.0001)$\\ 
 $\tau_{P}$ & 0.1 \\ 
 ${\mu}_\sigma$ & $\log(100)$  \\ 
 $\tau_{\sigma}$ & 0.1 \\ 
 ${\mu}_\nu$ & $\log(10000)$  \\ 
 $\tau_{\nu}$ & 0.1 \\ 
    \hline
    \end{tabular}
    \npnoround
\end{table}
\clearpage






SDE_values.tex
\begin{table}[h!]
\caption{\label{app:sde_val_fur}Bayesian state space model parameter posterior means, standard deviations and true values for Figure \ref{simlog}, \ref{sim} and \ref{real}. True values for the simulated data used for Figure \ref{4nonu}, \ref{simlog} and \ref{sim} are also given.}
\centering
\resizebox{\columnwidth}{!}{%
\npdecimalsign{.}
\nprounddigits{3}
\begin{tabular}{c c n{2}{5} n{1}{3} n{1}{3} n{1}{3} n{1}{3} n{1}{3} n{1}{3} n{1}{3} n{1}{3} n{1}{3}}
\hline
\emph{Panel} & \emph{Model} & \multicolumn{2}{c}{\emph{$\hat{K}$}}  & \multicolumn{2}{c}{\emph{$\hat{r}$}}  & \multicolumn{2}{c}{\emph{$\hat{P}$}} & \multicolumn{2}{c}{\emph{$\hat{\nu}$}} & \multicolumn{2}{c}{\emph{$\hat{\zeta}$}}  \\ 
\hline
\noalign{\vskip 0.4mm} 
\multicolumn{4}{l}{\emph{Figure \ref{simlog}, SLGM with lognormal error}}&&&&&&&&\\
\noalign{\vskip 0.2mm} 
A & SGLM+L & 0.1495369660 &( 0.001)&
2.9822427445 &(0.01350133)&
1.001567e-04 &(1.111598e-06)&
3.8597809e-03 &(2.126679e-03)&
0.0173228498 &(0.005037563)\\ 
B & RRTR & 0.149999335 &(0.00341)&
2.99015285 &(0.01092)&
9.93119e-05 &(1.06863e-06)&
5.6843212256e-03 &(2.36e-03)&
0.011547265557 &(0.00618)\\ 
C & LNAM & 0.149614787 &(0.001)&
2.98773783 &(0.01318)&
9.97960e-05 &(1.124219e-06)&
4.1396496436e-03 &(2.179554e-03)&
0.016364920984 &(0.00517)\\
D & LNAA & 0.14953488 &(0.0005)&
3.00514737 &(0.01995)&
9.64713e-05 &(2.945503e-06)&
3.098832e-05 &(2.53443e-05)&
0.0194997429 &(0.00344)\\ 
E & SGLM+L & 0.1096583 &(0.0007301528)&
3.974696 &(0.04724553)&
5.054189e-05 &(1.568196e-06)&
6.158735e-03 &(5.527275e-03)&
0.05052818 &(0.01393793)\\ 
F  & RRTR & 0.1086907   &(0.006533453)&
3.983695  &(0.03462259)&
5.046032e-05  &(1.137058e-06)&
5.928172e-03  &(4.596061e-03)&
0.03704644 &(0.00901516)\\ 
G  & LNAM & 0.1097371  &(0.0007490187)&
 3.985208  &(0.04606133)&
 5.043313e-05&(1.580017e-06)&
 6.187716e-03 &(5.190529e-03)&
  0.05156039 &(0.0126579)\\
H  & LNAA  &  0.1097912  &(0.0007867848)&
 3.959346  &(0.06704208)&
 5.207434e-05  &(4.309876e-06)&
 4.539877e-05   &(4.394698e-05)&
 0.05912248 &(0.00978808)\\ 
I
& SGLM+L & 0.3003784 &(0.0009748994)&
  5.996904  &(0.02851697)&
 1.962434e-05  &(4.0414e-07)&
 9.543053e-03 &(4.03492e-03)&
0.02406215 &(0.01543261)\\ 
J & RRTR & 0.3005997   &(0.0044928)&
6.015265  &(0.01697153)&
 1.942819e-05 &(2.835413e-07)&
 1.240917e-02 &(2.306528e-03)&
 0.007745952  &(0.006357571)\\ 
K & LNAM &  0.3004755 &(0.0009411441)&
  6.014510 &(0.03065974)&
1.953061e-05  &(4.202024e-07)&
8.942722e-03  &(4.251644e-03)&
0.02690633 &(0.01575695)\\
L & LNAA  &  0.3004586  &(0.001146215)&
 6.037006   &( 0.06747846)&
 1.895225e-05 &(1.502188e-06)&
  8.122229e-05  &(1.595971e-04)&
 0.04733552  &(0.007512705)\\ 
\noalign{\vskip 1mm} 
\multicolumn{4}{l}{\emph{Figure \ref{sim}, SLGM with normal error}}&&&&&&&&\\
\noalign{\vskip 0.2mm} 
A
& SLGM+N & 0.1502746 &(0.002206958)&
3.099127 &(0.08534504)&
9.298825e-05 &(7.304899e-06)&
5.326174e-03 &(1.008844e-03)&
0.05948689 &(0.02986999)\\ 
B & RRTR & 0.2127014398 &(0.12323)&
1.367943557 &(0.26339)&
4.55242737e-03 &(2.1184933e-03)&
2.5392518626e-01 &(1.0969e-01)&
0.41885683478 &(0.12943)\\ 
C  & LNAM & 0.1714519199 &(0.03280)&
1.579950397 &(0.27107)&
5.2409643e-03 &(2.0479954e-03)&
2.05421122203e-01 &(7.805e-02)&
0.4730700594  &(0.05100)\\ 
D & LNAA & 0.150452393 &(0.00222)&
2.9899791 &(0.26239)&
1.189e-04 &(7.098716e-05)&
5.49023564e-03 &(1.06e-03)&
0.052829431599  &(0.03326)\\
E
& SLGM+N &  0.1094272  &(0.0007586429)&
4.183277 &(0.07356718)&
 4.389757e-05 &(4.128667e-06)&
9.679088e-04  &(2.806267e-04)&
0.05662970 &(0.01165919)\\ 
F & RRTR & 0.1574587528 &(0.08749706)&
2.6313040010 &(0.3367076)&
4.398197e-04 &(1.677808e-04)&
1.040098786e-01 &(1.008875e-01)&
0.3735663206 &(0.1616813)\\ 
G  & LNAM & 0.1159608439 &(0.009422523)&
3.0185447100 &(0.373545)&
 4.967272e-04 &(1.396969e-04)&
3.34648648e-02 &(4.308955e-02)&
0.4753192595 &(0.04357369)\\ 
H & LNAA &  0.1095810  &(0.0007913119)&
4.009893 &( 0.1581109)&
5.011655e-05 &( 1.442577e-05)&
1.093103e-03 &(3.638245e-04)&
0.05275491 &(0.01322628)\\ 
I
& SLGM+N & 0.3052235554 &(0.00282216)&
5.2670916442 &(0.1249732)&
3.263387e-04 &(3.407446e-05)&
1.11938047e-02 &(1.974384e-03)&
0.0446024422 &(0.0310588)\\ 
J & RRTR & 0.313643513 &(0.05677519)&
3.029806490 &(0.2325609)&
1.307166e-03 &(2.897121e-04)&
2.22843445e-01 &(3.707941e-02)&
0.074919987 &(0.08628192)\\ 
K  & LNAM  & 0.312629722 &(0.02045301)&
3.391999960 &(0.4300475)&
1.118432e-03 &(3.2685e-04)&
1.17641652e-01 &(8.435387e-02)&
0.360004398 &(0.1647072)\\ 
L & LNAA &  0.3022532 &(0.002361439)&
 5.862218 &(0.5228873)&
2.889783e-05 &(2.598836e-05)&
8.773735e-03 &(1.466338e-03)&
0.04084204  &(0.02844386)\\
\noalign{\vskip 1mm} 
\multicolumn{4}{l}{\emph{ Figure \ref{real}, observed yeast data}}&&&&&&&&\\
\noalign{\vskip 0.2mm} 
A & SLGM+L &  0.1096916 &(0.007456264)&
4.098380 &(0.2994734)&
7.603410e-06 &(3.205816e-06)&
3.456518e-01 &(5.318567e-02)&
0.1128021 &(0.109331)\\ 
B & SLGM+N &  0.1098586 &(0.003273399)&
3.905401 &(0.1725181)&
1.043704e-05 &(3.085561e-06)&
1.851678e-04 &( 7.460388e-05)&
0.1666784 &(0.02814871)\\ 
C & RRTR &  0.1143094 &(0.02602507)&
3.763728 &(0.2006525)&
1.079002e-05 &(3.154911e-06)&
3.378525e-01 &(4.839998e-02)&
0.07771955 &(0.07732242)\\ 
D & LNAM & 0.1103702 &(0.01113706)&
3.776766 &(0.216272)&
1.077129e-05 &(3.277218e-06)&
3.36179e-01 &(5.137415e-02)&
0.1036413 &(0.1079923)\\
E & LNAA & 0.109176 &(0.003306315)&
3.832318 &(0.1984635)&
1.068951e-05 &(3.680477e-06)&
1.768693e-04 &(6.606697e-05)&
0.1637506 &(0.03259627)\\
\\ 
\hline
\multicolumn{2}{c}{\emph{True values}} & \multicolumn{2}{c}{\emph{K}}  & \multicolumn{2}{c}{\emph{r}}  & \multicolumn{2}{c}{\emph{P}} & \multicolumn{2}{c}{\emph{$\nu$}} & \multicolumn{2}{c}{\emph{$\sigma$}}   \\
\hline
\multicolumn{2}{c}{Figure~\ref{4nonu}, panels A, B, C and D} & \multicolumn{2}{c}{0.11} & \multicolumn{2}{c}{4} & \multicolumn{2}{c}{0.00005} & \multicolumn{2}{c}{N/A} & \multicolumn{2}{c}{0.05}  \\
\multicolumn{2}{c}{Figure~\ref{simlog} and \ref{sim}, panels A, B, C \& D}  & \multicolumn{2}{c}{0.15} & \multicolumn{2}{c}{3} & \multicolumn{2}{c}{0.0001} & \multicolumn{2}{c}{0.005} & \multicolumn{2}{c}{0.01} \\ 
\multicolumn{2}{c}{Figure~\ref{simlog} and \ref{sim}, panels E, F, G and H} & 
\multicolumn{2}{c}{0.11} & \multicolumn{2}{c}{4} & \multicolumn{2}{c}{0.00005} & \multicolumn{2}{c}{0.001} & \multicolumn{2}{c}{0.05}  \\
\multicolumn{2}{c}{Figure~\ref{simlog} and \ref{sim}, panels I, J, K and L} &
 \multicolumn{2}{c}{0.3} & \multicolumn{2}{c}{6} & \multicolumn{2}{c}{0.0002} & \multicolumn{2}{c}{0.01} & \multicolumn{2}{c}{0.02}  \\
\hline
\end{tabular}
\npnoround%
}
\end{table}




\begin{comment}
samp<-read.table("test.txt",header=F)
colMeans(samp)
sd(samp[,1])
sd(samp[,2])
sd(samp[,3])
sd(samp[,4])
sd(samp[,5])


N
> colMeans(samp)
           V1            V2            V3            V4            V5
 1.095810e-01  4.009893e+00  5.011655e-05  1.093103e-03  5.275491e-02
           V6
-9.939633e+00
> sd(samp[,1])
[1] 0.0007913119
> sd(samp[,2])
[1] 0.1581109
> sd(samp[,3])
[1] 1.442577e-05
> sd(samp[,4])
[1] 0.0003638245
> sd(samp[,5])
[1] 0.01322628

L
> samp<-read.table("test.txt",header=F)
> colMeans(samp)
          V1           V2           V3           V4           V5
0.1159608439 3.0185447100 0.0004967272 0.0334648648 0.4753192595
> sd(samp[,1])
[1] 0.009422523
> sd(samp[,2])
[1] 0.373545
> sd(samp[,3])
[1] 0.0001396969
> sd(samp[,4])
[1] 0.04308955
> sd(samp[,5])
[1] 0.04357369
>

R
> samp<-read.table("test.txt",header=F)
> colMeans(samp)
          V1           V2           V3           V4           V5
0.1574587528 2.6313000000 0.0004398197 0.1040098786 0.3735663206
> sd(samp[,1])
[1] 0.08749706
> sd(samp[,2])
[1] 0.3367076
> sd(samp[,3])
[1] 0.0001677808
> sd(samp[,4])
[1] 0.1008875
> sd(samp[,5])
[1] 0.1616813


S
Means(sam> colMeans(samp)
          V1           V2           V3           V4           V5
1.094272e-01 4.183277e+00 4.389757e-05 9.679088e-04 5.662970e-02
> sd(samp[,1])
[1] 0.0007586429
> sd(samp[,2])
[1] 0.07356718
> sd(samp[,3])
[1] 4.128667e-06
> sd(samp[,4])
[1] 0.0002806267
> sd(samp[,5])
[1] 0.01165919





sim3
n


> samp<-read.table("test.txt",header=F)
(> colMeans(samp)
           V1            V2            V3            V4            V5
 3.022532e-01  5.862218e+00  2.889783e-05  8.773735e-03  4.084204e-02
           V6
-1.079402e+01
> sd(samp[,1])
[1] 0.002361439
> sd(samp[,2])
[1] 0.5228873
> sd(samp[,3])
[1] 2.598836e-05
> sd(samp[,4])
[1] 0.001466338
> sd(samp[,5])
[1] 0.02844386
L

> samp<-read.table("test.txt",header=F)
> colMeans(samp)
         V1          V2          V3          V4          V5
0.312629722 3.391999960 0.001118432 0.117641652 0.360004398
> sd(samp[,1])
[1] 0.02045301
> sd(samp[,2])
[1] 0.4300475
> sd(samp[,3])
[1] 0.00032685
> sd(samp[,4])
[1] 0.08435387
> sd(samp[,5])
[1] 0.1647072

R

> samp<-read.table("test.txt",header=F)
> colMeans(samp)
         V1          V2          V3          V4          V5
0.313643513 3.029806490 0.001307166 0.222843445 0.074919987
> sd(samp[,1])
[1] 0.05677519
> sd(samp[,2])
[1] 0.2325609
> sd(samp[,3])
[1] 0.0002897121
> sd(samp[,4])
[1] 0.03707941
> sd(samp[,5])
[1] 0.08628192
>


S
eans(samp> colMeans(samp)
          V1           V2           V3           V4           V5
0.3052235554 4.2800916442 0.0003263387 0.0111938047 0.0446024422
> sd(samp[,1])
[1] 0.00282216
> sd(samp[,2])
[1] 0.1249732
> sd(samp[,3])
[1] 3.407446e-05
> sd(samp[,4])
[1] 0.001974384
> sd(samp[,5])
[1] 0.0310588


\end{comment}


paper.bbl
\begin{thebibliography}{}

\bibitem[\protect\citeauthoryear{Addinall, Holstein, Lawless, Yu, Chapman,
  Banks, Ngo, Maringele, Taschuk, Young, Ciesiolka, Lister, Wipat, Wilkinson,
  and Lydall}{Addinall et~al.}{2011}]{QFA1}
Addinall, S.~G., E.-M. Holstein, C.~Lawless, M.~Yu, K.~Chapman, A.~P. Banks,
  H.-P. Ngo, L.~Maringele, M.~Taschuk, A.~Young, A.~Ciesiolka, A.~L. Lister,
  A.~Wipat, D.~J. Wilkinson, and D.~Lydall (2011).
\newblock {Quantitative Fitness Analysis Shows That NMD Proteins and Many Other
  Protein Complexes Suppress or Enhance Distinct Telomere Cap Defects}.
\newblock {\em PLoS Genet\/}~{\em 7\/}(4), e1001362.

\bibitem[\protect\citeauthoryear{Banks, Lawless, and Lydall}{Banks
  et~al.}{2012}]{jove}
Banks, A., C.~Lawless, and D.~Lydall (2012).
\newblock {A Quantitative Fitness Analysis Workflow}.
\newblock {\em J. Vis. Exp\/}~{\em 66}, e4018.

\bibitem[\protect\citeauthoryear{Campillo, Joannides, and
  Larramendy-Valverde}{Campillo et~al.}{2013}]{stolog}
Campillo, F., M.~Joannides, and I.~Larramendy-Valverde (2013).
\newblock {Estimation of the parameters of a stochastic logistic growth model}.
\newblock In submission, (available
  at,\:\:\url{http://arxiv.org/abs/1307.2217}).

\bibitem[\protect\citeauthoryear{Carletti}{Carletti}{2006}]{embook}
Carletti, M. (2006).
\newblock Numerical solution of stochastic differential problems in the
  biosciences.
\newblock {\em J. Comput. Appl. Math.\/}~{\em 185\/}(2), 422--440.

\bibitem[\protect\citeauthoryear{Chen, Lawless, Gillespie, Wu, Boys, and
  Wilkinson}{Chen et~al.}{2010}]{calibayes}
Chen, Y., C.~Lawless, C.~S. Gillespie, J.~Wu, R.~J. Boys, and D.~J. Wilkinson
  (2010).
\newblock Calibayes and basis: integrated tools for the calibration, simulation
  and storage of biological simulation models.
\newblock {\em Briefings in Bioinformatics\/}~{\em 11\/}(3), 278--289.

\bibitem[\protect\citeauthoryear{Gamerman and Lopes}{Gamerman and
  Lopes}{2006}]{gamerman}
Gamerman, D. and H.~Lopes (2006).
\newblock {\em Markov Chain Monte Carlo: Stochastic Simulation for Bayesian
  Inference}.
\newblock Chapman and Hall/CRC Texts in Statistical Science Series. Taylor \&
  Francis.

\bibitem[\protect\citeauthoryear{Golightly and Wilkinson}{Golightly and
  Wilkinson}{2005}]{darren2005}
Golightly, A. and D.~J. Wilkinson (2005).
\newblock Bayesian inference for stochastic kinetic models using a diffusion
  approximation.
\newblock {\em Biometrics\/}~{\em 61\/}(3), 781--788.

\bibitem[\protect\citeauthoryear{Guti\'{e}rrez, Rico, Rom\'{a}n-Rom\'{a}n, and
  Torres-Ruiz}{Guti\'{e}rrez et~al.}{2006}]{gutierrez}
Guti\'{e}rrez, R., N.~Rico, P.~Rom\'{a}n-Rom\'{a}n, and F.~Torres-Ruiz (2006).
\newblock {Approximate and generalized confidence bands for some parametric
  functions of the lognormal diffusion process with exogenous factors}.
\newblock {\em Scientiae Mathematicae Japonicae\/}~{\em 64\/}(2), 313--330.

\bibitem[\protect\citeauthoryear{Heidelberger and Welch}{Heidelberger and
  Welch}{1981}]{Heidelberger}
Heidelberger, P. and P.~D. Welch (1981).
\newblock A spectral method for confidence interval generation and run length
  control in simulations.
\newblock {\em Commun. ACM\/}~{\em 24\/}(4), 233--245.

\bibitem[\protect\citeauthoryear{Heydari, Lawless, Lydall, and
  Wilkinson}{Heydari et~al.}{2012}]{heydari}
Heydari, J.~J., C.~Lawless, D.~A. Lydall, and D.~J. Wilkinson (2012).
\newblock Bayesian hierarchical modelling for inferring genetic interactions in
  yeast.
\newblock {\em {In submission}\/}.

\bibitem[\protect\citeauthoryear{Hurn, Jeisman, and Lindsay}{Hurn
  et~al.}{2007}]{woodtrees}
Hurn, A.~S., J.~I. Jeisman, and K.~A. Lindsay (2007).
\newblock {Seeing the Wood for the Trees: A Critical Evaluation of Methods to
  Estimate the Parameters of Stochastic Differential Equations}.
\newblock {\em Journal of Financial Econometrics\/}~{\em 5\/}(3), 390--455.

\bibitem[\protect\citeauthoryear{Kalman}{Kalman}{1960}]{kalmanoriginal}
Kalman, R.~E. (1960).
\newblock A new approach to linear filtering and prediction problems.
\newblock {\em Transactions of the ASME--Journal of Basic Engineering\/}~{\em
  82\/}(Series D), 35--45.

\bibitem[\protect\citeauthoryear{Kijima}{Kijima}{2013}]{stochfinance}
Kijima, M. (2013).
\newblock {\em Stochastic Processes with Applications to Finance, Second
  Edition}.
\newblock Chapman \& Hall/CRC financial mathematics series. CRC Press.

\bibitem[\protect\citeauthoryear{Kloeden and Platen}{Kloeden and
  Platen}{1992}]{euler_maruyama}
Kloeden, P. and E.~Platen (1992).
\newblock {\em Numerical Solution of Stochastic Differential Equations}.
\newblock Applications of mathematics : stochastic modelling and applied
  probability. Springer.

\bibitem[\protect\citeauthoryear{Koller}{Koller}{2012}]{stochfinance2}
Koller, M. (2012).
\newblock {\em Stochastic Models in Life Insurance}.
\newblock EAA Series. Springer.

\bibitem[\protect\citeauthoryear{Komorowski, Finkenstadt, Harper, and
  Rand}{Komorowski et~al.}{2009}]{komorowski}
Komorowski, M., B.~Finkenstadt, C.~V. Harper, and D.~A. Rand (2009).
\newblock {Bayesian inference of biochemical kinetic parameters using the
  linear noise approximation}.
\newblock {\em BMC Bioinformatics\/}~{\em 10}, 343.

\bibitem[\protect\citeauthoryear{Lawless, Wilkinson, Young, Addinall, and
  Lydall}{Lawless et~al.}{2010}]{Colonyzer}
Lawless, C., D.~J. Wilkinson, A.~Young, S.~G. Addinall, and D.~A. Lydall
  (2010).
\newblock {Colonyzer: automated quantification of micro-organism growth
  characteristics on solid agar}.
\newblock {\em BMC Bioinformatics\/}~{\em 11}, 287.

\bibitem[\protect\citeauthoryear{Li, Wang, and Su}{Li
  et~al.}{2011}]{optimal_harvesting}
Li, W., K.~Wang, and H.~Su (2011).
\newblock Optimal harvesting policy for stochastic logistic population model.
\newblock {\em Applied Mathematics and Computation\/}~{\em 218\/}(1), 157 --
  162.

\bibitem[\protect\citeauthoryear{Peleg, Corradini, and Normand}{Peleg
  et~al.}{2007}]{logisticrevisited}
Peleg, M., M.~G. Corradini, and M.~D. Normand (2007).
\newblock {The logistic (Verhulst) model for sigmoid microbial growth curves
  revisited}.
\newblock {\em Food Research International\/}~{\em 40\/}(7), 808 -- 818.

\bibitem[\protect\citeauthoryear{Plummer}{Plummer}{2010}]{rjags}
Plummer, M. (2010).
\newblock {\em {rjags: Bayesian graphical models using MCMC}}.
\newblock R package version 2.1.0-10.

\bibitem[\protect\citeauthoryear{Plummer, Best, Cowles, and Vines}{Plummer
  et~al.}{2006}]{coda}
Plummer, M., N.~Best, K.~Cowles, and K.~Vines (2006).
\newblock {CODA: Convergence Diagnosis and Output Analysis for MCMC}.
\newblock {\em R News\/}~{\em 6\/}(1), 7--11.

\bibitem[\protect\citeauthoryear{Rom\'{a}n-Rom\'{a}n and
  Torres-Ruiz}{Rom\'{a}n-Rom\'{a}n and Torres-Ruiz}{2012}]{roman}
Rom\'{a}n-Rom\'{a}n, P. and F.~Torres-Ruiz (2012).
\newblock {Modelling logistic growth by a new diffusion process: Application to
  biological systems}.
\newblock {\em Biosystems\/}~{\em 110\/}(1), 9--21.

\bibitem[\protect\citeauthoryear{Tsoularis and Wallace}{Tsoularis and
  Wallace}{2002}]{analysisoflogistic}
Tsoularis, A. and J.~Wallace (2002).
\newblock Analysis of logistic growth models.
\newblock {\em Mathematical Biosciences\/}~{\em 179\/}(1), 21 -- 55.

\bibitem[\protect\citeauthoryear{Turner~Jr., Bradley~Jr., Kirk, and
  Pruitt}{Turner~Jr. et~al.}{1976}]{theoryoflogisticgro}
Turner~Jr., M.~E., E.~L. Bradley~Jr., K.~A. Kirk, and K.~M. Pruitt (1976).
\newblock A theory of growth.
\newblock {\em Mathematical Biosciences\/}~{\em 29\/}(3–4), 367 -- 373.

\bibitem[\protect\citeauthoryear{Verhulst}{Verhulst}{1845}]{Verhulst1847}
Verhulst, P.~F. (1845).
\newblock {Recherches math\'{e}matiques sur la loi d'accroissement de la
  population.}
\newblock {\em Nouveaux m\'{e}moires de l'Academie Royale des Science et
  Belles-Lettres de Bruxelles\/}~{\em 18}, 1--41.

\bibitem[\protect\citeauthoryear{Wallace}{Wallace}{2010}]{LNA}
Wallace, E. W.~J. (2010).
\newblock {A simplified derivation of the Linear Noise Approximation}.
\newblock In submission,( available
  at,\:\:\url{http://arxiv.org/abs/1004.4280}).

\bibitem[\protect\citeauthoryear{West and Harrison}{West and
  Harrison}{1997}]{dynamicmodels}
West, M. and J.~Harrison (1997).
\newblock {\em {Bayesian Forecasting and Dynamic Models}\/} (Second ed.).
\newblock Springer Series in Statistics. New York: Springer-Verlag.

\bibitem[\protect\citeauthoryear{Wilkinson}{Wilkinson}{2009}]{wilkinson_nature}
Wilkinson, D.~J. (2009).
\newblock {Stochastic modelling for quantitative description of heterogeneous
  biological systems}.
\newblock {\em Nature Reviews Genetics\/}~{\em 10\/}(2), 122--133.

\bibitem[\protect\citeauthoryear{Wilkinson}{Wilkinson}{2012}]{wilkinson2012stochastic}
Wilkinson, D.~J. (2012).
\newblock {\em Stochastic Modelling for Systems Biology}.
\newblock Chapman \& Hall / CRC Mathematical and computational biology series.
  CRC Press.

\end{thebibliography}
