Image Credit: Climate Data Blog
By Richard Linsley Hood – Edited by Just The Facts
The goal of this crowdsourcing thread is to present a 12 month/365 day Cascaded Triple Running Mean (CTRM) filter, inform readers of its basis and value, and gather your input on how I can improve and develop it. A 12 month/365 day CTRM filter completely removes the annual ‘cycle’, as the CTRM is a near Gaussian low pass filter. In fact it is slightly better than Gaussian in that it completely removes the 12 month ‘cycle’ whereas true Gaussian leaves a small residual of that still in the data. This new tool is an attempt to produce a more accurate treatment of climate data and see what new perspectives, if any, it uncovers. This tool builds on the good work by Greg Goodman, with Vaughan Pratt’s valuable input, on this thread on Climate Etc.
Before we get too far into this, let me explain some of the terminology that will be used in this article:
—————-
Filter:
“In signal processing, a filter is a device or process that removes from a signal some unwanted component or feature. Filtering is a class of signal processing, the defining feature of filters being the complete or partial suppression of some aspect of the signal. Most often, this means removing some frequencies and not others in order to suppress interfering signals and reduce background noise.” Wikipedia.
Gaussian Filter:
A Gaussian Filter is probably the ideal filter in time domain terms. That is, if you consider the graphs you are looking at are like the ones displayed on an oscilloscope, then a Gaussian filter is the one that adds the least amount of distortions to the signal.
Full Kernel Filter:
Indicates that the output of the filter will not change when new data is added (except to extend the existing plot). It does not extend up to the ends of the data available, because the output is in the centre of the input range. This is its biggest limitation.
Low Pass Filter:
A low pass filter is one which removes the high frequency components in a signal. One of its most common usages is in anti-aliasing filters for conditioning signals prior to analog-to-digital conversion. Daily, Monthly and Annual averages are low pass filters also.
Cascaded:
A cascade is where you feed the output of the first stage into the input of the next stage and so on. In a spreadsheet implementation of a CTRM you can produce a single average column in the normal way and then use that column as an input to create the next output column and so on. The value of the inter-stage multiplier/divider is very important. It should be set to 1.2067. This is the precise value that makes the CTRM into a near Gaussian filter. It gives values of 12, 10 and 8 months for the three stages in an Annual filter for example.
Triple Running Mean:
The simplest method to remove high frequencies or smooth data is to use moving averages, also referred to as running means. A running mean filter is the standard ‘average’ that is most commonly used in Climate work. On its own it is a very bad form of filter and produces a lot of arithmetic artefacts. Adding three of those ‘back to back’ in a cascade, however, allows for a much higher quality filter that is also very easy to implement. It just needs two more stages than are normally used.
—————
With all of this in mind, a CTRM filter, used either at 365 days (if we have that resolution of data available) or 12 months in length with the most common data sets, will completely remove the Annual cycle while retaining the underlying monthly sampling frequency in the output. In fact it is even better than that, as it does not matter if the data used has been normalised already or not. A CTRM filter will produce the same output on either raw or normalised data, with only a small offset in order to address whatever the ‘Normal’ period chosen by the data provider. There are no added distortions of any sort from the filter.
Let’s take a look at at what this generates in practice.The following are UAH Anomalies from 1979 to Present with an Annual CTRM applied:
Fig 1: UAH data with an Annual CTRM filter
Note that I have just plotted the data points. The CTRM filter has removed the ‘visual noise’ that a month to month variability causes. This is very similar to the 12 or 13 month single running mean that is often used, however it is more accurate as the mathematical errors produced by those simple running means are removed. Additionally, the higher frequencies are completely removed while all the lower frequencies are left completely intact.
The following are HadCRUT4 Anomalies from 1850 to Present with an Annual CTRM applied:
Fig 2: HadCRUT4 data with an Annual CTRM filter
Note again that all the higher frequencies have been removed and the lower frequencies are all displayed without distortions or noise.
There is a small issue with these CTRM filters in that CTRMs are ‘full kernel’ filters as mentioned above, meaning their outputs will not change when new data is added (except to extend the existing plot). However, because the output is in the middle of the input data, they do not extend up to the ends of the data available as can be seen above. In order to overcome this issue, some additional work will be required.
The basic principles of filters work over all timescales, thus we do not need to constrain ourselves to an Annual filter. We are, after all, trying to determine how this complex load that is the Earth reacts to the constantly varying surface input and surface reflection/absorption with very long timescale storage and release systems including phase change, mass transport and the like. If this were some giant mechanical structure slowly vibrating away we would run low pass filters with much longer time constants to see what was down in the sub-harmonics. So let’s do just that for Climate.
When I applied a standard time/energy low pass filter sweep against the data I noticed that there is a sweet spot around 12-20 years where the output changes very little. This looks like it may well be a good stop/pass band binary chop point. So I choose 15 years as the roll off point to see what happens. Remember this is a standard low pass/band-pass filter, similar to the one that splits telephone from broadband to connect to the Internet. Using this approach, all frequencies of any period above 15 years are fully preserved in the output and all frequencies below that point are completely removed.
The following are HadCRUT4 Anomalies from 1850 to Present with a 15 CTRM and a 75 year single mean applied:
Fig 3: HadCRUT4 with additional greater than 15 year low pass. Greater than 75 year low pass filter included to remove the red trace discovered by the first pass.
Now, when reviewing the plot above some have claimed that this is a curve fitting or a ‘cycle mania’ exercise. However, the data hasn’t been fit to anything, I just applied a filter. Then out pops some wriggle in that plot which the data draws all on its own at around ~60 years. It’s the data what done it – not me! If you see any ‘cycle’ in graph, then that’s your perception. What you can’t do is say the wriggle is not there. That’s what the DATA says is there.
Note that the extra ‘greater than 75 years’ single running mean is included to remove the discovered ~60 year line, as one would normally do to get whatever residual is left. Only a single stage running mean can be used as the data available is too short for a full triple cascaded set. The UAH and RSS data series are too short to run a full greater than 15 year triple cascade pass on them, but it is possible to do a greater than 7.5 year which I’ll leave for a future exercise.
And that Full Kernel problem? We can add a Savitzky-Golay filter to the set, which is the Engineering equivalent of LOWESS in Statistics, so should not meet too much resistance from statisticians (want to bet?).
Fig 4: HadCRUT4 with additional S-G projections to observe near term future trends
We can verify that the parameters chosen are correct because the line closely follows the full kernel filter if that is used as a training/verification guide. The latest part of the line should not be considered an absolute guide to the future. Like LOWESS, S-G will ‘whip’ around on new data like a caterpillar searching for a new leaf. However, it tends to follow a similar trajectory, at least until it runs into a tree. While this only a basic predictive tool, which estimates that the future will be like the recent past, the tool estimates that we are over a local peak and headed downwards…
And there we have it. A simple data treatment for the various temperature data sets, a high quality filter that removes the noise and helps us to see the bigger picture. Something to test the various claims made as to how the climate system works. Want to compare it against CO2. Go for it. Want to check SO2. Again fine. Volcanoes? Be my guest. Here is a spreadsheet containing UAH and a Annual CTRM and R code for a simple RSS graph. Please just don’t complain if the results from the data don’t meet your expectations. This is just data and summaries of the data. Occam’s Razor for a temperature series. Very simple, but it should be very revealing.
Now the question is how I can improve it. Do you see any flaws in the methodology or tool I’ve developed? Do you know how I can make it more accurate, more effective or more accessible? What other data sets do you think might be good candidates for a CTRM filter? Are there any particular combinations of data sets that you would like to see? You may have noted the 15 year CTRM combining UAH, RSS, HadCRUT and GISS at the head of this article. I have been developing various options at my new Climate Data Blog and based upon your input on this thread, I am planning a follow up article that will delve into some combinations of data sets, some of their similarities and some of their differences.
About the Author: Richard Linsley Hood holds an MSc in System Design and has been working as a ‘Practicing Logician’ (aka Computer Geek) to look at signals, images and the modelling of things in general inside computers for over 40 years now. This is his first venture into Climate Science and temperature analysis.
Willis Eschenbach says:
March 17, 2014 at 4:45 pm
“So you are telling us that you can flawlessly reconstruct the original input trace, given only the result of a gaussian smoothing … but only if you have access to the original input trace??”
I did not say that if you re-read my post. What I said was that the signal is split accurately into two parts, only one of which is plotted above. Please look carefully before you comment.
“In practice we would probably truncate with a Hamming window rather than with a rectangle”
The idea is to use the central lobe of the second sinc in the way you would use a Hamming or other windowing function, to ensure it fades to zero at the end of the kernel. This is times to be at zero crossing of the of the first sinc so as to have the smoothest fade out possible. In practice the means you usually have 3,5 or 7 etc lobes. Outside that it’s zero so the windowing sinc never goes negative. That is , in effect, a further rectangular truncation.
If you think replacing the rectangular truncation by a second bell shaped function would be better I’d like to see it, though I suspect it would do more damage to the nice trapezoidal form than it would help by reducing the ripple. For a start it would mess up the nice flat section of the pass band.
I’m always ready to be proved wrong and if you can improve on Lanczos I’d be happy to be wrong and probably adopt your filter.
BTW I usually use lanczos if I need to remove annual before doing FFT, since it does less damage to the bit I want to keep. I can remove 12mo and lower yet still have 24mo and above pretty much intact.
cd says:
March 17, 2014 at 4:40 pm
“You cannot apply a Gaussian filter to an image, save the image then pass it to someone without the original image, then given just the kernel design reconfigure reproduce the original image. ”
No-one – certainly not me anyway – is claiming any such thing.
What you can say is that the low pass part is an accurate part of the signal below 15 years in period. The high part of the signal can be obtained by subtracting that output from the original input (possible 1 month phase delays being correctly accounted for).
Do you know how a band pass splitter circuit work? In most cases they are cascaded low pass filters with each stage subtracted from the previous stage to reduce component count and so as not have to construct full band pass, two sided, filters for each stage.
This is the digital equivalent of such a first stage splitter.
Willis Eschenbach says:
March 17, 2014 at 4:50 pm
“The idea that we can squint at a smoothed temperature record and say what has caused the swings, however, is a bridge far too far.”
Well you have a fairly good centre frequency with a fairly tight bandwidth around it to look for potential clients in 🙂
There has to be something doing it and I don’t buy convenient volcanos, SO2 and CO2 as being the cause.
cd “Any filter will result in loss of information ”
No, nothing is “lost”. If you have box fruit and put all the apples in one bowl and all the oranges in another, you have not “lost” the apples. You still have them. What you have is two bowls or fruit. This may help you determine the weight of apples that you have bought,. for example. You would not be able to do this with all the fruit in one box.
Once you have analysed the weight of the apples, the oranges and assorted ‘noise’ of fruit remaining in the box you can either put them all back to together, sell the oranges or eat the apples.
Nothing is “lost”
cd says:
March 17, 2014 at 4:51 pm
“I think you are playing fast with the rules here and being quite “obfuscatious” with your interpretation of Willis’ comments.
Any filter will result in loss of information”
No I was being accurate. There is no loss, the ‘missing’ data is in the other ‘stop’ band as it must be.
The two parts added together MUST be the input. There is a possibility that some frequencies around the corner value will ‘leak’ into both but the sum of the parts….
cd says:
March 17, 2014 at 4:51 pm
CD – exactly so, and Richard said the same thing to me at least twice above. Your guess!
It would be like filing your income tax by submitting half the data and then argue that nothing was actually missing because you yourself still have a full copy! That should work!
“There has to be something doing it and I don’t buy convenient volcanos, SO2 and CO2 as being the cause.”
I agree about what it’s not but you need to be careful about seeing it as a definite cycle which “must” have a cause, when there are barely two full cycles and a lot of what is there is “adjusted” to buggery.
3RM is close to gaussian and it is a trap of gaussian filters to make everything look like nice smooth sine waves. Be careful not to jump to too firm conclusions.
Richard
I am disappointed to say the least, I find myself on Willis’ side whose dealing with you has been far from fair – or at least I thought.
What you can say is that the low pass part is an accurate part of the signal below 15 years in period.
Right, that is a mantra, you’ve been repeating that throughout even when it is irrelevant to the question being asked. I agree with your work but not to you answer to Willis’ rather proper question.
Do you know how a band pass splitter circuit work?
Sorry you didn’t understand me this obfuscation, it is what people do when they tie themselves in knots. Do you understand sample density prevalence in discrete series and on and on we go. I think Willis’ statement was clear your answer far from it.
In most cases they are cascaded low pass filters…
You see a cascade filter IS akin to a Basis led method such as B-Spline. Yet when asked you said it wasn’t, I really doubt that! It certainly sound like one.
cd says:
March 17, 2014 at 5:12 pm
“I am disappointed to say the least, I find myself on Willis’ side whose dealing with you has been far from fair – or at least I thought.”
I am sorry if I am not explaining myself better then. My fault.
“Right, that is a mantra, you’ve been repeating that throughout even when it is irrelevant to the question being asked. I agree with your work but not to you answer to Willis’ rather proper question.”
Which is what? That filters MSUT lose information? Pray tell, where does that loss go? It must go somewhere. Where do you think it is and what is its magnitude?
“Sorry you didn’t understand me this obfuscation, it is what people do when they tie themselves in knots. Do you understand sample density prevalence in discrete series and on and on we go. I think Willis’ statement was clear your answer far from it. ”
I am sorry, I’ll try again.
“You see a cascade filter IS akin to a Basis led method such as B-Spline. Yet when asked you said it wasn’t, I really doubt that! It certainly sound like one.”
A cascaded filter set as described allows for appropriate correction of the previous stage until you get to the point where any errors that are left are below the input digitisation ones.
Greg says:
March 17, 2014 at 5:11 pm
“3RM is close to gaussian and it is a trap of gaussian filters to make everything look like nice smooth sine waves. Be careful not to jump to too firm conclusions.”
The only conclusion I have jumped to is that there is a strong ~60 wriggle to the data that needs a good, physics based, explanation. Chance and convenience are the last choices items on my list.
Greg
Once you have analysed the weight of the apples, the oranges and assorted ‘noise’ of fruit remaining in the box you can either put them all back to together, sell the oranges or eat the apples.
One should judge the weight of someones understanding by the appropriateness and clarity of their analogy. This however, is rather pathetic. And if you wonder why I am so terse and to the point, it is because you’re just as dismissive with others. If you demand such high standards with others, can I suggest you either lower your standards or stop pretending coming from such high ground to begin with.
Nothing is “lost”
So give us a workflow using a freeware image processing package (e.g. gimp). If you can do it then you’ve just got the Noble prize for physics cause I know a whole populous of astronomers (including astrophysicist) who need your brilliance right now.
Sorry Nobel…
You seem the type of chap to make a big deal of this Greg.
Bernie Hutchins says:
March 17, 2014 at 5:10 pm
“It would be like filing your income tax by submitting half the data and then argue that nothing was actually missing because you yourself still have a full copy! That should work!”
It would be more like taking some money out of your account and still having the full value at your disposal (provided you don’t spend any of it anyway). It is just in two places now, not one. The total will be the same.
No matter what low-pass filter is used to reveal to the eye the low-frequency signal components, it doesn’t change the fact that no sound basis for predicting their future is provided thereby. Unlike the tides, temperature variations are NOT simply the superposition of discrete sinusoids (line spectra); they are stochastic processes whose spectrum may contain fairly narrow-band components, but encompasses a wide band CONTINUUM of frequencies nevertheless. Even in the multidecadal range, GISP2 data show significant peaks at periods considerably different than ~60yrs. Only amateurs would contend that the precise nature of the LP filter makes a critical difference in the practical problem. Ironically enough, the great side-lobe leakage of SG filters is blithely ignored here, while arguing the “superiority” of cascaded boxcar or gaussian filters.
cd: What is it that I am not explaining well?
A simple low pass filter (Gaussian or whatever) is a binary chop into two bands, high and low pass. Together they sum to the original series. There is no loss (assuming that we do not get rounding errors in the calculations anyway which in this case is not even going to be close to being a factor).
This ‘filters must lose information’ is a myth. Please tell me where you think the loss is going?
1sky1 says:
March 17, 2014 at 5:27 pm
“Ironically enough, the great side-lobe leakage of SG filters is blithely ignored here, while arguing the “superiority” of cascaded boxcar or gaussian filters.”
Well as the output of the S-G can be directly compared to the CTRM where they overlap then any such error should be easy to see. 🙂
cd, Bernie et al.
Information is only lost if you throw it away. Filtering does not throw anything away, it divides it. Richard has separated out the long term components to examine that part. Like separating out the apples to weight them. I don’t see any suggestion from him that this should now be the new surface temperature record and we should throw away the rest.
The way I’d look at this would be to try to isolate the supposed 60 y periodicity and subtract it and see what’s left. Then do a similar thing again. If there’s something around 34 years 😉 subtract that. By this time the data should be long enough to do a reliable FFT on what remains.
Now that would not be possible if the information was “lost.”
Enough of this rather pointless argument, please.
RichardLH
Pray tell, where does that loss go? It must go somewhere.
Now if ever there was red flag that it is it. It is transferred, depending to the coefficients (and even if dynamic) information is transferred to one part of the kernel window (data point) to another (typically THE kernel cell).
Do you know how convolution work? No honestly!
A cascaded filter set as described allows for appropriate correction of the previous stage until you get to the point where any errors that are left are below the input digitisation ones.
That is the same principle of a B-Spline, and is based in the prescribed level of the spline.
Greg
I appreciate your response: but information is lost?! All you’re doing is squirming, you now seem to be moving toward a transform argument; that we move toward a kernel with a series of operations with widths 1 to n-1 then we can reconstruct the signal. I don’t think this is even right, perhaps I misunderstand but remember if I got to implement it in code then waffle wont do. Give me an algorithm in pseudocode!
I actually like Richard’s argument, but your response to Willis’s tight statement (and yes in my rather modest understanding (but I’m no novice either), shows that you guys are not so sure of your own turf but rather than coming clean you obfuscate – the ultimate condescension.
cd “One should judge the weight of someones understanding by the appropriateness and clarity of their analogy. This however, is rather pathetic.”
Well you say it’s pathetic and start getting all philosophical about the worth of a man from his analogies etc. but you don’t have any counter to the ‘pathetic’ analogy.
Nothing is lost by filtering.
Quite what you think I should do with gimp or what bearing this has on whether something is lost when using a low pass filter to split a dataset into h.f and l.f bands, is not clear.
Enough of this rather pointless argument. I’m out;.
cd says:
March 17, 2014 at 5:32 pm
“Now if ever there was red flag that it is it. It is transferred, depending to the coefficients (and even if dynamic) information is transferred to one part of the kernel window (data point) to another (typically THE kernel cell).”
As I pointed out above and more than once. There will be some frequencies around the corner value where some portion will end up in the wrong ‘bin’. Not all, just some. And only for a small range of frequencies.
As I deliberately chose a portion of the spectrum where there is little ‘power’ (i.e. did the sweep needed to check the output sensitivity to corner values) then any such values are going to be very small in any case.
In any case Gaussian is the least likely to cause any in-band re-assignments. It is considered to be the ‘ideal’ time domain filter which is what is required.
So I ask again, what are the losses and where are they going? For a 15 year corner Gaussian (or near Gaussian) filter specifically.
Greg says:
March 17, 2014 at 5:31 pm
cd, Bernie et al.
“Information is only lost if you throw it away. Filtering does not throw anything away, it divides it. Richard has separated out the long term components to examine that part. Like separating out the apples to weight them. I don’t see any suggestion from him that this should now be the new surface temperature record and we should throw away the rest.
The way I’d look at this would be to try to isolate the supposed 60 y periodicity and subtract it and see what’s left. Then do a similar thing again. If there’s something around 34 years 😉 subtract that. By this time the data should be long enough to do a reliable FFT on what remains.
Now that would not be possible if the information was “lost.”
Enough of this rather pointless argument, please.”
Yes – let’s stop being silly. It’s not lost because it WAS lost but you bring it back when you subtract.
Recall all the students back to class. Tell them we have been lying for decades. Always a good idea to keep them off balance by changing the rules.
Greg
Nothing is lost by filtering.
OK. And again I ask you:
1) Apply filter to image.
2) With only the filtered image and the filter kernel design.
3) Give me the original image.
I use Gimp cause others can then try your suggestions (without expense), after all no information is lost so it should be simple to do. Please tell us how to get go through steps 1 to 3 without any knowledge of other steps other than the filter design?