%Poster for ICSB, October 2005
% template for latex A0 posters
\documentclass[portrait]{a0poster}
\usepackage{multicol,graphicx,times}%,mathptm}
\usepackage{pstricks,amsmath,amsfonts,url}
\usepackage{verbatim,psfrag}
\usepackage{wrapfig}
\usepackage{prepbib}

% my colours
\newcmykcolor{myyellow}{0 0 0.2 0}
\newrgbcolor{mygreen}{0 0.5 0}
\newrgbcolor{myred}{0.9 0 0.1}
\newrgbcolor{myblue}{0 0.1 0.9}
\newrgbcolor{mymagenta}{0.9 0.0 0.9}

% set number of columns here
\newcommand{\numcols}{3}
\newcommand{\oonumcols}{0.333} % one over number of columns
%\newcommand{\oonumcols}{0.3} % one over number of columns

%added by eryk 2005
\newcommand{\cb}{\textsc{CaliBayes}}
\newcommand{\note}[1]{\framebox{\texttt{#1}}}
\newcommand{\pnote}[1]{\par\framebox{\parbox{15cm}{\texttt{#1}}}\par}

\renewcommand\ref{\hangindent=5mm\hangafter=1\noindent}

\renewcommand\refname{}

% DAH had to add this to make it go A4
\setlength{\textwidth}{0.95\textwidth}

% define a figure environment
\newlength{\colwidth}
\setlength{\colwidth}{\oonumcols\textwidth}
\newlength{\figwidth}
\setlength{\figwidth}{0.85\colwidth}
\newcounter{figscount}
\setcounter{figscount}{0}
\newcommand{\figcap}{blah}
\newenvironment{myfig}[1]%
{\begin{minipage}{\figwidth}\addtocounter{figscount}{1}\vspace*{1cm}%
\renewcommand{\figcap}{#1}}%
{\vspace*{0.2cm}\newline{\small\mygreen Figure~\arabic{figscount}: %
\figcap}\end{minipage}}

% alter sectioning commands
\makeatletter
\newcommand{\mysection}{\@startsection{section}{1}{0mm}{-\baselineskip}{0.5\baselineskip}{\sf\bfseries\Large}}
\newcommand{\mysubsection}{\@startsection{subsection}{2}{0mm}{-\baselineskip}{0.5\baselineskip}{\sf\bfseries\large}}
\makeatother
\newcommand{\sectionpost}[1]{\mysection{{\myblue #1}}}
%\renewcommand{\section}[1]{\mysection{{\myblue #1}}}
\renewcommand{\subsection}[1]{\mysubsection{{\myred #1}}}
\renewcommand{\abstractname}{\mygreen Summary}
\newcommand{\sectionnonumber}[1]{\myblue #1}

\renewcommand{\vec}[1]{\boldsymbol{#1}}

\newcommand{\display}[1]{\par\centerline{
\psframebox[linewidth=2mm,linecolor=myred]{\psframebox*[fillcolor=myyellow]{
\myblue #1}}}
}



\raggedright

% now for the poster itself...
\begin{document}

% Begin code to draw border round poster 
\psset{linewidth=0.3cm}
\newlength{\frameleft}
\newlength{\frameright}
\newlength{\frametop}
\newlength{\framebottom}
\setlength{\frameleft}{-1.8cm}
\setlength{\frameright}{\textwidth}
%\addtolength{\frameright}{-\frameleft}
\addtolength{\frameright}{8mm}
\setlength{\frametop}{2cm}
\setlength{\framebottom}{-\textheight}
\addtolength{\framebottom}{-\frametop}
%\addtolength{\framebottom}{6cm} % extra lift at bottom
\addtolength{\framebottom}{4cm} % extra lift at bottom
\psframe[linecolor=myblue,cornersize=absolute,linearc=2]
(\frameleft,\framebottom)(\frameright,\frametop)
% Border finished 

% I prefer a large sans serif font throughout
%\large\sffamily
\sffamily

\centerline{\myblue\huge\bfseries\url{http://www.calibayes.ncl.ac.uk/}}

\vspace*{1cm}
\hrule
\vspace*{1cm}


% title stuff
\begin{minipage}{0.840\textwidth}
\Huge \bfseries
{\myred {\veryHuge Bayesian calibration of biological simulation models}}

\vspace*{1cm}
\large \bfseries
{\myblue
{\LARGE D.A.Henderson, W.E.Wolski, R.J.Boys,
  T.B.L.Kirkwood and 
  D.J.Wilkinson}  \\ \\
Schools of Mathematics \& Statistics and Clinical Medical Sciences \\
University of Newcastle upon Tyne, UK
}
\vspace{1cm}

%{\Large\mygreen $^{\dagger}$\texttt{d.j.wilkinson@ncl.ac.uk}} \qquad 
%{\Large\mygreen $^{\ddagger}$\texttt{w.e.wolski@ncl.ac.uk}}
\end{minipage}
\begin{minipage}{0.15\textwidth}
\centerline{\includegraphics[width=8cm]{crest_colour}}
\end{minipage}


\vspace*{0.6cm}
\hrule
\vspace*{0.6cm}
% end of title stuff

% now go into multi-column mode for the poster content
\columnsep=100pt
\columnseprule=3pt
\begin{multicols}{\numcols}

%\begin{abstract}
\centerline{\mygreen Summary}
\vspace{1cm}
\hfill
\begin{minipage}{0.8\colwidth}
\raggedright\normalsize

Modern Systems Biology is greatly concerned with using post-genomic
experimental data sources in order to ``fit'' dynamic and predictive
computer models of complex non-linear biological processes. If this is
successful, scientists will be able to conduct ``in silico''
biological experiments that would be impossible to carry out ``in
vivo'' and these could, in turn, transform our understanding of
molecular biology. Unfortunately there are a great many obstacles in
the way of this grand vision, and many are fundamentally statistical
in nature.  Here,
an outline will be presented of the key issues involved in model
calibration, and the advantages that Bayesian calibration offers over
other commonly used approaches. The methodology is applied to a
case study concerning the cell-cycle in frog eggs.





%Therefore current research focuses on
%computational~\shortcite{Moles2003,Klipp2004} and
%statistical~\shortcite{BoysWK-u,GolightlyW-u} methods to 
%infer parameter values from experimental high throughput data~\cite{Weckwerth2002}.

%\vspace{1cm}


%\end{abstract}
\end{minipage}
\hfill\mbox{}
% start of poster proper

\sectionpost{Introduction}

\begin{itemize}
\item Large, complex, quantitative biological models are gaining
increasing significance for research hypothesis formulation in
biological and medical science \shortcite{Kirkwood2003,Proctor2005}. 

\item These
models are typically quantified using a set of parameters, such as
kinetic rate constants. 

\item Determining appropriate values for these
parameters directly from experiments can be a challenging task. 

\item The
process of determining parameter values based on a comparison of
experimental data with output from a computer model is known as
calibration \cite{KennedyO01}. 

\item Calibration of complex computer codes
is a well-studied problem in a range of different scientific
disciplines. It turns out that the calibration problem is particularly
well-suited to a Bayesian non-parametric statistical analysis.  

\end{itemize}

\begin{comment}
The primary aim of the \cb\ project is to capitalise on the
development of GRID-based modelling and simulation resources such as
BASIS ({\mygreen\texttt{www.basis.ncl.ac.uk}}) by building a higher level
computational GRID facility designed for inferring parameters of
simulation models.  The goal is to find a combination of
parameters that give the best match to available experimental data,
together with some idea of the uncertainty in these parameter values.
%The main \cb\ service carefully selects parameter
%combinations at which the simulator will be run. 
%These choices are
%made sequentially in order to make efficient use of information from
%previous runs. 
The inference
engine is implemented in the statistical
computing language \textbf{R}~\cite{Rlang} with optimised routines written in C and
C++, and communicates with the BASIS simulator using
\texttt{SSOAP}~\cite{SSOAP}.  

%The inference engine
%manipulates models specified in the Systems Biology
%Markup Language (SBML)~\shortcite{Hucka2003}, submits the SBML
%model description to BASIS and retrieves the simulation results using
%standard web-service (WS) protocols.


\end{comment}

\begin{comment}
\sectionpost{Computing architecture}


The \cb \ project will adopt a service oriented architecture (SOA)
using  \texttt{RSOAP}~\cite{Warnes04} as it provides a simple API for managing and communicating with multiple \textbf{R} processes using SOAP.
We are working to develop a standard WS~API for biological simulators as well as one to interface with online databases of experimental post-genomic data. Further important research tasks of the \cb\ project are to provide protocols and methods to define inference tasks.

\end{comment}

\sectionpost{Example: cell cycle in  frog eggs}

We demonstrate a Bayesian calibration approach by analysing a model
of the cell cycle in frog eggs \cite{ZwolakTW05}. 
\begin{itemize}
\item The model consists
of eight Michaelis-Menten type kinetic equations describing four
reversible reactions. 

\item The total number of parameters in the model is
$12$, out of which four parameters are Michaelis $K_m$-constants which
were preset. The collection of the eight parameters is denoted $\vec{\theta}=(\theta_1,\ldots,\theta_8)$.

\item Using the experimental data, made up of eight 2-4
timepoint short time series ($n=25$ datapoints in all), we were interested
in determining values for the eight parameters.

\item The computer model (simulator) for this system is a set of
  ordinary differential equations (ODEs), for which a fast numerical
  solver was available.

\end{itemize}

\subsection{Bayesian calibration}
\label{sec:bayescal}

\begin{itemize}
\item The computer model  is fast to run and so can be  evaluated at many sets of parameter values.

\item It can therefore be treated as a deterministic (non-linear)
  function, $f(\vec{x};\vec{\theta})$, of the parameter values $\vec{\theta}$ and any covariate
  information $\vec{x}$, such as initial quantities, experimental conditions, etc.

\item This deterministic function is embedded into a stochastic model,
  which relates the experimental data $\vec{y}=(y_1,y_2,\ldots,y_n)$ to the computer model via the
  regression relationship:
\[
y_i = f(x_i;\vec{\theta})+\sigma\epsilon_i.
\]
We take $\epsilon_i\sim t_\nu$ (independently)
for $i=1,2,\ldots,n$, that is, the stochastic component of the model
is independent, zero mean Student $t$--distributed noise. We present
results for $\nu=3$ degrees of freedom -- a distribution with much heavier tails than a Gaussian distribution.

\item With the above stochastic model, the likelihood
  $L(\vec{\theta},\sigma;\vec{y},\vec{x})$ is easy to write down.

\item Prior uncertainty about parameter values is expressed through uniform $U(-10,10)$ distributions for
  $\log(\theta_i)$, $i=1,2,\ldots,8$, and for $\log(\sigma)$.

\item Standard Markov chain Monte Carlo (MCMC) methods can be used to sample from the
  (analytically intractable) 
  posterior distribution of the unknown parameters $\pi(\vec{\theta},\sigma | \vec{y},\vec{x})$.

\end{itemize}

\subsection{Results}

Results for two of the eight model parameters are shown in
Figure~1.

\begin{myfig}{
Scatterplot of sampled values from the joint posterior distribution of
parameters $\theta_1\equiv v'_d$ and $\theta_2\equiv v''_d$, together with histograms
representing their marginal
posterior densities. The blue dots represent the values estimated
in \shortciteN{ZwolakTW05}.}
\includegraphics[width=0.75\figwidth,angle=-90]{parameter-post-nlr-t3}
\end{myfig}

\bigskip 

The marginal posterior distributions highlight the uncertainty about
the values of the parameters after observing the data in a way that is
not easily obtained from a non-Bayesian analysis.

\medskip

Our inference for the parameters $v'_d$ and $v''_d$ is largely
consistent with that of \shortciteN{ZwolakTW05}. This is true for
most, but not all of the eight parameters.  

\medskip

However, our main focus is on predicting the output of the model, whilst
accounting for the posterior uncertainty in the parameters.  This
allows us to perform \emph{in silico} experiments; investigating
the range of plausible outputs for a given set of experimental conditions.

\bigskip

Figure~2 illustrates how posterior uncertainty in the $\vec{\theta}$
values impacts on the output of the deterministic model
$f(\vec{x};\vec{\theta})$, for experimental conditions C and D from
\shortciteN{ZwolakTW05}.

\begin{myfig}{ The red dots are the experimental data. The blue line
    is the simulator output corresponding to the set of parameters
    which gave the best match to the data (in terms of likelihood).
    The red dashed line is the output corresponding to the estimated
    parameters from \shortciteN{ZwolakTW05}. The green line is the
    point-wise posterior mean of $f(\vec{\tilde{x}};\vec{\theta})$, and the blue
    shading gives point-wise 95\% posterior probability intervals for
    $f(\vec{\tilde{x}};\vec{\theta})$.}
\includegraphics[width=0.75\figwidth,angle=-90]{predictive-nlr-t3}
\end{myfig}

\bigskip 

As well as uncertainty, a point estimate of the parameter values can
be obtained from the parameter values giving the best fit to the
experimental data over the course of the MCMC simulation.  The values
we obtain give a better fit to the data (in terms of loglikelihood) than the estimated parameters from \shortciteN{ZwolakTW05}.

\medskip

Overall, this simple Bayesian approach seems to work well, even in this
case-study involving sparse data, from different sources. The Bayesian
approach allows the uncertainty in the parameter values and model
output to be quantified in an appropriate way.


\sectionpost{Inferential challenges}

\begin{itemize}
\item The Bayesian calibration approach described in Section~2 relied
  on the fact that the computer model was fast to run, and so could be
  evaluated many times. For large and/or complex biological models,
  running the simulator at many input points can be prohibitive.

\item Biological processes are intrinsically \emph{stochastic}, and as
  such deterministic models can not always be used for making
  inferences.  Perhaps the main obstacle in the way of calibration of
  stochastic simulators is the fact they are likely to be too slow to
  be used in a calibration process where many stochastic simulations
  are required.
\end{itemize}

When a (deterministic) simulator is slow to evaluate, one possible
solution is to emulate it using a  Gaussian process (GP) model.


\subsection{Emulation via Gaussian processes}

\begin{itemize}
\item Gaussian processes provide a flexible model for the output from
  a deterministic computer simulation model that can be fitted on relatively
  few data points.

\item The GP can then be used as a (fast) emulator of the computer
  model.
 
\item To illustrate this, we fit a Gaussian process to simulated data
from the frog model using the software
\texttt{PErK}~\shortcite{SantnerWN03}.  
%The $650$ input points were
%generated by a Latin hypercube design and cover a space in the range
%from $e^{-10}$ to $e^{4}$ for each parameter . 

\item The scatterplot in Figure~3 indicates that based on only 650
evaluations of the computer model, the GP can be used to emulate the
computer model reasonably well.  

\item In the case of the frog model,
evaluating the GP model at a set of parameter values $\vec{\theta}$ is approximately 4 times faster than
running the fast ODE solver. The speed improvement will become larger
for more complex models.

\item Once estimated, the GP can be used as a surrogate for the computer model in
  methods to estimate model parameters. 
\end{itemize}

\begin{myfig}{
Scatterplot of predicted against true agreement with experimental data for a set of 1500 input parameters sampled uniformly from the parameter space.}
\mbox{}\centerline{\includegraphics[width=0.7\figwidth,angle=-90]{scatter}}\mbox{}
\end{myfig}

\begin{comment}
\sectionpost{Search of the best parameters.}


We predicted the output of $500000$ input parameter tuples. The prediction using the Gaussian process is approximately 4 times faster than running the simulator even if in case of the Frog model it is a fast ODE solver. The speed improvement will become larger for more complex models or if using stochastic simulators. 

The histogram is showing the distribution of the predicted distance to the experimental data. The red vertical dashed line visualizes the threshold used to select regions of plausible parameters.

\begin{myfig}{Histogram of the Manhattan distances, to the experimental data, as predicted by the GP for $200000$ design points. The red vertical line indicates the threshold used afterwards to select candidate parameter region for further parameter search.}
\mbox{}\centerline{\includegraphics[width=\figwidth,angle=-90]{hist}}\mbox{}
\end{myfig}

The pair plot below  parameter values which have generated small values of the objective function. 

\end{comment}
 
\sectionpost{Future work}

\begin{itemize}
\item The ultimate goal of the CaliBayes project is to develop
  \emph{fully} Bayesian
methods which allow us to quantify all sources of uncertainty in the
calibration process. 

%\item For larger, more complex models, the computational advantages of using
%Gaussian processes to emulate the computer model will become more important.  

%\item We aim to employ more sophisticated Gaussian process modelling
%  techniques \cite{KennedyO01} in order to quantify
%  uncertainty in the values of the model parameters.

%\item Such methods allow predictions to take into account this
%  uncertainty in the parameter values.

%\item However, the problem of adapting these methods to deal with multivariate  outputs (as in the frog cell-cycle model) needs to be addressed.

\item The Bayesian calibration of \emph{stochastic} models is of
  particular interest due to the inherent variability in biological
  processes which cannot be represented by a deterministic model.

\item Progress has been made using exact and approximate Bayesian
  models (e.g.\ \citeN{BoysWK-u,GolightlyW05}) but these methods
  are not likely to be applicable when data are sparse, from different
  sources and of different types, and when the models are large and
  contain a mixture of low and high copy-number species.

\item When a (black box) stochastic simulator is too slow to make many
  thousands of simulations, some form of emulation using a fast
  surrogate will almost certainly be required.

%\item Once the statistical methods are in place our aim is to develop
%  a GRID-based inference engine that will utilise both GRID-based
%  simulators and GRID-based data repositories in order to accurately
%  estimate model parameters.  Intermediate web-services will be
%  required to \emph{transform} the simulator output into the same
%  format as the experimental data, and to \emph{compare} the
%  experimental and simulated data.

\end{itemize}

%\scriptsize
\footnotesize

\bibliographystyle{prepbib}  % Style BST file
\bibliography{abbrev,CaliBayes}  % Bibliography file (usually '*.bib' ) 


% end of poster proper
\end{multicols}

\vspace*{1cm}
\hrule

\begin{center}
\begin{minipage}{0.2\textwidth}
\rightline{\includegraphics[width=8cm]{bbsrc}}
\end{minipage}
\begin{minipage}{0.75\textwidth}
\leftline{\large\mygreen BBSRC Bioinformatics and e-Science
    Programme II, BBS/B/16550, January 2005 -- December 2007}
\end{minipage}
\end{center}

\end{document}

