\chapter{Multiple realizations}
\label{chap:random}


\mhead{Introduction}
In Chapter~\ref{chap:reg} and Appendix~\ref{chap:nonstat} I introduced
several different methods to characterize the covariance. 
In each case I presented  a single answer to the interpolation or
tomography problem.  In some ways this is exactly what we are looking for:
given a problem, give me `the answer'.  
The single solution approach
has a couple significant drawbacks.  First, the solution tends to 
have a low spatial frequency. Second, it does not provide information
on model variability or provide error bars on our model estimate.
Geostatisticians  have these
abilities  in their repertoire through
what they refer to as `multiple realizations'.
They introduce a random component, based on
properties of the model,   to their estimation procedure.  
Each realization's frequency content is more accurate and by
comparing and contrasting  the equi-probable realizations models,
variability can be assessed.
\par
In this appendix I present a method to achieve the same goal
using geophysical techniques.  I modify the model styling goal,
replacing the zero mean vector with a random vector.  I show  how
the resulting models have a more pleasing texture and can
provide information on variability.

\mhead{Multiple realizations}


The missing data problem is probably the simplest to understand and
interpret results.  As explained earlier, we solve the
problem in its preconditioned form using
\beqa
\data &\approx& \bf J \prec \pvar \\ \nonumber
\zero &\approx&  \iop \pvar ,
\eeqa
where:
\begin{description}
\item [$\data$] is a binned version of our known points
\item [$\bf J$] is our known data selector
\item [$\prec$] is our preconditioning operator that characterizes the
model's covariance (anything from a PEF to a steering filter to 
a non-stationary PEF)
\item [\pvar] is our preconditioned variable, which is equal to
$\reg \model$.
\end{description}
We assume
that  $\reg$ accurately and fully describes the model spectrum. 
This is rarely the case.  A steering filter is only able to handle
a single, smoothly varying dip.  A large enough PEF could, in theory,
perfectly describe the model's spectrum. Using a large PEF
is infeasible because our model is generally not stationary 
      and, in a missing  data problem, we
usually do not have enough valid equations to estimate it.
Estimating a  non-stationary PEF is also often impossible because
of shortage of equations.
\par
So what happens to  the model spectrum that our PEF fails
to describe?  I hypothesize that the $\zero$ in the model styling goal
results in creating a model which has zeros in the spectrum corresponding
to the zeros in $\reg$'s spectrum. 
On one hand, this seems like the most intelligent choice.  When we don't
have any information, don't guess.  An alternative  is to 
to try random values where don't know the spectrum.
We could think of applying a trick similar to the way I handled
smoothing slowness rather than the change in slowness in Chapter~\ref{chap:reg}.
We begin by   writing our model styling goal for the regularized problem as,
\beqa
\zero \pox \left( \bf B +  \bf A \right) \model \\ \nonumber
- \bf B \model \pox  \bf A \model ,
\eeqa
where $\bf B$ is a filter describing the spectrum not described by $\bf A$.
Of course, we don't know $\bf B$  or $\model$ {\it a priori}. We can simulate
different values of $\bf B \model$ by replacing our $\zero$ with the random
vector $\bf v$.
If we replace $\zero$ with  $\bf v$ we get,
\beqa
\bf d \pox  \bf J \prec \pvar  \nonumber \\ 
\bf v \pox \epsilon \bf I \pvar\label{eq:rand} .
\eeqa

The magnitude of $\bf v$ should be approximately the average of
the absolute values  of the  residual  using $\zero$ for $\bf v$.
The  larger the
value of $\epsilon$, the  smoother our model estimate.



If  we use a steering
filter as $\reg$ with the  fitting goals 
(\ref{eq:rand}) on the well log interpolation problem,
Figure~\ref{fig:well-logs}, we obtain Figure~\ref{fig:well.movie}.  Note
how each images frequency content is much closer to the original models
content. 

\plot{well.movie}{width=6.0in,height=3.0in}{
Three different  random realization results of interpolating Figure~\ref{fig:well-logs}.}



\mhead{Seabeam}

A more interesting interpolation  problem is the SeaBeam dataset\cite{gee}.
Figure~\ref{fig:sea.init} shows a day's worth of data
collected by SeaBeam, which measures
water depth under, and to the side, of a  ship. 

\sideplot{sea.init}{width=3.0in,height=3.0in}{Depth of the ocean under ship tracks.}


Figure~\ref{fig:sea.pef} shows the result of estimating a PEF from the
known data locations and then using it to interpolate the entire mesh.
Note how the solution has a lower spatial frequency as we move away from
the recorded data and the original tracks of the ship are still
clearly visible.  

\sideplot{sea.pef}{width=3.0in,height=3.0in}{Result of using a PEF to
interpolate Figure~\ref{fig:sea.init}, taken from GEE \cite{gee}.}

Figure~\ref{fig:sea.movie} shows nine different realizations
by applying (\ref{eq:rand}).  The texture of each realization 
matches better the texture of the known data (Figure~\ref{fig:sea.pef})
and we get a range of possible interpretations.  For example,
the top and  bottom right realizations vary significantly in
the way the handle the ridge going down the center of the image.

\plot{sea.movie}{width=6.0in,height=8.0in}{
Nine different realizations of the Seabeam interpolation problem. 
Note how the realizations  vary away from the known data points.  }

The left panel of Figure~\ref{fig:sea.var} provides a measure of
the variance of the various realizations.  Note, how as expected,
the larger the gap, the more variance in the estimations. The right
panel of Figure~\ref{fig:sea.var} shows the average of all of the realizations.
As expected the average image is comparable to the PEF result.


\plot{sea.var}{width=6.0in,height=3.0in}{
The left panel is the variance in the prediction for the different realizations
starting from Figure\ref{fig:sea.init}. The right panel is the average of
all of the realizations.  Note the similarity with Figure~\ref{fig:sea.pef}.}