Biomedical Engineering Reference
In-Depth Information
expectation of some function
f
(
x
) with a nonuniform probability
p
(
x
). This is just
an integral, ,
f
(
x
)
p
(
x
)
dx
, so we could sample points uniformly and compute
f
(
x
)
p
(
x
) for each one. But if some points have very low probability, so they only
make a small contribution to the integral, spending time evaluating the function
there is a bit of a waste. A better strategy would be to pick points according to
the actual probability distribution. This can sometimes be done directly, espe-
cially if
p
(
x
) is of a particularly nice form. A very general and clever indirect
scheme is as follows (14). We want a whole sequence of points,
x
1
,
x
2
, ...
x
n
. We
pick the first one however we like, and after that we pick successive points ac-
cording to some Markov chain: that is, the distribution of
x
i
+1
depends only on
x
i
,
according to some fixed function
q
(
x
i
,
x
i
+1
). Under some mild conditions,
17
the
distribution of
x
t
approaches a stationary distribution
q
*(
x
) at large times
t
. If we
could ensure that
q
*(
x
) =
p
(
x
), we would know that the Markov chain was con-
verging to our distribution, and then, by the ergodic theorem, averaging
f
(
x
)
along a trajectory would give the expected value of
f
(
x
). One way to ensure this
is to use the "detailed balance" condition of the invariant distribution, that the
total probability of going from
x
to
y
must equal the total probability of going
the other way:
p
(
x
)
q
(
x
,
y
) =
p
(
y
),
[39]
qxy
(, )
py
()
=
w
hxy
(, )
.
[40]
qyx
(, )
px
()
So now we just need to make sure that [40] is satisfied. One way to do this is to
set
q
(
x
,
y
) = min(1,
h
(
x
,
y
)); this was the original proposal of Metropolis et al.
(141). Another is
q
(
x
,
y
) = (
h
(
x
,
y
))/(1 +
h
(
x
,
y
)). This method is what physicists
usually mean by "Monte Carlo," but statisticians call it
Markov chain Monte
Carlo
, or "MCMC." While we can now estimate the properties of basically arbi-
trary distributions, we no longer have independent samples, so evaluating the
accuracy of our estimates is no longer a matter of
trivial
probability.
18
An im-
mense range of refinements have been developed over the last fifty years, ad-
dressing these and other points; see the further reading section for details.
Keep in mind that Monte Carlo is a stochastic simulation method only in a
special sense—it simulates the probability distribution
p
(
x
),
not
the mecha-
nism that generated that distribution. The dynamics of Markov chain Monte
Carlo, in particular, often bear no resemblance whatsoever to those of the
real system.
19
Since the point of Monte Carlo is to tell us about the properties of
p
(
x
) (what is the expectation value of this function? what is the probability of
configurations with this property? etc.), the actual trajectory of the Markov
chain is of no interest. This point sometimes confuses those more used to direct
simulation methods.
