Lumpy vs Model E
Over at Climate Audit, Steve Mosher has been showing readers who Model E hindcasts HadCrut data. He suggested I show how “Lumpy” hindcasts. And, as we’ve been discussing trying to tease out time constants using assumptions made by Schwartz in his 2007, I thought I’d show you the method I’d really prefer to estimate the time constant and climate sensitivity: that’s Lumpy.
Here’s how she well she hindcasts compared to Model E:
If you need metrics that describe the fit, the only ones I have are these: Standard error for Lumpy, from 1880-1999: σlumpy =0.16C. The standard error for Model E during the same time period: σmodelE =0.17C.
Of course, this isn’t quite fair, as Lumpy is literally fit to GISS data which tracks the HadCrut data well. However, she only has two parameters: Time constant and climate sensitivity. I won’t reveal those until after we’ve discussed the “splendidness” or “crumminess” of Lumpy. That way, you can decide if you love her or hate her without knowing too much about what she says about the future climate.
What is Lumpy?
“Lumpy” is the simplest possible lumped parameter model to describe the earth’s climate. It assumes the earth’s climate is isothermal (i.e. the same temperature, T, exists everywhere) and applies a linearized form of conservation of energy. This results in the same ordinary differential equation suggested by Schwartz (2007). Here I’ve re-organized to express constant, τ and the climate sensitivity S:
To estimate the time constant τ, and sensitivity, S, I “drove” the model using forcings used in the GISS Model E and fit the monthly Land/Ocean temperature data from 1880-1999. (I stopped in 1999 because I have no monthly forcings for later times.)
Why are Lumpy’s temperature predictions so smooth?
Lumpy is driven by smoothed forcing data. Also, because she’s such a simple model, she has no internal dynamics. So, her projections are smoothed. To obtain “noisy” predictions, I would either need “noisy” forcing or a more complex model. (Ideally, perfection requires both.)
What sort of noisy forcings might I apply to Lumpy to make things “noisy”? Well, for example, the GISS forcings show the Total Solar Irradiance as an more or less periodic smoothly varying function with a period of about 11 years. Measured irradiances look like this:
Figure 2: Image from ‘pmod wrc: Solar Constant’. In this figure we see forcing from the sun that few would describe as smoothly varying over a period of 11 years. Certainly, averaged over periods like a month, one sees more or less smooth variation, but over shorter periods, the Total Solar Irradiance (TSI) is sizable.
Based on known theoretical behavior of equation (1), the noise in the TSI would be sufficient to create visible noise in Lumpy’s temperature predictions. (It might not make them any better, but it would make them noisy. 
)
It’s worth noting that model E, which does contain physics to create some weather, does shown some “noisy” variability. However, it is also smoother than Hadcrut data. This could be because of measurement noise, is also due to averaging over 5 model runs, but some excess smoothness might be due to the smoothed forcing functions. These could include the smoothed solar irradiance, smoothed behavior or aerosols after eruptions, smoothing of the effect of minor volcanos, and pretty much anything smoothed in the input files.
So, which is better?
Well….
Obviously, Lumpy is not a sophisticated model. She doesn’t explain much. But, she does show how we can, in principle, use physics when curve fitting. Also, there is some glimmer of hope we can get a more empirically based estimate of the time constant and climate sensitivity with a tool like this. (Or not. But if we can’t, then we can’t. )
Unfortunately, I also didn’t do the fitting in a sophisticated way, so I don’t have a clue what the uncertainties on my parameters are. But, learning issues I need to make Lumpy better (and/or figure out the error bars on the parameters) is one of the real reasons I’m fiddling with the even simpler (Schwartz + Meaurement noise) model and trying to learn how to best estimate the uncertainty bounds (as well as figure out if there are ways to look at that data to make the uncertainties smaller.)
While some may prefer other courses, I actually prefer to learn how to use the various statistical tools on the simplest possible method where I don’t expect to get a precise final answer. Then, later, I can apply the methods to the version of the model that I think more likely to give better answers. (Though, of course, Lumpy may be totally hopeless too.)
Update: I should never write these before coffee. I couldn’t sleep, wrote as I drank coffee and hit “publish” as my husband was hurrying me out the door to meet the inlaws for breakfast.
We went an hour earlier than usual because the inlaws golf outing was canceled due to rain. All through breakfast I thought “Dang! Comments will be full of people reminding me Model E is an ensemble average.” Indeed they were. I nserted about 10 words to reflect this.)
Comments
Basil (Comment#2708) May 11th, 2008 at 6:54 am
Both GISS-E and Lumpy don’t show a fast enough rate of increase from ~1910 to ~1940. And both miss the downturn at ~1940. Any idea why?
But Lumpy does a little better than GISS-E. I’m impressed.
anon and on (Comment#2709) May 11th, 2008 at 7:01 am
The modelE data is an ensemble mean - that’s why it smoother than the real world.
steven mosher (Comment#2710) May 11th, 2008 at 7:40 am
Basil, Lumpy, like modelE depends on data for forcing. So for example there is data file of forcing due
to volcanos. Also, GHG forcing, solar forcing, aerosol forcing, etc. As one goes back in time these forcings become
more uncertain, so the hindcast is more uncertain. Also, there is the nasty beast of internal variability.
The common beleif is that internal variability evens out over the long term. However, folks like
bender, UC, and me point to things like 1/f ,systems with long memories, etc and so one could attribute
the mismatches to a.) data uncertainity in forcings prior to 1960. b.) observtion uncertainity in early years. c)
longer term cycles in internal variability. d.) model inefficency. or some combination of all 4.
steven mosher (Comment#2711) May 11th, 2008 at 7:44 am
lucia, when you fit lumpy to giss did you rebaseline giss to 1961-90? my modele and hadcru are both
a 1961-90 anomaly. To rebaseline giss simply take the mean anomlay of giss from 1961-90 and subtract
it from the entire series ( should be about .092C)
lucia (Comment#2712) May 11th, 2008 at 7:46 am
Fred & Steve–
I clicked publish, and went out to breakfast with the inlaws. While in line I was thinking “Oh shoot! I forgot to mention the averaging.” Ahh… for the days when no one read my blog and I could fix boo-boos hours later!
Yes. The Model E is smoother also because of the averaging– that should knock variability down sqrt(5)!
Tom Gray– Lumpy can do forecasts on the same basis as Model E. I need forcing data. Of course, Lumpy may be crummier….
lucia (Comment#2713) May 11th, 2008 at 8:08 am
steve– I did rebaseline to set the mean temperatures the same. The reason for this is that, based on the definitions of the anomolies and the forcings, I actually also find the best baseline that fits one to the other. So, I just found the average difference and slid.
This also gives Lumpy an unfair edge here of the St.devations of the uncertainties.
The other thing that gives it an unfair edge on st. deviations is Model E, even though averaged, still contains some “weather noise”, just less than otherwise. At least some of that noise is uncorrelated with the ensemble that is the earth’s trajectory of temperatures and so adds to the standard deviation.
At one point, I unintentionally computed the standard deviation between Lumpy and ModelE. It’s about 0.08C. So, they track each other much better than they track the weather– which is as it should be since both use the same forcing file and both should conserve energy overall.
basil-
No. Neither Lumpy nor model E do well during 1910-1940.
lucia (Comment#2714) May 11th, 2008 at 9:28 am
Basil–
I forgot the possible reasons “why” for differences:
1) The forcing functions are incorrect, even when seen as smoothed.
2) The deviation is due sort of “weather noise” we expect to cause deviations between mean behavior the behavior of an single realization of whether that can happen when the system is driven by perfectly smooth forcing. (This would be the “It’s Navier-Stokes like non-linearity” type explanations.)
3) The effect of random variations around the “slow varying” forcing functions could, hypothetically, result in excursions that are not due to the non-linear nature of the equations for climate but rather due to the noisy behavior of the driver. (That is: it’s not a Navier-Stokes like turbulence thing, but something we could see when we drive a simple “lump” with white noise.)
Reasons (1) and (2) are familiar to everyone. If the forcings are incorrect, your guess is as good as mine about possible reasons. I need to explore (3) using some order of magnitude arguments.
But to give you an idea of how much “white noise” can push a system think about this:
1) Suppose we look at the solar cycle, and subtract the 11 year periodic cycle. We’ll be left with noise. Eyeballing the graph, suppose that if we examine the “noise” each month we find the standard deviation of that TSI noise as experienced per square meter of earth is σTSI= 0.1 W/m^2 from month to month Now, suppose for the sake of a back of the envelop calculation
a) that TSI noise is white,
b) the climate really obeys equation 1 above and
c) the time scale is 10 years and d) the sensitivity to 1 W/m^2 results in 0.5 C variations in temperature, so that 0.1 W/m^2 results in a 0.05 C variability in “equilibrium temperature” each month. (The equilibrium temperature would be defined as (S F ) in the equation above. (BTW: All numbers pulled out of the air.)
We would expect the variability in climate temperature due to this to be: σT= 0.05C sqrt(month) * sqrt(10 years*12 months/year)/2 = you do the math…
Of course, it’s unlikely the TSI noise is white (and I haven’t looked). But, I think you can see that if the the external forcing is noisy and the time constant is long, we can get big excursion due to the response to white noise. Because of the nature of the system, these temperature excursions would not look white, they would look red. (If the external forcing isn’t white, the temperature excursions wouldn’t look red either.)
Note that these excursions are not due to the non-linear nature of the navier stokes! (Also, note that I need to down load files and see what the order of magnitude of different things are.)
steven mosher (Comment#2715) May 11th, 2008 at 10:08 am
Lumpy forecasts.
you have solar forcing data, so extrapolate ( add the 11 year cycle for kicks)
You have a volcanic forcing, average the past 100 years and straightline it forward.
Or model volacano as a poisson process and monte carlo
You have aerosols et al? same thing, straightline it.
GHG? you have the forcing curve right? Call that Business as Usual. Extrapolate.
Rough order guess.
Here would be the neat thing. use lumpy and some forcing estimate to predict the trend
we will see from 2001 to 2011.
Then offer Annan a bet
lucia (Comment#2716) May 11th, 2008 at 11:39 am
steven– If you think I’m going to bet on Lumpy, you have more confidence in Lumpy that I! But, I guess if I were going to bet, I could imitate all the climate bloggers who try to “prove” by betting, make up a ridiculously one sided bet and insist the other side would take up the challenge if only they had confidence in their own models.
On the extrapolation: It’s easy enough to do. I have annual average forcings for some years forward and can make up ’scenarios’ further into the future. I can come up with all sorts of fake error bars. Real error bars would require me to figure out the real uncertainty in forcing, real uncertainty in my estimates for time constant and sensitivity yada, yada. Obviously, I’m not going to have real uncertainty intervals next week. (Or possibly any time soon. That’s one of the reasons I want to fiddle with the Schwartz idea with the white noise.)
But, I have a question: Do I get to do this the IPCC way and just provide my projections ‘relative to 2000′ without stating what the temperature in 2001 was? This would permitting me to slide stuff around afterwards, match the future and hide the fact that, based on the brand new definition of the temperature in 2000, the hindcast no longer fits.
Or do you insist I make projections based on HadCrut’s defition of anomaly? And when my projections are totally wrong, do I get to say my uncertainty intervals don’t include ‘weather noise’?
Also, how many quatloos do I get if this fits 1999-2007 within how much? ( There’s a good chance it will be more or less ok. I did this in the past for annual average data, but that resulted in somewhat different values. At it happens, I haven’t looked at how the monthly data version of “Lumpy” fits post 1999 data — amazing as that may sound! )
fred (Comment#2720) May 12th, 2008 at 1:29 am
Jane Austen too was burdened by domestic chores. She wrote her novels all the way through, then rewrote them again, until satisfied. She wrote in the main living room. They were done on small scraps of paper which she would cover with sewing or other papers when interrupted by people coming in. But at least when stopped in mid sentence she could reflect they were not yet published for everyone to read!
pliny (Comment#2723) May 12th, 2008 at 6:05 am
Lucia,
Lumpy looks great - it certainly gets most of what ModelE gets.
My guesstimate - time constant 7 years, sensitivity 0.5 C/(W/m2)
Nick Stokes
steven mosher (Comment#2725) May 12th, 2008 at 7:25 am
gavin has a great post up on internal variability, have a look
lucia (Comment#2728) May 12th, 2008 at 8:57 am
Steve– Thanks for letting me know. I’ll be writing a response.
Francois O (Comment#2734) May 12th, 2008 at 4:46 pm
Lucia,
Great job! In a similar thread at CA, I suggested that one should also make runs covering the range of uncertainty in the forcings. You’d get the error bars from that. Aerosols, for example, have a huge uncertainty before 1990, and even worse before 1980. In Hansen’s model, the aerosol forcing is literally tuned to reflect the 1940-1975 cooling, but without any observational data to support it. So convenient! Maybe you could run without aerosol forcing just to see. And as you know, Leif Svalgaard has his own version of TSI, and gets angry when someone doesn’t use it.
About natural variability: maybe what you could do is figure out how much noise, and of what kind you need in the forcing to account for the observed variance. Hey, maybe that could be a way to identify which forcing is important: find which of the forcings has the right noise characteristics!
My opinion was that a 20-line program could do just as well as Model-E. You’ve obviously proven me right!
david (Comment#2737) May 12th, 2008 at 6:22 pm
lucia,
I think this result is much more significant than you think. What you have shown is that model-E is in fact a linear system. This means that to model the response (eg temperature) to an input (e CO2 levels) all you need to do is evaluate a convolution integral, which you can do at home with a calcuator. No more super computers needed. This will save $billions. I’ll get back with more details if you are interested.
David
david (Comment#2738) May 12th, 2008 at 11:20 pm
lucia,
Obviously the climate isn’t a strictly linear system, as there are lots of complex nonlinear processes involved. But complex systems can sometimes be approximated quite adequately as a linear system. The example I know about is the distribution of drugs in the body. A mathematical model of the processes involved, expressed in fundamental terms, would involve temperature and pressure gradients, fluid flow (laminar and turbulent depending where), diffusion … and complex boundary conditions (time dependent). The physical scale this needs to be modelled on would require something like a micron sized grid. So I think it is physical system of comparable complexity to climate.
However, if your interest is in predicting a response to a known input then in most cases (when the system is linear, which is usually the case in spite of the potential nonlinearities in the underlying processes)there is a much easier way. A linear system’s response P(t) to an input R(t) is P=Q*R where Q is a function characterising the system (generally called the unit impulse response) and * is the convolution operation. What I think you have shown is that Model-E has a unit impulse function (ie Q) which can be adequately approximated by a simple two parameter exponential function.
David
Sean Houlihane (Comment#2739) May 13th, 2008 at 2:11 am
Re your caption for fig.2 - I believe this chart of TSI is corrected to 1AU, i.e. it excludes the annual 90W/m2 variation. Whilst it may be valid to use this in the context of a linear model, I think your caption is misleading, since it is the TSI at earth which ought to be averaged once you are considering timescales of less than years. I guess this will affect the noise analysis.
lucia (Comment#2740) May 13th, 2008 at 5:26 am
Sean–
Yes. The annual variation has been removed from the TSI. That’s likely appropriate for the sorts of analysis everyone does in climate modeling because the annual variation in temperaure is removed from the temperaure anomalies. If the average temperaure obeys a linear process- as in the equation above — removing the annual variations works out ok mathematically. (Of course…. if it’s no linear, it does not. And then, that would be a big open question.)
David–
Yes. I know I can do this in phase space. But with the weird forcings for the earth, it’s easier to do this in the time domain!
Demesure (Comment#2748) May 13th, 2008 at 4:03 pm
Wow,
With Lumpy, you’ll be saving millions of tons of carbon wasted in super computing.
You deserve a (timeshared) Nobel peace prize, Lucia.
lucia (Comment#2749) May 13th, 2008 at 4:04 pm
Demesure–
Heh! Rest assured, Lumpy has problems. . .
Still, I’m going to project with Lumpy just for kicks.
Julian Flood (Comment#2755) May 14th, 2008 at 2:23 am
Re reply to basil-
quote No. Neither Lumpy nor model E do well during 1910-1940 unquote
How do they both match if one abandons the Folland and Parker bucket correction? The hadcrut temp in your graph above would have a lower slope from 1910 to 1940. Miss out the next five years and you’re back with a nice steady slope which, by eye, matches. Ish.
http://www.climateaudit.org/in.....ct3b44e627
JF
david (Comment#2796) May 14th, 2008 at 10:11 pm
lucia,
No, the convolution integral is in the time-domain. That is, P = Q*R = integral from 0 to t of {Q(t-s)R(s)ds}. The way predictions are done in the applications I know about is that a known R is introduced and P is measured. Then Q can be found by deconvolution techniques (various methods - normally in the time domain; you can try doing this in the Laplace domain, but these approaches generally do poorly in the presesnce of noise). The prediction is then done by convolution of the just calculated Q with the new R to give the predicted new P.
One of the problems of Lumpy is that as it’s a d.e. everyone will immediately interpret it in terms of an implied mechanism, and find it lacking. You can avoid that with a linear systems approach, because your upfront assumption is just linearity. Nothing that looks remotely like a mechanism. The “structure” that you have put into Lumpy emerges from the data and the deconvolution step, but you can proceed with a prediction without any attempt to establish what’s physically going on.
David
lucia (Comment#2804) May 15th, 2008 at 3:54 am
David–
I suspect your everyone doesn’t include people from continuum mechanics. In my area, people prefer mechanistic approaches to describing the behavior of physical systems. They find purely statistical ones unappealing!
The ODE is conservation of energy. The criticism, from a phenomenological standpoint is it may be over simplified. This could be true. But relaxingconservation of energy isn’t an improvement from a physical standpoint.
Also, I do integrate a function of the form you describe.
david (Comment#2870) May 16th, 2008 at 12:56 am
lucia,
Sorry, my everybody includes the kind of people who comment on climate modelling; in particular the kind of people who developed Model-E. Of course, a mechanistic approach is best when you have enough understanding to get it to work. But I thought that is what Model-E is about and that the general consensus among skeptics is that it isn’t close to working. But it seems that I’ve misunderstood your model. I assumed that it was a purely phenomenological model. The equation has the same form as a first-order chemical reaction with a varying input and it surprises me that temperature would follow the same time course, whatever the degree of simplification. Or are you saying something like dH/dt = -kH + R where H is heat content and k a constant, and then transforming assuming an average world heat capacity?
The systems approach I’m talking about is deterministic. All that’s involved is: P=W{R} where R is the input (forcing) to the system (world) represented by an operator W and the response (climate) is P. Model-E is a mechanistic approach to determine in detail the way W operates. The system approach tries to make less dramatic progress by trying to see what can be done with less information. For example, if W has the property W(R1+R2) = W(R1)+W(R2)then the operator can be expressed as a convolution integral and the prediction I discussed becomes possible. I agree that this is not specially interesting but when my standard of living is under extreme threat based on possibly more interesting models whose predictions I think are completely unreliable, I’ll cope with a bit of boredom.
David
Nick Stokes (Comment#2872) May 16th, 2008 at 2:55 am
David,
Lucia’s equation is a heat transfer equation, in fact, Newton’s Law of Cooling. It’s widely used, with the T multiplier called a heat transfer coefficient. And yes, your H version is equivalent, as you say.
I was at first attracted to your convolution idea, because it doesn’t require a single time constant, and would yield the sensitivity (what we’re really after) as the integral of the impulse response. And you can get it by dividing Laplace Transforms.
But then I remembered how dodgy the results of that could be. And the methods that are supposed to do better involve assumptions. And the difficulty of the finite time data range, and what to do about the ends. And having to integrate the impulse response to infinity, with slowly decaying tails. So I’m not sure.
lucia (Comment#2877) May 16th, 2008 at 5:56 am
David–
Yes. The equation is the form of a first order chemical reaction. It is also a simplified expression of conservation of energy.
If we know the forcing, the solution is a convolution integral of the form you describe. But, I’m restricting myself to a specific form that conserves energy. I could extend the model, and permit myself different forms– but I would always restrict myself to those that conserve energy.
So– it’s curve fitting, but restricted to the form that conserve energy.
Does that make sense?
david (Comment#2912) May 16th, 2008 at 6:40 pm
lucia:
Yes, that makes sense. My problem with your formulation has been with T, which I see as a secondary quantity which is not in its self subject to the kind of process that Lumpy implies (no law of conservartion of temperature that I’m aware of; BTW W(R1+R2) = W(R1)+W(R2) probably implies conservation of energy but I haven’t anything like an argument apart from “it’s obvious, isn’t it?”). But as a phlogiston devotee I can easily envision a bucket of H that leaks out at a rate proportional to H and if T is proportional to H as well all I can say is that I’m broadly supportive of proportionality. As for curve fitting the linear systems approach is just that but more than yours. By using a linear first order d.e. you constrain at the outset the unit impulse response to be an exponential function. The fact that Lumpy does quite a good job suggests that this doesn’t matter much, and that a linear systems approach would probably conclude the same. So from a mathematical (yes, you can express the solution of your d.e. as a convolution integral)or practical point of view it is much the same. All I am saying (so far) is that there is a strong defence (since Model-E has this property as an outcome, as far as I can see) against the criticism that Lumpy is a ridiculous oversimplification: that what you really have in mind is W(R1+R2) = W(R1) + W(R2) but that you thought it more congenial to put it in a form familiar to mechanics.
David
Nick Stokes (Comment#2913) May 16th, 2008 at 7:51 pm
David (and Lucia),
Actually, Lumpy is not a ridiculous over-simplification, and has very respectable ancestry. See for example, this 2005 Nature article, especially Box 1.
david (Comment#2914) May 16th, 2008 at 8:30 pm
Nick,
You are right, there are no ends in my equation. With P=Po at t=0 it should be P(t) = PoQ(t) + Q*R. But no particular problem (in the time domain) as t->infinity, as Q(t) -> 0 (and since we are all going to fry soon who’s interested in such an egotistical time frame anyway?). Yes, just a ratio of Laplace transforms (p(s)/r(s))gives the transform of the unit impulse response (q(s)), but that’s not the way to go in practice. The problem is that in the real world there is a +e(t) to be added and the Laplace transform of an unknown high frequency function of time is unknowable (and possibly very significant in the s-domain even if it is insignificant in the t-domain).
A method I’ve used with noisy data is to assume a functional form for Q(t), say F(par,t) where par stands for a vector of parameters (eg based on lucia’s model, F(par,t) = par1.exp(-par2.t)). Then evaluate P(t) (numerically, in general) at times ti for which observations are available, to construct: sum over all i of {Pi-P(ti))^2} where Pi is the observation at ti and P(ti) is the model prediction. This is then a standard nonlinear regression problem and you get the vector par’ which minimises the usual sum of squared differences between the observations and model (and the usual statistical properties, such as sd of the parameter values). Then P’ = PoQ’ + Q’*R where Q’ = F(par’,t) to generate predictions of P (and a range of P’ if you want, based on the sd of par’).
The drawback of this is that you assume a functional form for Q, which you don’t know (there are methods that don’t require this but they are not always stable). If you look at the GCMs you would conclude that this is impossible but it isn’t hard to get an intuition about it by recalling that Q is what P would be with a single short blast of forcing (unit impulse input). So with lucia’s model, based on H, a single exponential function for Q is plausible (heats up at time zero, then cools down at a rate that decreases as the heat content decreases). In the end, it’s like any curve fitting issue: choose a bad Q and you won’t get a good fit.
The additional problem here is that we don’t know what R is either. Lucia is using the forcings used by Model-E so I think she has a model of Model-E which might or might not be a climate model.
David
Arthur Smith (Comment#2915) May 16th, 2008 at 9:43 pm
Hi guys - interesting discussion, but I’m not sure how “conservation of energy” fits in here? I mean, as a whole, the energy of the planet is expected to increase steadily with any given forcing, until the net forcing is brought back down to zero through higher surface temperatures. If the pseudo temperature T increases, then total energy of the system has surely increased. What modifications to the relevant equations would not “conserve energy” to some larger degree than this? Or are you meaning something by the word “energy” other than what I would take it to mean? Thanks!
david (Comment#2916) May 17th, 2008 at 1:02 am
Nick,
Yes, same model, different application. I’m surprised that anyone from the Hadley Centre would be allowed to use such a simple model. I thought they would insist that only a GCM was up to the task, so that no-one else could check the results.
Lucia,
Have this paper or references 1-3 been discussed on climateaudit? They appear to be about how the sensitivity figures were obtained, which has been discussed there. (I can read but not write on CA)
David
david (Comment#2917) May 17th, 2008 at 1:21 am
Arthur,
Rate of change of heat content = rate in - rate out,
so it heats up if rate in> rate out, cools if rate out> rate in.
David
lucia (Comment#2938) May 18th, 2008 at 6:48 am
John A: Lumpy just is female. Real women have curves.
steven mosher (Comment#2706) May 11th, 2008 at 6:24 am
thanks that is great. you go lumpy! also note that modelE is the average of 5 runs