Schwartz & Scafetta Estimate Climate Time Scale
One sign of an important paper is to count the number of immediate responses. By that measure, Schwartz 2007 “Heat capacity, time constant and sensitivity of the Earth’s climate system ” appears to have made quite an impact.
The impact is not surprising as this paper represents the relatively rare attempts to estimate climate sensitivity empirically rather than using models. The empirical value turned out to be much lower than generally predicted by climate models.
It appears that three comments on Schwartz’s paper will soon appear in J. Geophys Res. along with Steven Schwartz’s has also published his response. The various papers are a) Schwartz’s response, Scafetta’s comment, Foster, Annan, Schmidt and Mann and Knutti (for which I have no link). I’ve looked at all but Knutti.
Since these things are likely to be discussed, and I have posted my initial thoughts in comments, I thought I’d comment on Scafetta’s very good paper, and Schwartz’s response to Scafetta’s comment. Later I will comment on FASM (though I have already said much of what I think in comments on Friday’s blog post.)
Comments on Scafetta
Scafetta’s approach to Schwartz’s paper was to use the idea in the paper as a springboard for further analysis. In some ways, I see the initial Schwartz paper as “the great idea”, that idea being “Let’s try to estimate the climate time constant based on data. Here’s my first shot!” I see Scafetta’s paper as one that takes the great idea and improves method of applying the idea to get a more accurate answer.
This is the way things work in science. Mistakes may be made along the way, but we need to start out with ideas that permit further analysis.
What’s the great idea?
So, Schwartz’s idea was that we could treat the earth’s climate as a simple “lumped parameter” with one time constant and futher assumed that this “climate lump” was being heated at a rate that, on average, increased linearly with time, but also varied randomly. The random variation in heating was assumed to be “white noise”.
Under this assumption, Schwartz concluded the random component of surface temperature of the planet would be something called “red noise”. Moreover importantly, if one detrend the average temperature, calculated an correlation, ρ of lagged residuals, the correlation, ρ would be a constant. Based on that constant, we could estimate the time constant of the planet.
So, basically, one does some calculations, plots, and read a number off the plot. Schwartz concluded the time constant of the climate was about 5 ±1 years.
But there were problems.. .
Scafetta’s paper noted that Schwartz’s analysis of the time series for Global Mean Surface Temperature is imperfect and might not give the correct result for the time constant for three reasons.
- One is that the autocorrelation (ρ) as a function of lag time in Schwartz 2007 is not a constant function of time, as would have been suggested if the temperature series obeyed the simplified equation for the earth’s climate system described in Schwartz.
- Scafetta also noted that assuming the underlying trend in temperature over time did not suggest a linear increase in temperature over time.
- Finally Scafetta also discussed the known biases that arise when estimating the correlation in lagged residuals from data.
To deal with these issues, Scafetta did three separate things:
- He treated the temperature series as having two time constants: one faster, one slower. This would make sense if, for example, we thought the air or ocean made up two systems. (Alternatively, one could image the air, ocean and deep ocean make up three interconnected systems, or break up the planet in more and more pieces.)
- He de-trended the temperature series two different ways. First, he assumed a linear trends, as did Schwartz. But, recognizing the shortcoming of this method, Scafetta also detrended using the estimate for the mean surface temperature as a function of time as predicted from GISS model results.
- He performed an analysis to estimate the possible bias the estimating the autocorrelation.
Scafetta’s results for the time constant are: a) τ1=0.4 ±0.1 years and τ2=8.7 ±2 years under linear detrending and b)) τ1=0.4 years and τ2=8.1 ±2 years detrending using GISS model data. Scafetta then considers the effect of bias in estimating the autocorrelation from short data samples, and concludes the bias could be 3 years. Therefor, accounting for bias the time scale for the climate may be 12 ± 2 years. (Note: the bias may be smaller than 3 years, as its magnitude depends on the magnitude of the actual time constant.)
All three extension by Scafetta appear to me to be improvements on Schwartz.
Afterwards, assuming Schwartz’s estimate for the heat capacity is correct, Scafetta estimates the time constant for the climate most consistent with the data is 1.7 C for a doubling of CO2. Accounting for the possibility of bias the sensitivity may be as large as 2.7C for a doubling of CO2.
How does Schwartz respond to Scafetta?
Schwartz seems to more or less endorse Scafetta’s observation.
Schwartz discusses the physical basis for the possibility of time scales of approximately months and accepts the two time scale system as both plausible and an improvement on Schwartz’s first analysis. After accepting this possibility, Schwartz then explain that, owing to the shortness of the 4 month time scale, we could assume the time scale varies as a function of lag as follows:
where ρ(Δt) is the autocorrelation measured as a function of lag time, δt, and ρo is the zero intercept at Δt = 0.
Rearranging this equation, and plotting as a function of lat time on semi-log paper, we would expect to see that at lag time large compared to 0.4 months, ρ should fall on a straight line. This is illustrated in figure 2a from Schwartz, reproduced below:
From Schwartz 2008: In the figure above, Schwartz shows series of data both before and after correcting for bias but uses a different method of bias correction than selected by Scafetta.
To obtain time constants, Schwartz plots the time constant τ calculated using equation (1) above as a function of lat time Δt. and then (I think) he eyeballs the data to obtain estimates of 8.8±2 years and 7.2±1.5 years.
Ultimately, Schwartz revises his estimate of the time constant to 8.5±2.5 years, which is noticeably larger than 5 ±1 year reported in the original paper.
How does this relate to my earlier blog post?
Those familiar with my first blog post discussing Schwartz will recall I observed the even if the earth’s climate obeyed the equation described by Schwartz, the data should obey the relationship discussed above and performed a linear fit on ln(ρ) vs time. I obtained 7.2 ± 2 years based on annual average GISS land ocean data and 8± 2 years based on annual average MetStation data. These magnitude are comparable to the values obtained by Schwartz– and well should be as the two methods are equivalent.
However, I did not account for bias in anyway. (This is because I was entirely unaware of the issue.)
Because I seem compelled to always look at data myself, I will be examining different methods of accounting for bias to see how they might affect my estimates and my uncertainty intervals.
However, as my analysis currently stands my estimates are comparable to those of Scafetta and Schwartz, and were obtained entirely independently. I’m reasonably sure neither read my blog back in December, (when the blog was experiencing 2 visits a day, mostly from me.) I had certainly not read their papers, which were posted on the web quite recently. So, it’s clear that people with different backgrounds are all making similar (though slightly different) arguments and coming forward with similar results.
What does a time scale of 8 years mean?
First, this estimate for the time constant of the earth and the climate sensitivity is much lower than for models. So, there is an inconsistency here.
All other things being equal, the lower time constants means:
- Less “heat in the pipeline” than predicted by models.
- Possibly lower total climate sensitivity– if we also believe Schwartz’ estimate of the climate heat capacity.
Obviously, the method used by Schwartz, and Scafetta involve a great deal of approximation, but approximation is common in science. However, despite the huge amount of simplification, when the very simple model in Schwartz first paper is corrected by recognizing either noise or very short time scale weather noise, the simple model seems to give fairly decent representations of the major features autocorrelation in lagged residuals for measurements of GMST.
Given the reasonable agreement with data, it is difficult to wave away the estimate for 8 years as simply “wrong”.
So, what are we to make of the mismatch between models and this estimate?
Steven Schwartz discussed four possible explanations for the mismatch: a) The 8 year estimate could be in error due to uncertainties in the measured data, b)The time record may be too short to use this method to estimate the time constant, c) The method used to infer the time constant could be inappropriate and, d) Climate models may be inaccurate.
Figuring out which of the four applies will require thought and work. Obviously, researchers will continue to strive to resolve the mis-match between models estimate for the time constant and Schwartz estimate. In the meantime, it appears a time constant of 8 years falls in the range of “plausible” based on 125 years of measurements of the Earth’s mean surface temperature.
Previous Post:
« Chileans Flee Volcanic Eruption.
Next Post:
What is the Climate Time Constant? Refining the Estimate I. »
54 Responses to “Schwartz & Scafetta Estimate Climate Time Scale”
You can leave a response, or trackback from your own site.
Phil. (13 comments.) May 7th, 2008 at 2:54 pm
Arthur, somewhere recently on one of the blogs someone showed that there was an ~1 year peak in the FT of the difference between UAH and RSS.
Andrew May 7th, 2008 at 5:00 pm
Scafetta’s use of a model to “detrend” the temperature record seemed to strange to me given his and Bruce West’s earlier work in which they concluded that models (like the GISS model) may underestimate the contribution of solar to 20th century warming.
pliny May 7th, 2008 at 6:14 pm
apsmith, yes, I think FFT’s are the best approach. But shouldn’t you be looking for peaks in the complex Fourier domain? These would correspond to the decay time constants that SS and co are looking for. The imaginary part would be the inverse time constant.
pliny May 7th, 2008 at 9:50 pm
Lucia, I proposed using synthetic data to test the algorithm for recovering relaxation time constants. I see now that it has been done. Tamino, at RC, did it even before the SS paper had appeared in JGR. He used timelags of 5 years, and shows that the autocorrelations look quite different to those from real temp data, but also the recovered relaxation times are substantially less than five years. In this case it reflects bias; I think that at 15 years, the effect of detrending would add significantly.
The comments to appear that you link to, by Foster et al (not unrelated, I know) have a more thorough working with AR(1) in their section 5. They show, for example, that if you put in a relaxation time of 30 years, you get back about 13 years. Detrending would hurt there.
The SS paper has been criticised for poorly representing the atmosphere. However, these criticisms do not relate to that; they say that even if the temperature really was AR(1), the analysis would give a badly wrong answer.
On Scafetta, I’ve realised what is wrong with adding a second AR(1). It isn’t initially obvious why sensitivity should be related to the relaxation time. The reason is explained (not well) in SS eqns (4-7), or better by Foster et al in their eq (4):
C dT/dt = F − λT + e
To relate the constants, you need both the small-time (initial) and large-time limits of this differential equation. The problem with Scafetta’s model is that the temperature with two AR(1)s does not satisfy this equation, or any simple modification. You can’t say that the short-term process can be ignored in the long-term limit, because the initial condition is also needed.
D.R. Williams May 8th, 2008 at 12:02 am
The second paragraph from above:
“The impact is not surprising as this paper represents the relatively rare attempts to estimate climate sensitivity empirically rather than using models. The empirical value turned out to be much lower than generally predicted by climate models.”
is followed later by:
“Steven Schwartz discussed four possible explanations for the mismatch: a) The 8 year estimate could be in error due to uncertainties in the measured data, b)The time record may be too short to use this method to estimate the time constant, c) The method used to infer the time constant could be inappropriate and, d) Climate models may be inaccurate.”
Surely part of the validation of the climate models is to test the validity of the prediction against the data?
Andrew May 8th, 2008 at 12:41 am
D. R. Williams, if Schwartz wants to stay in good graces with th modelers (and I’d say he was on thin ice already) the last thing he would suggest is that the reason for his empirical results mismatching with the models has anything to do the ~models~ being wrong. And so he quite literally suggested this last. I know I should be more disturbed by this, but I’ve come to expect this sort of thing. Consider it astounding that he acknowledged that possibility at all!
lucia May 8th, 2008 at 5:17 am
Pliny–
Using the synthetic data is a good idea.
However, the results Tamino got are wrong.
I addressed Tamino’s post almost as soon as he posted it. His analysis comes to the wrong conclusion because he makes a classic mistake (as does every freaking peer reviewed paper on this subject).
Measured data is never precisely the same as the underlying thing. It always includes noise.
If you assume the measurement, T, is the sum of “real” GMST, Θ and white noise, u, you find
ρT=ρΘ * { 1/ (1+<uu>/<ΘΘ>) }
Many people forget about measurement noise before launching off and doing analyses, but it’s always there!
When you add the noise, the data look almost precisely the way they should if they follow the IPCC equation. (Arthur’s looking at some peaks in the spectrum– but the peak Tamino saw was due to another feature– integration.)
I’m repeating those steps but adding the noise. I’ll be explaining more about the noise and the integration when I have the first “chunk” done. (Then, I’m hoping Arthur will tell me where he finds robust peaks!)
Cassanders May 8th, 2008 at 5:31 am
@Lucia,
Again from a numerically challenged ecologist.
When doing measurements I am used to deal with concepts like precicion (i.e. the inherent variability between replicates) and accuracy ( the trueness of the measurement WRT to the “real” value”), and the eventual distributions of these.
Could you elaborate a bit on “white noise”? How is it related to accuracy and precicion (if related at all)? Is it an amalgamate of the two?
Does the *white* signify that it is distribution-free?.
Cassanders
In Cod we trust
Douglas Hoyt May 8th, 2008 at 6:08 am
What are the results if you do the analysis for the Northern Hemisphere alone and for the Southern Hemisphere alone?
lucia May 8th, 2008 at 6:19 am
Cassanders,
With regard to individual temperature measurements, white noise mostly relates to precision. If we examined a daily temperature measurement, and for some reason, we expected the imprecision of the individual measurement to have a standard error of 0.1 C, then 0.1 C would be the magnitude of the standard error of the white noise.
But, the 0.1 C only tells you how big the white noise is. The feature that makes it white is that the error from 1 measurement is uncorrelated from the any other measurement. This happens, when, the errors on say “day 2″ have nothing to do with the errors made on “day 1″.
In contrast, if the errors on “day 2″ are somehow related to the errors on “day 1″, then the noise will not be “white”. It might be red, pink, blue, 2nd order blah, blah, blah. (For example, if you sample temperature once a second using a thermocouple with a 5 minute response time, the errors will look “red”. This is because it takes the thermometer 5 minutes to respond to a change in temperature. However, if you sample every 24 hours using this thermocouple, the noise will look “white”, because the error due to the response time is tiny and all you are left with is error due to other factors (like calibration error, electrical noise etc.)
Most decent experimentalists do their best to make sure the random component of errors in the data are white. This makes it much easier to deal with data for most down stream processing purposes.
Because most experimentalists do this, many downstream users then forget that even white noise doesn’t always average out. It averages out when you are calculating the mean. But, for example, in turbulence, the experimentalist is actually trying to measure standard deviations and correlations. When measuring these, you must never forget that the white noise doesn’t average out. You need to know how big it is and subtract it. (Better yet, make sure it’s small enough to not matter.)
lucia May 8th, 2008 at 6:21 am
Doug– I don’t know the answer to that question. I think Schwartz compared in his first paper and found they were similar in magnitude. However, I’m not sure.
pliny May 8th, 2008 at 6:24 am
Lucia,
Thanks to the pointers to your previous comments on Tamino. I’ve been reading them, and also the post by UC at CA.
However, the part of Tamino’s work that I was noting was where he generated AR(1) data with a 5 year lag, put it through the analysis, and came up with a shorter relaxation time. I can’t see that that would be affected by the noise of experimental measurement. It’s just a statistical processing issue. Incidentally, I note that UC seems to have done the same thing, but with sometimes longer relaxation times as a result. Neither T nor UC where very explicit about what they did at this stage (though UC gave R code for the GISS treatment).
It’s actually a bit unclear to me what you are suggesting. Is it that data should be synthesized with an AR(1) process, then have white noise added? Is it the difference from the pink noise of the AR(1) that matters?
Anyway, I’ll keep trying to get through your Dec post, UC’s posts and the CA discussion. Incidentally I noticed that Francois O and UC both brought up the detrending issue, and I haven’t seen a real answer to it.
lucia May 8th, 2008 at 6:57 am
Pliny–
I know what I’m saying is unclear– because I plan to document in more detail after I run the synthetic data with the white noise.
However, what I am saying is:
a) there is a bias.
b) the adding the white noise to the synthetic process may make a difference in the magnitude of the bias. (I don’t know how much though, but it’s not zero.)
d) adding white noise to the synthetic process may make a difference to the standard deviation in the population of time constant determined experimentally. (I don’t know if it makes it larger or smaller, and in particular, I don’t know if this will matter to my method or Schwartz’s.)
e) Tamino (FASM) didn’t determine the magnitude of bias or standard deviation with white noise added to the synthetic process process.
So, while there is a bias and there is a large standard devaition, I don’t currently know how large the bias or the standard deviation is/are!
I want to run my method with a bias correction. Also, for my method, I need to know the bias over a large range of ρ because I apply a linear fit to the ρ vs lag time data, so the information in FASM isn’t enough for me to go on. (And given what Schwartz said in his response, it appears there is no standard method. So, as long as I’m doing the synthetic data thing, I’m just going to estimate the biase with the synthetic data!)
I’m going to deal with detrending issues after looking at the problem with no detrending. I’ll look at the issue Arthur addresssed then too. The reason for this is that the bias/scatter issue exists both with and without detrending, and with/without the “spikes” in the spectrum. So, to deal with all of them, I need to do this first.
I’ll probably be posting today.
Arthur Smith May 8th, 2008 at 7:39 am
I’m feeling still quite new to this empirical temperature analysis business, so please forgive me if this is an obvious question - I’m also re-reading Schwartz’s 2007 paper to try to understand it all better. If somebody here has better references on the subject I’d love to read them! But anyway, here’s my question:
Isn’t the whole point of the relaxation time (tau) to be a measure of the response of temperature to forcing? I.e. don’t you somehow need to separate out the driving causes of temperature change from the temperature response? So how can we possibly hope to extract it from the temperature time series alone? I feel I’m missing some key insight here!
lucia May 8th, 2008 at 7:55 am
Yep!!! In the full problem, you need to know the spectral properties of both the random component and the mean component.
Schwartz gets around the need to know the forcing using the classic method: Making an assumption!
Schwartz assumes the random component of the external forcing for the temperature anomaly is white. This is certainly not exact, and may or may not be a good approximation. In principle, the definition of “anomaly” takes care of the annual periodicity (and this can be shown.) However, if the random component of the external component is not spectrally white (or nearly so), this will introduce errors.
In his Schwartz’s first paper, and in his reply to comments, Schwartz assumed the mean component is a linear ramp function that’s been going on “forever”. That’s the assumption that permits linear detrending.
So, Schwartz made two assumptions to deal with the exact issue you are pondering. As with all assumptions, either one or the other could be sufficiently wrong to result in inaccurate answers. (This is why showing the autocorrelation looks “right” is an important fiduciary step. If the data looked self-inconsistent with the assumptions in any way, that would cast serious doubt on the method.)
Scafetta tried to account for the non-linear nature of the mean component by by detrending using GISS model mean data but also did a linear detrend. This made very little difference to the answer. (This likely means the amount of ‘false noise’ introduce by detrending is relatively small compared to the amount of real, honest to-goodness random “noise” in the extrenal forcing to the system. )
Arthur Smith May 8th, 2008 at 9:09 am
Lucia - that makes sense - but the variations in the temperatures, with or without detrending, which is just the simple thing I’ve been looking at (and I thought that’s what Schwartz etc were looking at), are definitely not “white”! The periodic variations are surely from periodic forcings, whether external (year, solar cycle?) or internal oscillations in the climate system - so do you have to also fit the magnitudes of those, and subtract them out, to properly do Schwartz’s method? I notice Schwartz 2007 does mention using “deseasonalized” monthly temperature records, which I don’t think he explains - that would presumably eliminate the yearly cycles.
Is there a theoretical basis for the “white noise” assumption? I don’t see the term “white noise” in Schwartz 2007 at all…
lucia May 8th, 2008 at 9:48 am
Arthur — Schwartz cites underlying references for the white noise assumption.
The physical reasons I can think of are imperfect, but here the are:
1) If you examine the TSI charts, you” see loads of high frequency noise. That might tend towards “white”.
2) Over the course of the year, we have somewhat random things happen, which are not related specifically to weather. These include different people spewing out more or less particulates, those particulates falling on snow in different places, volcanoes erupting etc.
3) Internal behavior of the climate could add more noise to the “f” term.
So, collectively a bunch of things might tend to result in the noise being “white” or approximately so. (Or, the assumption could be bad.)
BTW. I’m adjusting what I did for the bias– still assuming the white noise is ok. I’m doing it in this order not because the stuff you are suggesting is irrelevant, but because I need to have a clue how to deal with bias in my method first.
I’m up to 15 years or so. I need to proof read though.
What is the Climate Time Constant? Refining the Estimate I. | The Blackboard May 8th, 2008 at 12:33 pm
[...] Comments: Schwartz & … [...]
Martin Ringo May 8th, 2008 at 2:58 pm
I have a question (and I apologize to all you who understand the physics of the Schwartz and respondents better than I — which should be just about everybody).
Everybody appears to be working with the empirical AutoCorrelation Function (ACF). That is the estimate of the autocorrelations at each lag from the detrended time series. Coming from a different discipline, economics, I almost never work with the ACF. Rather I work with the estimated ARMA model. Then once I have an estimate of the ARMA (that is stationary), I calculate the theoretical ACF from the ARMA coefficients. The reason for this is that just about all the coefficients of the empirical ACFs are full of the small but pervasion random autocorrelation one gets with the ACF of a purely random, non-autocorrelated series. (Look at the probability of the Q-statistic move up and down for the ACF of white noise, and you will get an idea of the problem.)
Further, if you are fitting a pure AR(…) model, i.e. no MA terms, then don’t you get the time constant (or is it time constants when there are more than one AR terms?) directly from the AR coefficient?
Any help (with explicit equations if possible) would be appreciated.
Thanks, Marty
lucia May 8th, 2008 at 3:15 pm
Martin-
I think the problem is that the physical scientists want to connect each term in the fit to a physical meaning or at least explain where they come from. So, we would either need some explanation of what the other terms in the ARMA mean, physically, or an explanation of why they are required to get need to get rid of some known feature that comes up in the process of handling data.
So, under Schwartz’s model, the time constant of the climate is claimed to have a physical basis and be meaningful. I can justify white noise based on measurement imprecision.
The time constant comes from saying the climate obeys
dT/dt = -T/τ + α F where F is forcing and τ is the time constant and α is just another parameter.
That’s basically conservation of energy.
The idea that we can just pull that out of the series the precise way we are arises if we say F is white noise. This because an AR(1) process if and only if F is white noise. Otherwise, it’s a different process.
But I can’t connect higher order terms to the physics. Or I could come up with a physical basis for higher order terms, but I then need to say that I believe the “forcing” to be red noise or something that creates the other terms in the ARMA. (And you can see Arthur is arguing F is probably not white because T isn’t red.)
Now, mind you, the forcing isn’t white. But scientiests will often accept the possiblity of white noise because it’s the simplest thing you can pick. But if you start advancing other types of noise, they are going to demand an physically sound reason for that particular noise being appropriate for the Forcing, F.
And the problem with that to get scientists to accept other terms in F, we need to explain what they might be. I can explain adding white noise to the Temperatures– that measurement uncertainty.
On the other hand, if there is a good statistical explanation for why the terms “just arise” in a string of data even when the real underlying process for temperature, T, is red, then that can fly. (Sort of like, we use N-1 in the denominator of a standard deviation because of what happens when we estimate the variance using the sample mean.)
I just don’t know the reason the other terms in the ARMA yet. (Doesn’t mean there isn’t one. I just don’t understand the physics of math.)
pliny May 8th, 2008 at 3:38 pm
Arthur,
I had trouble with that sensitivity-relaxation connection too, and Schwartz explains it poorly. But here is how I understand it. It follows from the linear system assumption. Suppose that initially T=0, and you add a new steady flux F. Then, with C=heat capacity:
dT/dt = a + bT (call this Eq 1 (eq 4 in FASM), it’s the general linear assumption, with a,b to be found)
Initially T=0, so you’re just steadily heating a body, so a=F/C
and with the assumption that the response to any perturbation has relaxation time τ, then b=-1/τ
therefore, when finally T has settled and dT/dt=0, then F/C=T/τ, or sensitivity T/F = τ/C
Note that it is essential that, for this logic, T actually satisfies Eq 1. That is what is quite wrong, I think, with Scafetta’s mixing of two AR(1) processes. You lose this essential requirement.
D.R. Williams May 8th, 2008 at 4:44 pm
Andrew, my question:
“Surely part of the validation of the climate models is to test the validity of the prediction against the data?”
is directed at the modeling process, not the research that provides an estimate by examining actual data.
I’m a bit mystified at the lack of such an elementary notion as model validation. If the models do indeed predict a single time constant, then the empirical methods that do likewise should not be subject to examination and discussion as to whether one, two, or more time constants.
Indeed, it seems that the cart has been put before the horse. Why was the empirical data fit not done in 1993? Why no previous discussion about the number of time constants necessary to describe the system. Would it not be appropriate to do the analysis and have the discussion before building a model? And are these model predictions themselves an end-product, or are they used as inputs elsewhere?
D.R. Williams May 8th, 2008 at 4:47 pm
sorry, a bit seems to have been lost in editing:
…or more time constants are needed.
Martin Ringo May 8th, 2008 at 4:51 pm
Lucia,
Thanks, I halfway there. And now let me cause more problems, albeit the same idea.
The physical basis for an AR(1) seems to me to me the same as that for an AR(3) process. (Let me be clear on notation here: AR(1) means a first order autoregressive process. AR1 is the first order coefficient of any autoregressive process, i.e. an AR(1), AR(6) and an ARMA(3,3) all have an AR1. A higher order AR processes, say an AR(24), might have only a few non-zero coefficients. I do not believe there is any standard notation, but I use AR(1, 3, 12, 24) to denote a 24th order AR process with non-zero coefficients at lags 1, 3, 12, and 24.) An AR(1) process has geometric decay every period upon which the processes is defined, e.g. on monthly data. An AR(4) process has geometric decay every 4th period. This is elementary. But with regard to the physical world the time unit is arbitrary. Why couldn’t the fundamental unit of time be gamma0*PI*seconds? Or maybe, within our conventional time scales, why isn’t the decay nested, i.e. a decay of a decaying process? I haven’t done any work with such series, but I suspect that they have some distinct differences from a pure AR process.
Or let me get a bit nuttier here, suppose our data were continuous (we can always fit the N observations perfectly with an N-1st degree polynomial in time). Now we don’t have to look at the “delta t” — sorry I am not certain of how to paste a Greek delta — as an interval that can continuously vary instead of varying only for the natural numbers. Of course as the nuclear physics used to say to me “Die ganze Zahl schuf der liebe Gott, alles Ubrige ist Menschenwerk.” (But they were just quote Leopold Kronecker, the tormentor of Cantor.) But maybe the positive half of the real line is a bit approximation for time than various natural numbers, and then what is fundamental delta? Not clear. We could solve for the delta which gives the largest first order autocorrelation subject to being greater than 2*(1/sqrt(N)) and say that best fits the exponential decay model. Or maybe better is keep varying delta then create the discrete time series and test the AR(1) model versus the AR(K) or ARMA(K,P) model in a likelihood ratio test. If we get a clear cluster of AR(1) non-rejections, we have a division of time consistent with an AR(1) view of the world. But right now, whether it is monthly or annual, Schwartz is fitting an AR(1) shaped model into a non-AR(1) shaped data set. (Oh, in case, I hadn’t mentioned it: the serial correlation Lagrangian Multiplier test — the Breusch-Godfrey Test — fails for lag inclusion of 3 through 12 with annual data, and an ARMA(1,1) model can’t be rejected in favor of an AR(1) one.)
I guess if wanted to do that, I think the way to proceed is to fit a trend plus AR(K?) or ARMA(K?,P?), then take the AR(1) process out of the fit and subtract the remaining part from the data. Then you have — in theory — an AR(1) process plus the residual from the “true” model. The AR1 then gives you your tau: -1 / ln(AR1). End of story. No empirical correlograms or what. Just an AR(1) simplification and an estimate consistent with that simplification.
pliny May 8th, 2008 at 5:24 pm
Martin,
I think you’re missing the primary purpose of this activity. It isn’t to find the best model for the time series. Schwartz’ model is a first order linear differential equation with constant coefficients - as simple as you could get. If you allow all those simplifications, the equation gives you a single time constant, and relates it to the climate sensitivity. The AR(1) model is a way of identifying this time constant. Other time series models don’t model this particular differential equation, and you lose the link to climate sensitivity, which is what it is all about.
Martin Ringo May 8th, 2008 at 11:03 pm
Pliny,
You may right in that I am missing the point. But I am trying. The Schwartz article is an empirical estimate at the time constant, and what I am looking at is how he is pulling a number out of the data. So let me try saying my point another way. Look at Schwartz Figures 5 and 6. I’m saying that those patterns are AR(1). We probably all are willing to accept that whatever is generating those patterns has a positive AR1 coefficient and I am not questioning the underlying model, but I am questioning the way Schwartz goes about pulling his tau — which is just -1/ln(AR1) if the model fit well — from the data. If you say OK I am going to force the AR(1) model on the detrended data, then why not do exactly that: force it, make an estimate of the AR1, which is around 0.6 for the annual data and gives a tau of about 2 years. But taking the correlogram of the series is neither consistent or inconsistent with the underlying theory, it is — and this only my opinion, albeit a somewhat informed opinion — bad statistics.
If you build yourself a little model of the DGP, you will see if you use an AR(1,4,10,14) model with the fitted coefficients, you will generate ACF patterns similar to Schwartz’s on a regular basis. But if you stick in the AR(1) model coefficients, you will usually get a much lower tau and seldom get the pattern. Thus, if you are using the correlogram from the data, you are implicitly contradicting the AR(1) assumption on the DGP. This is what you tell me not to sweat. OK, but by using the correlogram from the data you will also be estimating something other that the time constant of an AR(1) process. Hence, I believe that it is necessary to explicitly go about extracting, estimating or whatever the AR1 value underlying the model.
pliny May 9th, 2008 at 12:01 am
Martin, I’m sure you’re right that much can be done to improve the statistical identification of the time constant. My worry is that we have to be able to then to go to the differential equation, of which the AR(1) process is a kind of discretisation, and say that the time constant determines the climate sensitivity. If a more elaborate model was fitted, we may have a better estimate of the coefficient, but with model assumptions that are inappropriate for the de. In effect, it is a better estimate for a different number. This is why I was suggesting that you can’t divorce the statistical analysis from the subsequent differential equation interpretation.
I like your idea of forcing the AR(1) model and estimating AR1, because the error could then give an indication of how good is the de model relating sensitivity to AR1.
Martin Ringo May 10th, 2008 at 12:46 am
Pliny,
Let me ask a question here which refers to Lucia’s next post (”What is the Climate Time Constant? Refining the Estimate I”) and to her previous post on estimating the time constant. How should one force the AR(1) model on the data?
As a non-physical scientist, my instincts are to just go ahead and estimate the AR(1) model temp(detrended) = Constant + B*lag1(temp(detrended)). No need to even say AR1 in the model since there are no other explanatory variables. For the GISS data, the AR1, i.e. the coefficient B of the lagged, detrended temperature, for annual data is about 0.63. For the monthly it is about 0.70.
But I could go to the ACF of the detrended data and do a nonlinear fit of the equation ACF(lag i) = B^(lag i) plus error of course. Now because I am using the geometric form I don’t have to worry about negative ACF values and can use the whole series. If do that, the AR1 for the annual goes up to 0.78 and the monthly goes up to 0.95.
In over the lowest to the highest that is a range on tau of 2 to 20. That range has nothing to do with the theory arguing for the estimation, just the question of what is the best AR1 representation of the data. As an econometrician I am going to say the annual AR(1) model has best “information” criterion although I really think the model is either some higher order AR process, an ARMA process, or maybe a fractionally differenced (long memory) ARIMA process (ARFIMA — AutorRegressive Fractionally Integrated Moving Average) — or maybe we just have the wrong units for the time scale altogether.
I should also note that when one uses the whole ACF the effect of higher order AR coefficients in the process affect ones estimate of the AR1 value unless they are accounted for explicitly. Thus, I would again come back to the direct estimation which comes back to the exclude explanatory variable bias, e.g. the AR1 coefficient in an annual ARMA(3,3) model is around 0.35, a tau of less than a year.
Now what I don’t know is suppose we have an AR(1) process. How often will we reject that process as the true model in comparison to some higher order ARMA model, AND when we do reject the AR(1), what is the best way to estimate the AR1 value? If the answer to the first part is a large percent, say 20-30%, then the answer to the second part, which is in the realm of statistics, with give us guidance. If the answer to the first part is a trivial percent, then it is back to the physics and the physical explanation of the higher order terms in the time series.
Marty
pliny May 10th, 2008 at 4:39 am
Martin,
I look at this more from a numerical DE perspective. I’ve put an electrical analogy on the latest thread. We have a linear inhomogeneous DE:
dT/dt-T/τ = H where we need to think of two cases for H:
1. H=-F/C, where F is a steady (extra GHG) heat flux, and C is heat capacity, and
2. H=w(t), some stochastic process which has generated the observed values.
Case (1) gives us the sensitivity. In the latter case, to get a value of τ out of what is observed, it is assumed that w(t) is a white noise process, or as nearly so as we can make it, by varying τ.
Discretised, to first order we have (T(n+1)-T(n))/dt - T(n)/τ = W(n), which becomes an AR(1) process.
How to force it, ie make W(n) as near iid as possible? A reasonable requirement is that W(n) should have zero autocorrelation for lag 1, and I think, though I haven’t checked, that this is equivalent to what the existing process is trying to achieve. But you could adopt other limits on the A-C of W(n), and vary τ to achieve them.
In terms of going to higher order, as I said earlier, I don’t think that is useful except insofar as it improves the approximation of dT/dt, and I doubt if this is critical.
I have to say that hanging this whole analysis on the assumption that observed temperatures were generated by a white noise flux seems a slender reed, but I think that is what is done.
Schwartz und die Klimasensitivität « Climate Review May 13th, 2008 at 12:29 pm
[...] vom Blog “The Blackboard” hat nicht nur die Arbeit von Schwartz angeschaut, sondern auch die veröffentlichte Kritik. [...]