One of the comments on my post on WattsUp suggested I attempt to do an SSA on the Hadcrut3 data. Here are the results:
(Note: higher resolution plots may be viewed by clicking on the plots below)

The best signal-to noise ratio was obtained with a window-size (L) of 20 and an eigenspan of 4 (2 eigenpairs). Low frequency oscillatory behavior is evident in the spectral plot which show the 2 lowest eigenvectors significantly above the noise floor.

Looking at the FFT of the model, we see again similar harmonically related spectral peaks to those found via harmonic decomposition.


Proceeding as before, we find the frequency of the psd spectral peaks
w1 = ArgMax[psdMod, w];
w2 = ArgMax[{psdMod, w > 1.5 w1}, w];
2*Pi/{w1, w2}
{185.83453775994096, 60.87208526557042}

which are fairly close to the values found via HD. Note the ratio = 3.05. Assuming we have more accuracy on the higher frequency mode, and assuming a harmonic ratio of three we form the model:

AY1 sin(phiY1 + 0.0344065 t) + AY2 sin(phiY2 + 0.103219 t)

where the radian frequencies {0.103219, 0.0344065} are w2 and w2/3 from above, and fit the model to the data to solve for the phase and amplitudes.

fitParmsSSA =
FindFit[yearly, modelSSA, varsSSA, t, Method -> NMinimize]
fitSSA = Table[modelSSA /. %, {t, 1, 113}];
residualSSA = yearly - fitSSA;
DateListPlot[TemporalData[residualSSA, Automatic]["Path"],
Joined -> True]
ListPlot[{yearly, fitSSA}, Joined -> True]

{AY1 -> -0.261994, phiY1 -> 1.38923, AY2 -> -0.160245,
phiY2 -> 0.464517}

The figures below shows the hybrid model (that is an HD model derived from the SSA eigenpairs as described above) versus the monthly Hadcrut3 data, hindcasted 25 and 50 years (to 1875 and 1850), and in the second figure forecasted 50 years to 2063. The hindcast matches the data closely for 23 years or so. One interesting feature that the event of 1878 is co-incident with the data matching the model, and the super-el-nino event of 1997 is co-incident with divergence. One highly speculative explanation (besides mere coincidence, or events like these causing issues in the SSA algorithm)) is the possibility that climate disruptions triggered a new phase state in the climate dynamics.


As with HD, I find it surprising that a fit can be achieved with so few eigenpairs (or in the case of HD, harmonic pairs). The fact that the SSA modes are also harmonically related (there is nothing in the SSA method that requires them to be so), adds credence to the notion that harmonically related periodicities are part of the current climate dynamics.

There are forecasting methods which use only the SSA results. Once I figure them out and code the algorithms I’ll post an update.

Part 2 is available here.