To figure out if we are talking past each other because we are using the term “Long Term Persistent” differently, I want to ask people questions based on this synthetic temperature series; the axis shows number of years.
- Do you think it shows pseudo-periodic behavior. (ENSO? PDO?)
- Do you think it shows “Long Term Persistence”?
- Do you think it could be described by ARIMA? ARFIMA?
- Do you think it resembles earth’s temperature series?
Feel free to elaborate. For those who want then numbers: here they are.
To make this more fun, I added a second image. Please answer the same questions:
update
For people who want something that is absolutely, positively LTP:
PureFracdiff
100 years of monthly data looks like this:
Don’t know about LTP but “by eye” (!) I see an effect around the 45year mark that I’d want to look into to see if there is anything to it.
Is the series of numbers intended to be a composite representative of multiple values or a series of actual temps from a single device? Or anomalies from a single device? Or other?
Enjoying your current enquiries and the discussions thereof – at my limited understanding it seems as if you are on to a good systematic investigation.
looks persistent flat with some internal cyclic trends
curious
For the time being, this is a series of numbers that I am presenting to people to ask them whether they think this “looks” LTP. It’s 100 years worth of monthly numbers.
Seems to be a wave like component every 20 years or so, could be sun?
BTW looks very long persistent
Ok thanks – please can you point me to your source for your defintion of LTP? If I’ve gleaned it properly you are taking it as an exhibition of power law decay with frequency in noise residuals – ie big hitting power items of very low frequency? Sorry if you have given this and I’ve missed it/misunderstood.
FWIW I’ve googled around on it and the (IMO) tangentially relevant hit I came up with was a paper on ecological population interdependence where LTP appears to be defined as the characteristic where multiple interelated, and potentially competing, drivers work to produce a stable system. Using my still under development 50 year ish fuzzy logic neural processing algorithm I’d say your graphic looks like a stable system that took a hit around the 45 year mark and is back on track. Hmmm – another analogy just came to mind. Hey ho! :0
Curious– My understanding is LTP is power law decay near low frequencies. Other people might supply other definitions. If they exist, I’d like to learn them.
Power law decays at high frequencies is just garden variety random stuff you see in all sorts of ARIMA etc. Those are also stable. So, mere stability can’t be the definition of LTP. (Or if it is, then white noise is LTP. No one calls white noise LTP!)
Thanks again – makes sense 🙂
I vote no on #2. The value of d from armaFit was 0.024 so no LTP to speak of.
Title:
ARIMA Modelling
Call:
armaFit(formula = x ~ arfima(1, 1), data = lucia)
Model:
FRACDIFF(1,1) with method: mle
Coefficient(s):
d ar1 ma1
0.02449 0.85870 0.28556
I haven’t tried looking at models with more coefficients.
Second data set coming soon. It would have been easier if you left off the row names. I’m sure there’s an easy way to get rid of them but I don’t know it.
“1.Do you think …”
I dont think you can underestimate the importance of understanding the origin of the data before going on to try to identify features within it.
I call this “step 1”. Without step 1, further steps are stabbing wildly into the dark.
Duh! Paste into Excel and use text to columns.
Title:
ARIMA Modelling
Call:
armaFit(formula = x ~ arfima(1, 1), data = lucia2)
Model:
FRACDIFF(1,1) with method: mle
Coefficient(s):
d ar1 ma1
0.03134 0.88212 0.38661
Still no (significant) LTP.
Dewitt-
You are correct d=0. However, oddly enough, the data was NOT generated with ARIMA.sim(1,0,1). It was generated with ARIMA(2,0,0).
What package holds armaFIT? Is it R? There seems to be an explosion of packages many of which are mostly wrappers around other packages. I can’t begin to keep track of all of them.
The fft power spectrum of the first data set looks like reddish noise. The second data set, OTOH, has a big peak at the lowest frequency.
Re: lucia (Apr 1 17:03),
The package is fArma. It’s probably just a wrapper for FRACDIFF. As I said, I didn’t look for the best model.
Title:
ARIMA Modelling
Call:
armaFit(formula = x ~ arfima(2, 0), data = lucia)
Model:
FRACDIFF(2,0) with method: mle
Coefficient(s):
d ar1 ar2
0.01637 0.56736 0.23127 0.00000
Are those anywhere close?
Function pacf shows the first two AR coefficients as being significant for both data sets.
DeWitt–Yes. d=0, but very close for the first.
Yes. The 2nd one rises near the lowest frequency.
To create the 2nd one, I used the same ARIMA.sim as the first, and then added a 2nd series created using fradiff with d=0.45, ar=0, ma=0. I scaled so there isn’t much fradiff power there but you can see it boosts the energy at the lowest frequency.
I don’t know what the 2nd once *should* give us when put through fradiff. I can generate a longer string and we can see later. But I did add enough fradiff to make it “LTP”.
What I want people who are telling me they just know something is LTP because of El Nino to see is that both of these have elements that look pseudo-periodic. Even the “not-LTP” AR2 case does mimic ENSO. It looks quite a bit like “weather”.
Its been a long time since I have been in a statistics class but I have fallen asleep in quite a few.
It is going to be very, very hard to prove LTP in the climate. No matter what test is applied and what significant stat is arrived at, it is not going to be convincing to most people. The eyeball method is better for most people.
In regards to the ENSO, most do not believe that one-year ENSO events (which are by far the most common) demonstrate persistence but a series of them in a row and/or multi-year events can leave limited multi-year temperature signals behind.
The AMO, however, has 30 year up cycles and 30 year down cycles in the data available since 1850. In the historical reconstructions, the cycles are less regular and can be even longer and bigger than today so one should not put much stock in it being a regular 30 year up and down cycle.
In my napping, I do remember one particular method that was the most convincing to me in time series analysis and this was just a preliminary step to see if you should carry on with further analysis, and it is the lagged scatterplot.
Given we are talking about the climate which has a seasonal signal and we are concerned with multi-year persistence rather than one month persistence, I plotted one, two and three year lags on your dataset2.
It does not demonstrate persistence. It demonstrates trends that last less that one year and no multi-year components are obvious. So no LTP in dataset2.
http://img638.imageshack.us/img638/1343/testingpersistence.png
On the other hand, HadSST2 for the southern hemisphere which I brought up yesterday demonstrates significant persistence on one, two and three year timeframes. It is not a random scatter (which should have a trend up given the warming) but there is significant density differences at different temperature levels – not at a really solid level, but good enough). So I conclude persistence for the HadSST2 SH at least up to three years.
http://img577.imageshack.us/img577/5138/haddsst2shpersistence.png
Re: curious (Apr 1 15:39),
Curious,
This paper talks a lot about LTP in the climate context. This one is more tangential but more readable.
Re: lucia (Apr 1 18:32),
If the LTP component is buried in non-LTP noise, it’s going to take much longer to detect it. I think I’m going to try some moving averages to get rid of the high frequency noise. If there is LTP, I’m pretty sure that low pass filtering isn’t going to get rid of it.
Bill–
The examples ares synthetic. So I know what they are. What I’m trying to see is what people get.
Construction of set 2 included a component with LTP. It’s small.
But do you call this LTP? What I’m trying to get to is an understanding of you people are using terms. LTP used in Cohn and Lins has a very specific definition.
DeWitt-
I agree that LTP buried in nonLTP noise is going to be hard to detect. It’s interesting to see what people do try.
Bill–
It took me a bit to figure out precisely what you were doing. Those lag plots are a way to visualize the lag-N autocorrelation. Oddly enough, the way LTP is defined, the long lag autocorrelation can be vanishingly small so the lag plots will look like bullseyes. One key feature though is the rate of decay in the lag with time is sufficiently slow that the tails of the autocorrelation still contribute to the sum of the area under the integral (or sum) of hte autocorrelation plot.
So those graphs are good ways to detect positive autocorrelation, but not particularly good ways to detect LTP. Also, I’ve fiddled with a bunch of these and it’s worth knowing that autocorrelation at larger lags is very biased. So,that’s not a very good way to detect the LTP. I’m not sure what is the most effective way, but the lagged plots seem to be pretty ineffective.
Nick — The link to the one you call more tangential doesn’t seem to work.
Re: lucia (Apr 1 21:36),
“Nick — The link to the one you call more tangential doesn’t seem to work.”
It’s a pdf file. Could it be just downloading it for you to some place?
Nick Comment#72593 – thanks, they both look good references.
They’ll take me a (long!) while to absorb but having looked at the Zorita et al paper abstract the comment they made about the difference between individual station records and averaged records caught my eye. I’d welcome any comments you have on this issue – it is an area I remain uncomfortable with.
This also caught my eye:
“The parameters of the autoregressive models were estimated from the observed records in the period up to 1960 to limit the influence of the anthropogenic forcing.”
and again I’d welcome any comments on the difference (if any) using the whole record for parameter estimation would/could make to the analysis.
I acknowledge they reject their null hypothesis with a very high level of confidence and I don’t find this surprising. My simple view is that looking at a range of historical records shows we have periodic interludes of warm and cold spells, so there is historical precedent for clustering of high and low temperatures. Are there any studies you know of which take a similar look at other 100 year ish periods of long station data? Say 1800 – 1900 and 1850 – 1950? I’ve had a read of the Rybksi paper which deals with the longer term (and, yes, it is tougher going!) but it doesn’t seem to provide this equivalent overlapped period analysis? Fully accept that to someone who knows their stats this may be deduceable from the paper (or not relevant?). Also suspect that an “overlapping century analyisis” would be pretty straightforward for one of the whizzes…
Thanks in advance for any comments.
lucia – Offered in the spirit of fun enquiry: Not sure if this would be worth anything but I wonder if you could generate (say) 10 separate pseudo series plots using each of the parameter sets you used for the three plots in the post? I’m thinking putting the resulting sets side by side might be an interesting visual training for the fuzzy logicians amongst us. Then maybe do 3×10 more and put them in random sequence and see if us “fls” can group them? On reflection I guess this sort of thing will have been tackled in psychology?
Nick–
The link sends me to a page that says the document I see does not exist. But weird… it worked this time. (I’d clicked from the email I get for each comment. Weird…)
Hi Lucia
Some help needed. Please check out the comments on this Bishop Hill thread and this page.
.
The issue is the funkiness of a temperature reconstruction graph published in Nature Climate Change (link). More specifically, look at the y-axis.
.
How can a graph such as this, come to have ‘0.6’ in the place of ‘0.0’? I want to know for my edification. The code line simply seems to generate ticks at every 0.2 interval and nothing more. I am assuming this line of Python code below is the one that is relevant. Re:
.
Nick–
Believe it or not…. The third figure I added in the update is at exactly the fracdiff Zorita used with d=0.45. That’s what trendless “weather” would look like under the assumptions of their tests with fractional difference noise. I plotted 100 years. (I computed 500 in the attachment.) I’ll repeat the graph:
As you can see, there are “runs” that might look like a trend for a while. But they aren’t horribly long.
To make the open circles in their figure 1, they must have just picked different values of d and then generated one huge number of these and then computed the number of records in the final 100 years and then created a probability distribution function, for numbers of records.
Of course, their results could change if they used even more complicated noise structures. Cohn and Lins did use more complicated FARIMA. (Is it FARIMA? Or ARFIMA? I seem to see both.)
Why don’t you ask whoever made the graph?
You won’t figure it out?
.
It is a question, it is a puzzle.
.
It has appeared in an article extolling the virtues of open access data and software code, in climate science.
.
Thanks for your help.
Shub–
No. Because I don’t really care how the 6 appeared where the 0 probably should have appeared. As a general rule, I’m not all that into discovering the underlying reason for pesky errors like how typos arise on figures.
Out of curiosity, if you think the puzzle is worth looking into, why don’t you set aside this glorious early spring afternoon in March to dive into the code and find out how a “6” might appear where a “0” ought to appear? I’m spending the afternoon getting ready to grill chicken, breathing spring air etc.
Lucia,
In addition to all those things you’re doing today, you ought to add “turning the page on my wall calendar”. 🙂
Thanks again.
.
One never knows what one learns!
Did these with FFTs, and I’m a little rusty, so someone might want to check with a no-kidding autocorr routine.
Try 2 on the image:
Re: lucia (Apr 2 08:30),
They would be if you defined your time step as a year or a decade rather than a month.
How do I embed an image? I tried standard anchor and image tags, but wordpress seems to have eaten some of it.
Anyway, here’s the link to a picture of the autocorrelations:
https://picasaweb.google.com/lh/photo/OMyEq5wT8mu0J6GMuJAcnhz9Qy9fw-1iTio99rKT-Zs?feat=directlink
Why so many data points in the third set? Is it repeated in the text file?
Re: jstults (Apr 2 17:03),
You can’t. You have to have editor privileges.
lucia,
fracdiff says that d=0.45 on the pure fracdiff set, ar and ma coefficients weren’t significant. arfima, which is a front end to fracdiff that is supposed to fit all the coefficients also failed miserably as it produced large ar and ma coefficients. I tried a whole raft of Hurst coefficient estimators. They mostly fail miserably on the first two data sets but do ok on the third set.
Dewitt– The pure fracdiff is d=0.45 with ar=0 and ma=0.
The middle one is a weighted sum of the other two. It’s (first + 0.05* second). So, it’s sort of a bastard. The spectrum looks like the sum of the first. I’m not sure what it “should” come out to. I was thinking of creating more later and seeing what I get if I had thousands of years of data vs. short numbers of years.
I know it has to be LTP because one element is LTP. So that part will continue to be exhibit correlation long after the other bit has faded away.
Did you only use fracdiff on the middle one? Or was there some other tool?
Sure. And since this is all self-similar scaling, I guess the hurst coefficients are identical for monthly averaged and annual averaged.
jstultz–
Wordpress does protect the blog by preventing visitors from posting images. It’s mostly a good thing because… well… who knows what someone might post at some random blog. The reason images appear here is when I notice an image link, I edit and add the image tags.
I tried to edit the comment to make the image appear, but those are link to a page, not the image itself.
Lucia,
I enjoyed reading your posts on LTP and the many thoughtful responses. Thank you!
I do have a couple of comments: First, given only a single time series, it is not easy to tell whether the underlying stochastic process possesses LTP. That point is made in “Naturally Trendy” (“NT”). Second, the motivation in NT for assuming LTP in the temperature time series was not based on statistical analysis of the time series at hand. Lins and I discussed this in NT:
Of course there are physical systems where LTP arises as a consequence of governing differential equations. Koutsoyiannis has explored such situations and has written a number of thoughtful papers on this topic. But that work does not necessarily address the question at hand.
Finally, as I think is clear from the discussions here and elsewhere, the concept of LTP is simple but its implications can be subtle.
I agree. Particularly if the series is short.
I’m not under the misimpression that your reason for diagnosing this is analyzing the NH thermometer record. I’m ok with people considering the possibility that the noise is LTP and think there are good reasons to suspect it’s existence.
Rather my concern is that not-withstanding your reasons for assuming LTP, when feeding your NH into fracdiff that module cannot distinguish between “LTP” and apparent LTP arising from non-linearity in the long term trend. So, whatever “d”‘s you get (or for that matter, the total power (σ), as well as p,q in (pdq) are affected by the deterministic trend.
These computed parameters are then used to test whether the trend is statistically significant.
I think if you might better understand my concern as follows:
Consider the analysis you did for figures 2 and 3. Now repeat the types of analysis summarized in figure 1, but when creating synthetic data, include a quadratic in the deterministic function. That is, for type I tests do this:
generate data with
y=bt^2+u with u = some noise and I suggest -T/2 < t < T/2 so you should get no net trend change in y over that span.
Then fit to
y=mt+u and test whether ‘m=0’ as your null.
And test whether you accept or reject ‘m’ with various types of noise. (Just as you did in figure 2.) See how various magnitudes of ‘b’ affect your results. You should find the existence of the quadratic that is neglected when computing the properties of the noise u will result in a smaller fraction of false rejections. The degree to which the quadratic will reduce your false rejections will depend on the magnitude of ‘b’ and off hand, I think the dimensionless parameter bT^2/σ where σ is a standard deviation of ‘u’ should be useful for evaluating the impact of the neglected non-linearity.
Now for the power tests do this:
y=mt+bt^2 +u
where u is some sort of noise. Here, in the span -T/2 <t <T/2, m is in a sense the deterministic trend — in the sense that it is repeatable and in the sense that you want to detect whether or not y is increasing with time.
Once again, see how ‘b’ affects this.
What you will find is that even in the case of OLS, the non-linearity in the deterministic trends will result in lower type I error than you intend, and will reduce power to find any deterministic trend. The problem will get worse for AR1. I suspect, but I am not sure, the problem will get even worse for LTP type noise.
After that consider the fact that during the period you analyze, the conventional view is that the deterministic trend is non-linear.
Given that, how is one to interpret your finding if they tend to believe the deterministic view about the evolution of the deterministic trend. Other people doing different analyses find the NH temperature rise is inconsistent with the properties of the random noise– and they find this even if they — like you– consider the possibility the noise is LTP.
So, my criticism isn’t “You should consider the possibility the noise is LTP”. It is that “Even if you assume the noise is LTP, you need to consider the effect of non-linearity in the deterministic trend when testing a hypothesis that m=0”.
Note that if you had rejected m=0 in this way, the non-linearity issue I discuss wouldn’t matter. We would know that the rejection would persist when you consider the non-linearity. But since your novel result is “fail to reject m=0”, while other reject, how is a reader to know whether the reason for your results is merely that you fed “fracdiff” a time series that contains a deterministic trend?
Of course, I’m bringing this up as a qualitative issue. It’s more work to deal with it quantitatively– but it’s an issue and I think it’s one that you might consider worth addressing.
Lucia,
I think I agree with what you wrote, but some clarification might help anyway.
The question NT addressed was whether the evidence for a deterministic trend in the temperature record was particularly compelling. NT concluded that the presence of LTP in climate processes (almost no one dismisses this) can explain the observed temperature trend (almost no one dismisses this, either). As a result, there is no need to invoke a deterministic trend to explain the observed record. Restated, the observed record could “easily” arise as a realization of a stationary LTP stochastic process.
I don’t think we disagree about any of this.
However, I think I also agree with your position (if I understand it correctly). Let me try to restate it so that it conforms to my way of thinking: If one reverses the null hypothesis (that is, assumes a deterministic quadratic trend with white or slightly off-white noise and then tests for evidence of LTP), the statistical evidence is going to be weak.
In response to your question — “how is one to interpret your finding if [she tends] to believe the deterministic view about the evolution of the deterministic trend”? — one _could_ reach an _apparently_ different conclusion: The observed trend is “not inconsistent” with a quadratic deterministic trend plus low-order ARMA noise.
However, this is a distinction without a difference. Either way you set up the null hypothesis you get an insignificant p-value. Which means, of course, the data do not reveal the validity of the hypothesis.
So what is _the_ truth? I don’t know the answer. But I’m pretty sure about one thing: Those who claim to “know” the answer are not fully comprehending the complexity of the situation.
Tim–
My concern is doing a test in a way that has lower power that necessary under the prevailing hypothesis. (The prevailing hypothesis: the ‘true’ trend is non-linear and looks sort of like the model mean trends from climate models).
I’m going to be thinking about this a bit more tomorrow when I can figure out how to meld what I’m thinking into your code. I’m going to email a rather stupid question about σ in figure 3 privately.
Partial autocorrelation function suggests AR(2) for the first series. In that sense it looks ‘most synthetic’, for the others one gets AR(4) or more. For the real HadCrut, PACF seems to behave more badly than any of these.
Things would be easier if we knew the non-linear deterministic trend (and would be able to predict it given the future forcings). I agree with Tim Cohn, no need to invoke a deterministic trend – unless the predictions based on that are more successful.
UC–
But the physics suggests that the forced component of temperature (i.e. deterministic portion) should evolve over time and increase. So the difficulty is that to be more convincing it would be nice to know what the residuals are relative to the forced trend, and then show something.
For sure, we should find those residuals. But can we figure out the forced trend without fitting it to the observed temperature series? If we have to do a fit, we have a problem: in the case of long memory errors the polynomial trend is estimated with a lousy precision when compared to iid errors. One can fit more and more exogenous factors, until the residuals look white and then conclude that the fit is good (because he then can assume iid errors). The other says the fit is bad because the errors were of long memory type in the first place. If someone derives the forced trend solely based on physics we don’t have this kind of problem.