by Richard Saumarez
A recent debate at Climate Etc. has involved autocorrelation and trends. Essentially, the argument is that the temperature signal is perturbed by a random signal, and the influence of this random fluctuation persists for an unknown length of time and so the temperature at any particular instant is influenced by its past history.
Why does persistence produce trends?
It may not be intuitively obvious why trends should occur in persistence models. Figure 1 shows a resting point (stage 0) at zero. A random, normally distributed value is added to it and it moves to point 1.
Figure 1 An Illustration of how persistence leads to trends. Upper trace shows how persistence stops alters the probability of the temperature returning to baseline, creating a trend. The bottom trace shows how a lack of persistence does not bias the probability of the next step.
The value at point decays slowly, it has persistence, so that after a time interval it reaches point 1a. A second random value is then added to it. Notice that probability of it returning to zero or below is the area of the distribution curve that is below zero, i.e.: very small. This illustrates why trends occur with persistence.
As the value is increased, note that the rate of decay increases – compare the decay rate between points1 and 1a with that between 5 and 5a. Therefore the observed value will remain above zero until there has been a run of negative disturbances that will return it.
The lower trace shows a non-persistence model. In this case, the value has returned to zero before the next disturbance and the observed values will have a mean of zero.
This is the basis of the argument that “persistence” will cause trends in temperature. If an impulse of radiation is stored as heat, it carries over so that at the next impulse, the temperature will be slightly raised and become raised still further, resulting in a temperature trend. The controversy is how significant is this effect and can observed temperature trends be explained on the basis of a simple heat reservoir and a random change in radiative fluxes?
What models of Heat storage or Persistence can we use?
There are two basic models, one is a power law model in which the temperature decays by a power of time and is often called the Hurst law, who defined in terms of Hydrology.
T = (Ht)g
(JC formatting note: the g superscript should be gamma).
(t is time, H is a scaling parameter and g, the exponent, describes how the process varies with time – g is negative).
Using this idea, Figure 2 shows a model, consisting of a leaky reservoir, which is surrounded by a dam that is more porous at the top than at the bottom. If it is filled suddenly with water at a high temperature, water will escape rapidly from the top and as the level subsides, there will be a slower leak through the bottom levels. Thus the level will decay at a slower rate as the reservoir empties, not only because the pressure head is smaller, but also because the rate of transfer is reduced by passing through a smaller area of less porous material.
Figure 2 A model of a power law system with the response shown (black). The responses for g=0.6-0.8 are shown in blue and the responses to the reservoir being filled to different levels are shown in red.
I have solved this model to give the decay curve and corresponds to a gamma of 0.7. It is important that the response to an impulse of water entering the reservoir will decay differently depending on how full the reservoir is (red curves) and therefore this system is non-linear. When computing a response, the system can be treated as quasi-linear by assuming that the actual level of the water does not vary much from its mean level. This allows linear methods of processing to be applied. If however, there is a wide range of levels, because the range of inputs is large compared to volume of the system, the behaviour of the system must be solved explicitly through its governing differential equations. Whether one can really treat the temperature system as quasi-linear, or not, depends on measuring the inputs (if we knew what they were) and the responses of a series of temperature records and determining the coherence between them.
The alternative model is a multi-compartment model that has exponentially decaying compartments as shown in figure 3. It is well known that temperature persistence does not conform to a simple, single exponential decay, but a system with a fast component and a large volume slow component can produce much longer levels of persistence. As shown in figure 3, this model produces a decay that is very similar to a power law decay and, if communication between the compartments can occur, the decay curve becomes almost indistinguishable from a power law model.
Figure 3 The response of a multi-compartment model. The power law responses for g=0.7 and 0.8 are shown in blue. Note that a 3 compartment model has a very similar response to a power law model.
There is a very important assumption in the use of persistence models to calculate global temperature trends.
This assumption allows one to assume an “average model” for persistence and to drive it in order to produce a simulated temperature record.
Determination of the model and its relationship to trends
This was done in a number of steps:
1) Temperature records were extracted from the HADCRUT III, UTAH and GHCNM databases and all those at least 100 years long with less than 5% of the data missing were selected. The missing data was interpolated linearly across the gaps.
2) The yearly (“seasonal”) fluctuations were removed.
3) This data was used to generate an autocorrelation and signal power model.
4) A system response was created using either a bi-exponential or a Hurst model.
5) Gaussian noise was used to create a sequence 1364 years long that have the same variance as the residual signals. A one hundred year record was randomly selected from each sequence, in order to account for measured temperature trends to start at a random time with respect to any trend in temperature. The amplitudes and duration of the trends were extracted from 100000 sequences.
The ACF is a measure of how a sample in the record is correlated to its neighbouring samples and gives a measure of the persistence as it demonstrates the effects of past samples in determining the current sample. In principle, a copy of the signal is made and is shifted in time by one sample. The correlation coefficient between the shifted copy and the original version determines the relationship between the signal and that at one sample apart. This is repeated for all possible shifts of the copy of signal and the correlation coefficient is plotted against the value of the shift (Figure 5). In practice this can be done efficiently through the discrete Fourier Transform of the signal.
There are a number of problems in calculating the ACF. The first is that as the signal copy is shifted and each sample multiplied with the sample at the same position in the original signal, the samples at the end of the shifted signal “drop off” the end and so do not contribute to the correlation at that shift, making the estimate less reliable. Therefore only shifts of about 1/5th of the total signal length should be used, in our case about 20 years.
The method used by L and LLE is slightly different from the autocorrelation function, in that they calculate a cumulative sum of the temperature excursion and manipulate this sum as a function of time by fitting a second order polynomial as a smoothing and de-trending step. They then derive a persistence factor a =1-g/2. Hence an alpha, which they quote, of, say, 0.6 corresponds to a g of 0.8.
Extracting the ACF from temperature records is quite difficult and may lead to artefacts.
The underlying temperature signal is dominated by the yearly variation as can be seen from the amplitude spectrum of the signal shown in figure 6, and this has to be removed to estimate the long-term persistence. The method by which this done can lead to quite different results for the ACF. One approach might be to simply low-pass filter the signal with a cut-off point well above 1 year-1 as shown by the purple filter response in figure 6. This has a huge effect on the ACF as shown by the purple broken line in figure 5. A selective band-pass filter can be used but this appears to produce “ringing” in the short term ACF and to alter its persistence in the tail. The method, used here, which causes least disturbance to the ACF is to interpolate across the spectral bands in question and the effect of this filter on random noise is shown as the green line in figure 5.
Figure 5 The Mean ACF of the temperature records shown as blue squares with a 2 compartment model (black) and a power law model (red) fitted to the data. The dashed purple line shows the effect of using a low pass filter to eliminate seasonal trends (i.e.: disastrous!) and the green line shows the effect of the seasonal de-trending method having very little effect on the ACF
Once the data has the yearly variation removed, it is normally distributed.
The Yearly ACF values are plotted as blue squares in figure 5 and the best-fit power law (red g=0.875) and bi-exponential (black) models are shown. In this data set, he bi-exponential model appears to fit slightly better, however, when the log ACF is plotted against the log delay, it is clear that there are considerable errors involved and it is not possible to say how closely the data conforms to a power law model.
Figure 6 Mean amplitude spectrum of the temperature records with peaks at 1,2,3 year-1. The red line shows a selective band pass filter and the purple line a low pass filter designed to remove seasonal trends. Interpolation across these peaks gives least artefact in the ACF.
Figure 7 Mean ACF plotted on a log-log plot. The vertical blue lines show the 95% limits of the distribution and the red dashed lines show a bi-compartment model.
It is not possible with this data set to be certain of the model to describe persistence, but one can use a power law or multi-compartment model with different assumed values and determine their abilities to produce trends.
Figure 8 Simulated record (blue) with low pass filtered version (black) and trends superimposed in red.
Figure 8 shows the residual temperature (light blue) using a power-law model with g=0.875. [This has been calculated using the ODE, not through the assumption of linearity and convolution]. The magnitude of disturbance has been drawn from a distribution so that the power in the model signal is drawn from the same distribution as in the observed seasonally detrended signals. This is repeated 100000 times to give a distribution of the expected residual signals. This has been repeated for gamma of 0.7.
The global temperature signal is considered to be the sum of a number of temperature signals drawn from stations that are perturbed with a highly correlated disturbance. Since we have only one long-term temperature record, we cannot say anything about the statistics of this signal. But, we can determine the ACF of the simulated signals and these show considerable variability as shown in figure 9. The mean ACF and its 95% limits are shown for a power law model of g=0.7 (red) and for a 2 compartment model which approximates g=0.9. Irrespective of how persistence is extracted from these signals, there is a considerable range of g fitted to the signal and the power law and 2 compartment model are indistinguishable given a single temperature record. The 95% limits for a true g =0.7 is 0.58-0.87.
Figure 9 Mean ACFs and their 95% error limits for a power law model (Red-g=0.7) and a 2 compartment model (Black)with an approximate g of 0.9.
Estimation of trends in the simulated signal.
The black line in Figure 8 shows a low-pass filtered version of the signal with trend lines plotted in red. There are clearly shorter-term trends within the signal and, if one decomposes a signal into very short-term trends, it is not possible to distinguish a long-term trend. The approach used here, is to low pass filter the signal heavily, so as to eliminate the rapid fluctuations and use the maxima and minima of this signal to determine the ends of trends. The trend line is then fitted to the data between these intervals.
This is shown more clearly in figure 10 where the low pass filtered signal is shown, and the trends extracted, using a filter with cut off points at 1/30 years (solid red line) and at /100 years (dotted red line). In general, the magnitude of the trend is less as longer trends are identified and there will be more variability around the trend.
Figure 11 Histograms of trend lengths and excursions. Power law processes with , A, g=0.9 (upper) and, B, g=0.7 (lower). The trends > 100 years are shaded red and the region between 60 and 100 years is also shaded
Using this approach, the distribution of trends for power law gammas of 0.9 and 0.7 (i.e.: a=0.55 and 0.65) can be calculated and are shown in figure 11. According to these calculations the probability of a 100 year trend in either direction with g=0.9 is 0.02 with a mean value of 0.27o C. When g =0.7, the probability of a 100 year trend is 0.09 with a mean value of 0.4o , while the probability of a >100 year trend of >0.7oC is 0.02.
The purpose of the analysis presented by Ludecke, Link and Ewert is to characterise the Earth’s temperature as a system, in particular the lag in temperature with energy changes.
Any system that involves lag can produce trends. Whether the lag is imposed via a Hurst process, a multi-compartment system or anything else is a secondary consideration. However, the problem is to characterise the process sufficiently accurately to make predictions about the magnitudes of trends is another matter. The data available is strictly limited and one cannot easily distinguish between different physical models that may underlie the trends. Satellite temperature and radiation flux data will eventually lead to much better temperature series but this data would be better incorporated into physical models that incorporate the mechanisms of energy storage and would render the approach taken here unnecessary.
One of the major problems in this analysis is that temperature signals at different points on the Earth are not independent. If heat is stored in the ocean and transferred by currents to another point, or by some other mechanism, the temperature rise at this point is not an “external trend” as stated, which has to be eliminated, but is an important part of the mechanism of the creation of regional trends within the system as a whole.
Figure 12 A model of energy transfer between different regions with peristence
The addition of energy into a region, as shown in figure 12, will give a false value for the persistence. However as system similar to that in figure 12 is analysable, in principle, provided the system is linear. If one is using a power law representation of the persistence, the problem is far more difficult. Essentially, one writes the energy flux equation for each component in the system as a set of coupled differential equations, or in the complex frequency domain, and solves then to obtain the parameters of the persistence in each element. This can be done via the ACF of the temperature at each point. In practice, this is practically impossible with the data available and the solution is highly unstable and very sensitive to small errors in data, even for a restricted subset of the problem. However, LL&E have attempted to address this problem by de-trending the data, which is an arbitrary approach and smoothes the data, which may not necessarily give a better value for the calculated persistence.
A major problem is the quality of the data on which a model is constructed. As I have pointed out, using correlation methods to determine persistence, long-term records are required and the number of “clean” long term records, especially in the tropics, are strictly limited and, to mind, global conclusions based on this data should treated sceptically.
 See Bendat JS & Piersol AG. “Random Data: Measurement and Analysis Techniques”. Chapter 7 – Multiple Input – Output Models, Conditional Spectral Density Functions.
Bio note: As Professor Curry asked me to give some biographical detail, I should explain that after medical school, I did a PhD in biomedical engineering, which before BME became an academic heavy industry, was in an electrical engineering department.
Richard Saumarez’ previous posts at Climate Etc.:
- Climate, control theory, feedback: does it make sense?
- Does the aliasing beast feed the uncertainty monster?
JC note: I received this post via email from RS. I posted this at Climate Etc. based on the interest generated by his previous posts and also by the Ludecke papers, and the fact that this is a topic I want to learn more about. This is a guest post, which I am very appreciative of, but the content reflects the perspective of RS and not myself.