#### Howdy, Stranger!

It looks like you're new here. If you want to get involved, click one of these buttons!

Options

# Monte Carlo Methods in Climate Science

One way the Azimuth Project can help save the planet is to get bright young students interested in ecology, climate science, green technology, and stuff like that. So, as planned out a while back, David Tweed and I have been writing an article for Math Horizons, an American magazine for undergraduate math majors.

Here's a draft:

It's been really a delight for me to write this, since David Tweed has put a lot of work into things I can't do very well - writing software illustrating the concepts - and this has let me focus on what I'm good at, namely writing. But so far he's done most of the work, and now it's time for me to polish the writing a few more times.

We'd really like to hear your comments! Personally, the kind of comments I don't want are ones asking for more technical details. The kinds I do want are ones where you sincerely don't understand something, or feel that students would have trouble understanding something. We're trying to cover a few big ideas rather rapidly, and that's hard to do well.

• Options
1.

Note that some freely downloadable examples can be found here if anyone wants to see the kind of level these things are pitched at. The latest article, the one about the mathematics of choice, inspired the section on "Exploring the topic" with its section of questions about the topic of the article, but we still have to come up with a good but simple environmental question to supplement the two fully mathematical ones! Any ideas gratefully accepted.

Comment Source:Note that some freely downloadable examples can be found [here](http://www.maa.org/publications/periodicals/math-horizons/math-horizons-sample-articles) if anyone wants to see the kind of level these things are pitched at. The latest article, the one about the mathematics of choice, inspired the section on "Exploring the topic" with its section of questions about the topic of the article, but we still have to come up with a good but simple environmental question to supplement the two fully mathematical ones! Any ideas gratefully accepted.
• Options
2.

David: to speed up the process of getting suggestions, I'm going to throw our draft onto the blog in the form of a blog article, and solicit suggestions there.

Comment Source:David: to speed up the process of getting suggestions, I'm going to throw our draft onto the blog in the form of a blog article, and solicit suggestions there.
• Options
3.

Hi! I'm in the process of reading it! It's well written...

Comment Source:Hi! I'm in the process of reading it! It's well written...
• Options
4.

Finally, P(M = m) is independent of our choice of settings. So, we can use Bayes' rule to compute P(S = s|M = m) up to a constant factor. And since probabilities must sum to 1, we can figure out this constant.

The important thing about P(M = m) is that it is constant (at least for the duration of the analysis) - because m is the "hard disc full of observed data" which does not change. And while the constant factor can be figured out in principle it is not often done in practice (it would take too long). A nice thing about MCMC is that it can explore an unnormalized density and you can extract lots of useful information from the sampled values without needing to know the constant.

The key to making this work is that at each step on the walk a proposed modification s' to the current settings s is generated randomly - but it may be rejected if it does not seem to improve the estimates. The essence of the rule is: The modification s to s' is randomly accepted with a probability equal to the ratio P(M = m|S = s')/P(M = m|S = s) Otherwise the walk stays at the current position. If the modification is better, so that the ratio is greater than 1, the new state is always accepted.

I suggest replacing this with:

At each step on the walk a proposed modification s' to the current settings s is generated randomly which is then either accepted or rejected. If P(M = m|S = s') is bigger than P(M = m|S = s), then the modification is accepted and the walk steps to s'. Otherwise, the modification s to s' is randomly accepted with a probability equal to the ratio P(M = m|S = s')/P(M = m|S = s). If the proposed modification is rejected the walk stays at the current position s. It is not obvious that this will make a walk with the right statistical properties, but this can be proved for lots of kinds of proposed modifications [Metropolis et al, 1953].

Figure 4 is very nice!

Comment Source:> Finally, P(M = m) is independent of our choice of settings. So, we can use Bayes' rule to compute P(S = s|M = m) up to a constant factor. And since probabilities must sum to 1, we can figure out this constant. The important thing about P(M = m) is that it is constant (at least for the duration of the analysis) - because m is the "hard disc full of observed data" which does not change. And while the constant factor can be figured out in principle it is not often done in practice (it would take too long). A nice thing about MCMC is that it can explore an unnormalized density and you can extract lots of useful information from the sampled values without needing to know the constant. > The key to making this work is that at each step on the walk a proposed modification s' to the current settings s is generated randomly - but it may be rejected if it does not seem to improve the estimates. The essence of the rule is: The modification s to s' is randomly accepted with a probability equal to the ratio P(M = m|S = s')/P(M = m|S = s) Otherwise the walk stays at the current position. If the modification is better, so that the ratio is greater than 1, the new state is always accepted. I suggest replacing this with: At each step on the walk a proposed modification s' to the current settings s is generated randomly which is then either accepted or rejected. If P(M = m|S = s') is bigger than P(M = m|S = s), then the modification is accepted and the walk steps to s'. Otherwise, the modification s to s' is randomly accepted with a probability equal to the ratio P(M = m|S = s')/P(M = m|S = s). If the proposed modification is rejected the walk stays at the current position s. It is not obvious that this will make a walk with the right statistical properties, but this can be proved for lots of kinds of proposed modifications [Metropolis et al, 1953]. Figure 4 is very nice!
• Options
5.

This a good book for further reading for the seriously interested.

Jun S Liu, Monte Carlo Strategies in Scientific Computing, Springer.

Comment Source:This a good book for further reading for the seriously interested. Jun S Liu, [Monte Carlo Strategies in Scientific Computing](http://www.springer.com/statistics/statistical+theory+and+methods/book/978-0-387-76369-9), Springer.
• Options
6.

Nice paper. Here are some trivial remarks you probably already found:

Repeated 'that' which I'm not sure you wanted (para 3 pg. 1)

The most popular theory is that that a ...

Last sentence para 3, extra 'i'

so when iit ...

When you open the Box models section, I think there should be a comma after climate.

Para 2 Box models section

containing related strongly related things

Space needed, top of page 2

ofCO$_2$ leaving

I think there's an extra - before but that should not be there.

Para 2, page 2,

Second, it turn out

Not necessarily, but after you say “future hight”, I think it might be nice to add “of the person in the swing” or something like that.

That's up to the section, Probability modeling.

Comment Source:Nice paper. Here are some trivial remarks you probably already found: Repeated 'that' which I'm not sure you wanted (para 3 pg. 1) > The most popular theory is that that a ... Last sentence para 3, extra 'i' > so when iit ... When you open the *Box models* section, I think there should be a comma after climate. Para 2 Box models section > containing related strongly related things Space needed, top of page 2 > ofCO$_2$ leaving I think there's an extra - before but that should not be there. Para 2, page 2, > Second, it turn out Not necessarily, but after you say “future hight”, I think it might be nice to add “of the person in the swing” or something like that. That's up to the section, *Probability modeling*.
• Options
7.

In para 1,

are element of some finite set

Shouldn't it be a multi-set and element(s)?

Sorry, I had to throw the multiset in there since this was brought up recently!

Missing bracket

$P(M = m$

Repeated words, before the section “An example”

to know to know that with 95%

Comment Source:In para 1, > are element of some finite set Shouldn't it be a multi-set and element(s)? Sorry, I had to throw the multiset in there since this was brought up recently! Missing bracket > $P(M = m$ Repeated words, before the section “An example” > to know to know that with 95%
• Options
8.
edited July 2013

Graham: firstly, you're right that the most important point is that $P(M=m)$ isn't needed for the MCMC acceptance ratio. I'm interested in why you think the definition of the rule with two cases (ratio>1, otherwise) is better? Ah, perhaps the second sentence ought to be run into the definition of the rule, the point just to make explicit the implicit poitn that the case where the ratio is greater than one then every random probability will be less than the ratio, and hence be accepted. (We've phrased it as is partly because we're admitting that we don't have the space to actually justify some of the stuff, in which case that distinction used in the proof doesn't get used. But if people think the other presentation is better we perhaps ought to change it.)

Jacob: I don't think we mean multiset; what we mean is that an individual measurement comes from a finite set of possibilities, with the "time series" for measurements forming a sequence or possibly, making an outrageous conflation of mathematical structures, a vector where the index corresponds to a timestep.

Everyone: thanks for pointing out the other typos and cases of phrasing that could be improved.

Comment Source:Hi, Thanks for all the helpful comments. Some quick replies: Graham: firstly, you're right that the most important point is that $P(M=m)$ isn't needed for the MCMC acceptance ratio. I'm interested in why you think the definition of the rule with two cases (ratio>1, otherwise) is better? Ah, perhaps the second sentence ought to be run into the definition of the rule, the point just to make explicit the implicit poitn that the case where the ratio is greater than one then every random probability will be less than the ratio, and hence be accepted. (We've phrased it as is partly because we're admitting that we don't have the space to actually justify some of the stuff, in which case that distinction used in the proof doesn't get used. But if people think the other presentation is better we perhaps ought to change it.) Jacob: I don't think we mean multiset; what we mean is that an individual measurement comes from a finite set of possibilities, with the "time series" for measurements forming a sequence or possibly, making an outrageous conflation of mathematical structures, a vector where the index corresponds to a timestep. Everyone: thanks for pointing out the other typos and cases of phrasing that could be improved.
• Options
9.

a probability equal to the ratio P(M = m|S = s’)/P(M = m|S = s)

because probabilities can't be greater than one (despite what the latest blog says ;-)).

Comment Source:David, it was this bit I liked least > a probability equal to the ratio P(M = m|S = s’)/P(M = m|S = s) because probabilities can't be greater than one (despite what the latest blog says ;-)).
• Options
10.

Hi David! I'll look closer at this, but think collection would be a better term than set. It's an interesting article and should get some readers interested!

Comment Source:Hi David! I'll look closer at this, but think collection would be a better term than set. It's an interesting article and should get some readers interested!
• Options
11.

I think we're talking at cross purposes: the "finite set" bit just refers to a single scalar measurement, saying that it's not the infinite set of values that you'd get if it was the real numbers. To the extent I understand these things, the differences between a collection and a set don't seem relevant. (There's also the question of whether we need a term to describe a time series of values, but that should be a separate sentence.)

Comment Source:I think we're talking at cross purposes: the "finite set" bit just refers to a single scalar measurement, saying that it's not the infinite set of values that you'd get if it was the real numbers. To the extent I understand these things, the differences between a collection and a set don't seem relevant. (There's also the question of whether we need a term to describe a time series of values, but that should be a separate sentence.)
• Options
12.

Hi David, I will read the paper again today and see if I have any other ideas/constructive feedbacks. I will look closer at the set business, but probably it's all OK and just a misunderstanding on my part. Cheers!

Comment Source:Hi David, I will read the paper again today and see if I have any other ideas/constructive feedbacks. I will look closer at the set business, but probably it's all OK and just a misunderstanding on my part. Cheers!