intro.tex

\section{Introduction}

Since the beginning of recorded time humans wondered whether we are
alone in the universe or whether there are other planets that might
support life. The scientific quest for extrasolar or exoplanets
(planets beyond our own solar system) began in the mid 19th century,
but it was only as recent as 1992 that astronomers were  able to
confirm the existence of other planetary systems.  Since then
\uhoh{describe resources, telescopes, etc} have been devoted to
identifying stars with orbiting planets, and in particular, those that
might host earth-like planets. As of June 2010 over 450 exoplanets
have been detected, the majority using techniques that rely on
measurements of a star's radial velocity.  The presence of a planet
orbiting a star will lead to a periodic ``wobble'' in the the star's
radial velocity, while if there are no planets present the expected
velocity is constant although measurements fluctuate because of
stellar jitter and other sources.  Based on a short sequence
of noisy measurements, the challenge is to determine whether there are
no planets or one, two or more planets, and if there are planets
present to characterize their orbital parameters.  If the data do not
provide strong evidence either in favor or against the presence of
planets, additional observations may be needed to provide
confirmation.  Determining which stars to observe and when to schedule
scarce telescope time optimally is critical given limited resources;
Bayesian methods in particular are ideally suited for updating
information/evidence and decision making in such sequential sampling
problems.


The primary question of whether there are no planets or one or more
may be poised as a Bayesian model selection problem, where the
collection of models under consideration are the zero-planet model
$\M_0$, one-planet $\M_1$, two-planet $\M_2$, etc.  In Bayesian
statistics, the posterior probability of a model $\M_j$ requires
accurately calculating the integrated or marginal likelihood (see
\citet{clyd:geor:2004} for a general overview or
\citet{clyd:etal:2007,ford2006bms} for some of the issues that arise in
astronomy).  In the exoplanet application, each model $\M_p$ has
$2+5p$ parameters where $p$ is the number of planets, so even the
simplest single planet model requires a seven dimensional integral
over the parameter space to obtain the marginal likelihood.  Because
of the nonlinear nature of the models, analytic solutions do not
exist.  While low in dimension compared to many modern problems of
model selection, the models are highly nonlinear making Laplace
approximations untenable given the modest sample sizes and multi-modal
posterior distribution; constructing good proposal distributions for
Reversible Jump-Markov Chain Monte Carlo (RJ-MCMC) algorithms is near
impossible; just designing an efficient MCMC algorithm for sampling
from posterior distributions for a model with one planet is
challenging on its own due to the multi-modal nature of the posterior
distribution \citep{ford2006bms}.   Even if an efficient MCMC
algorithm may be designed, estimates of the  marginal likelihood from
the output are typically unstable.


Out of available methods for numerical integration, importance
sampling (IS) generally produces the best estimates of marginal
likelihoods in terms of bias and variance \uhoh{reference Christian's
  comparisons?}. The key to building an efficient IS algorithm is to
design an importance function that mimics the integrand.  Building
such a function can be quite challenging even in lower to modest
dimensional settings; with the multi-modal posterior distributions
that arise in the exoplanet setting the task is all the more
ambitious!  In this paper, we propose a novel adaptive method for
constructing an importance function using a mixture of Student-$t$
distributions that gradually evolves to the target posterior
distribution.  This algorithm falls in the class of adaptive  Monte
Carlo methods \citep[see][for an overview and papers appearing in the
special issue on sequential Monte Carlo]{fern:2008}.  This adaptive
importance sampler permits accurate determination of posterior model
probabilities,  distributions of parameters of interest, as well as
sequential updating as data acrue.


The paper is organized as follows.  We begin by presenting an overview
of the radial velocity models for exoplanet detection in Section
\ref{sec:exo} and discuss the scientific issues that inspired this
work. In Section \ref{sec:IS} we review importance sampling in the
context of calculating marginal likelihoods for Bayesian model
selection. In Section \ref{sec:AAIS} we present the Adaptive Annealed
Importance Sampling method developed for the exoplanet problem. This
method integrates several key features of other methods: mixture
models, annealing, and adaptation via a sequence of distributions to
build a flexible importance distribution that mimics the posterior
distribution, but incorporates ideas successfully used in
trans-dimensional methods to adapt the number of kernels in the
mixture by utilizing split and merge moves.  In Section
\ref{sec:simulation} we present two simulation studies that
demonstrate the efficiency of the method compared to other approaches.
The first involves a challenging flared helix in three dimensions,
while the second seven dimensional problem that was selected to reflex
the target function in the exoplanet problem.  In Section
\ref{sec:app} we illustrate the method in the analysis of several real
star systems with one and two planets. Finally, we conclude with
discussion of possible extensions of the method in Section
\ref{sec:disc}.