3.2 Spectral methods

“...the use of any mathematical algorithm to derive hidden periodicities from the data always entails the question as to whether the resulting cycles are not introduced either by the particular numerical method used or by the time interval analyzed.”
(de Meyer, 1981)

Spectral analysis of the sunspot number record is used for prediction under the assumption that the main reason of variability in the solar cycle is a long-term modulation due to one or more periods.

The usual approach to the problem is the purely formal one of representing the sunspot record with the superposition of eigenfunctions forming an orthogonal basis. From a technical point of view, spectral methods are a complicated form of linear regression. The analysis can be performed by any of the widely used means of harmonic analysis:

(1) Least squares (LS) frequency analysis (sometimes called “Lomb–Scargle periodogram”) consists in finding by trial and error the best fitting sine curve to the data using the least squares method, subtracting it (“prewhitening”), then repeating the procedure until the residuals become indistinguishable from white noise. The first serious attempt at sunspot cycle prediction, due to Kimura (1913), belonged to this group. The analysis resulted in a large number of peaks with dubious physical significance. The prediction given for the upcoming cycle 15 failed, the forecasted amplitude being  60 while the cycle actually peaked at 105. However, it is interesting to note that Kimura correctly predicted the long term strengthening of solar activity during the first half of the 20th century! LS frequency analysis on sunspot data was also performed by Lomb and Andersen (1980), with similar results for the spectrum.

(2) Fourier analysis is probably the most commonly used method of spectral decomposition in science. It has been applied to sunspot data from the beginning of the 20th century (Turner, 1913a,bMichelson, 1913). Vitinsky (1973) judges Fourier-based forecasts even less reliable than LS periodogram methods. Indeed, for instance Cole (1973) predicted cycle 21 to have a peak amplitude of 60, while the real value proved to be nearly twice that.

(3) The maximum entropy method (MEM) relies on the Wiener–Khinchin theorem that the power spectrum is the Fourier transform of the autocorrelation function. Calculating the autocorrelation of a time series for points and extrapolating it further in time in a particular way to ensure maximal entropy can yield a spectrum that extends to arbitrarily low frequencies despite the shortness of the data segment considered, and also has the property of being able to reproduce sharp spectral features (if such are present in the data in the first place). A good description of the method is given by Ables (1974), accompanied with some propaganda for it – see Press et al. (1992) for a more balanced account of its pros and cons. The use of MEM for sunspot number prediction was pioneered by Currie (1973). Using maximum entropy method combined with multiple regression analysis (MRA) to estimate the amplitudes and phases, Kane (2007) arrived at a prediction of 80 to 101 for the maximum amplitude of cycle 24. It should be noted that the same method yielded a prediction (Kane, 1999) for cycle 23 that was far off the mark.

(4) Singular spectrum analysis (SSA) is a relatively novel method for the orthogonal decomposition of a time series. While in the methods discussed above the base was fixed (the trigonometric functions), SSA allows for the identification of a set of othogonal eigenfunctions that are most suitable for the problem. This is done by a principal component analysis of the covariance matrix . SSA was first applied to the sunspot record by Rangarajan (1998) who only used this method for pre-filtering before the application of MEM. Loskutov et al. (2001) who also give a good description of the method, already made a prediction for cycle 24: a peak amplitude of 117. More recently, the forecast has been corrected slightly downwards to 106 (Kuzanyan et al., 2008).

The dismal performance of spectral predictions with the methods (1) – (3) indicates that the sunpot number series cannot be well represented by the superposition of a limited number of fixed periodic components. Instead,

• the periods may be time dependent,
• the system may be quasiperiodic, with a significant finite width of the periodic peaks (esp. the 11-year peak),
• there may be non-periodic (i.e., chaotic or stochastic) components in the behaviour of the system, manifested as a continuous background in the spectrum.

In practice, all three effects suggested above may play some part. The first mentioned effect, time dependence, may in fact be studied within the framework of spectral analysis. MEM and SSA are intrinsically capable of detecting or representing time dependence in the spectrum, while LS and Fourier analysis can study time dependence by sliding an appropriate data window across the period covered by observations. If the window is Gaussian with a width proportional to the frequency we arrive at the popular wavelet analysis. This method was applied to the sunspot number series by Ochadlick Jr et al. (1993), Vigouroux and Delachie (1994), Frick et al. (1997), Fligge et al. (1999), and Li et al. (2005) who could confirm the existence and slight variation of the 11-year cycle and the Gleissberg-cycle. Recently, Kolláth and Oláh (2009) called attention to a variety of other generalized time dependent spectral analysis methods, of which the pseudo-Wigner transform yields especially clear details (see Figure 9). The time varying character of the basic periods makes it difficult to use these results for prediction purposes but they are able to shed some light on the variation as well as the presistent or intermittent nature of the periods determining solar activity.

In summary, it is fair to say that forecasts based on harmonic analysis are notoriously unreliable. The secular variation of the basic periods, obeying as yet unknown rules, would render harmonic analysis practically useless for the prediction of solar cycles even if solar activity could indeed be described by a superposition of periodic functions. Although they may be potentially useful for very long term prediction (on centennial scales), when it comes to cycle-to-cycle forecasts the best we can hope from spectral studies is apparently an indirect contribution, by constraining dynamo models with the inambiguously detected periodicities.

In what remains from this subsection, we briefly review what these apparently physically real periods are and what impact they may have on solar cycle prediction.

3.2.1 The 11-year cycle and its harmonics

As an example of the period spectrum obtained by these methods, in Figure 8 we present the FFT based power spectrum estimate of the smoothed sunspot number record. Three main features are immediately noticed:

• The dominant 11-year peak, with its sidelobes and its 5.5-year harmonic.
• The 22-year subharmonic, representing the even–odd rule.
• The significant power present at periods longer than 50 years, associated with the Gleissberg cycle.

The dominant peak in the power spectrum is at  11 years. Significant power is also present at the first harmonic of this period, at 5.5 years. This is hardly surprising as the sunspot number cycles, as presented in Figure 3, have a markedly asymmetrical profile. It is a characteristic of Fourier decomposition that in any periodic series of cycles where the profiles of individual cycles are non-sinusoidal, all harmonics of the base period will appear in the spectrum.

Indeed, were it not for the 13-month smoothing, higher harmonics could also be expected to appear in the power spectrum. It has been proposed (Krivova and Solanki, 2002) that these harmonics are detected in the sunspot record and that they may be related to the periodicities of  1.3 years intermittently observed in solar wind speed (Richardson et al., 1994Paularena et al., 1995Szabo et al., 1995Mursula and Zieger, 2000Lockwood, 2001) and in the internal rotation velocity of the Sun (Howe, 2009, Sect. 10.1). An analoguous intermittent 2.5 year variation in the solar neutrino flux (Shirai, 2004) may also belong to this group of phenomena. It may be worth noting that, from the other end of the period spectrum, the 154-day Rieger period in solar flare occurrence (Rieger et al., 1984Bai and Cliver, 1990) has also been tentatively linked to the 1.3-year periodicity. Unusually strong excitation of such high harmonics of the Schwabe cycle may possibly be explained by excitation due to unstable Rossby waves in the tachocline (Zaqarashvili et al., 2010).

The 11-year peak in the power spectrum has substantial width, related to the rather wide variation in cycle lengths in the range 9 – 13 years. Yet Figure 8 seems to suggest the presence of a well detached second peak in the spectrum at a period of  14 years. The presence of a distinct peak at the first harmonic and even at the subharmonic of this period seems to support its reality. Indeed, peaks at around 14 and 7 years were already found by other researchers (e.g., Kimura, 1913Currie, 1973) who suggested that these may be real secondary periods of sunspot activity.

The situation is, however, more prosaic. Constraining the time interval considered to data more recent than 1850, from which time the sunspot number series is considered to be more reliable, the 14.5-year secondary peak and its harmonics completely disappear. On the other hand, the power spectrum for the years 1783 – 1835 indicates that the appearance of the 14.5-year secondary peak in the complete series is almost entirely due to the strong predominance of this period (and its harmonic) in that interval. This interval covers the unusually long cycle 4 and the Dalton minimum, consisting of three consecutive unusually weak cycles, when the “normal” 11-year mode of operation was completely suppressed.

As pointed out by Petrovay (2010), this probably does not imply that the Sun was operating in a different mode during the Dalton minimum, the cycle length being 14.5 years instead of the usual 11 years. Instead, the effect may be explained by the well known inverse correlation between cycle length and amplitude, which in turn is the consequence of the strong inverse correlation between rise rate and cycle amplitude (Waldmeier effect), combined with a much weaker or nonexistent correlation between decay rate and amplitude (see Section 1.3.3). The cycles around the Dalton minimum, then, seem to lie at the low amplitude (or long period) end of a continuum representing the well known cycle length–amplitude relation, ultimately explained by the Waldmeier effect.

A major consequence of this is that the detailed distribution of peaks varies significantly depending on the interval of time considered. Indeed, Kolláth and Oláh (2009) recently applied time dependent harmonic analysis to the sunspot number series and found that the dominant periods have shown systematic secular changes during the past 300 years (Figure 9). For instance, the basic period seems to have shortened from 11 years to 10 years between 1850 and 1950, with some moderate increase in the last 50 years. (This is consistent with the known anticorrelation between cycle length and amplitude, cf. Section 1.3.3.)

3.2.2 The even–odd (a.k.a. Gnevyshev–Ohl) rule

A cursory look at Figure 3 shows that solar cycles often follow an alternating pattern of higher and lower maxima. In this apparent pattern, already noticed by the early observers (e.g, Turner, 1913c), odd cycles have been typically stronger than even cycles in the last two centuries.

This even–odd rule can be given two interpretations: a “weak” one of a general tendency of alternation between even and odd cycles in amplitude, or a “strong” one of a specific numerical relation between the amplitudes of consecutive cycles.

Let us first consider the rule in its weak interpretation. At first sight the rule admits many exceptions, but the amplitude of solar cycles depends on the particular measuring method used. Exceptions from the even–odd alternation rule become less common if a long term trend (calculated by applying a 12221 or 121 filter, see Section 1.3.1) is subtracted from the data (Charbonneau, 2001), or if integrated cycle amplitudes (sums of annual mean sunspot numbers during the cycle) are used (Gnevyshev and Ohl, 1948).

In fact, as evident from, e.g., the work of Mursula et al. (2001) where cycle amplitudes are based on group sunspot numbers and the amplitude of a cycle is defined as the sum of the annual GSN value over the course of the cycle, the odd–even alternation may be considered as strictly valid with only four exceptions:

• In the pairs 7 – 8 and 17 – 18, odd cycles are followed by stronger even cycles at the end of Dalton minimum and at the beginning of the Modern Maximum. These exceptions could be made to disappear by the subtraction of the long term trend as suggested by Charbonneau (2001).
• The pair 22 – 23 represents another apparent break of the weak even–odd rule which is not easily explained away, even though the relative difference is smaller if the Kislovodsk sunspot number series is used (Nagovitsyn et al., 2009). The possibility is obviously there that the subtraction of the long term trend may resolve the problem but we have no way to tell in the near future.
• Prior to cycle 5, the phase of the alternation was opposite, even cycles being stronger than odd cycles. As cycle 4 is known to have been anomalously long anyway (the so-called “phase catastrophe” in the solar cycle, Vitinsky et al., 1986) and its decaying phase is not well covered by observations (Vaquero, 2007), this gave rise to the suggestion of a “lost solar cycle” between cycles 4 and 5 (Usoskin et al., 2001Usoskin and Mursula, 2003). This cycle, however, would have been even more anomalous than cycle 4 and despite intensive searches in historic data the evidence is still not quite conclusive (Krivova et al., 2002; see, however, Usoskin et al., 2009a).

The issue whether the even–odd rule can go through phase jumps or not is important with respect to its possible origin. One plausible possibility is that the alternation is due to the superposition of a steady primordial magnetic field component on the oscillatory magnetic field generated by the dynamo (Levy and Boyer, 1982). In this case, any phase jump in the Gnevyshev–Ohl rule should imply a phase jump in Hale’s polarity rules, too. Alternatively, persistent even–odd alternation may also arise in nonlinear dynamos as a period–2 limit cycle (Durney, 2000); with a stochastic forcing occasional phase jumps are also possible (Charbonneau, 2001Charbonneau et al., 2007).

While we have no information on this from the 18th century phase jump, we can be certain that there was no such phase jump in polarities in the last two decades, even though the even–odd rule seems to have been broken again. It will be interesting to see when (and if) the even–odd rule settles in again, whether it will have done so with a phase jump or not. For instance, if cycle 25 will again exceed cycle 24 it would seem that no phase jump occurred and both theoretical options are still open. But if cycle 25 will represent a further weakening from cycle 24, followed by a stronger cycle 26, a phase jump will have occurred, which may exclude the primordial field origin of the rule if Hale’s polarity rules remain unchanged.

Let us now discuss the stronger interpretation of the even–odd rule. In the first quantitative study of the relative amplitudes of consecutive cycles, Gnevyshev and Ohl (1948) found a rather tight correlation between the time integrated amplitudes of even and subsequent odd cycles, while the correlation between odd cycles and subsequent even cycles was found to be much less strong. This gave rise to the notion that solar cycles come in “two-packs” as even–odd pairs. Nagovitsyn et al. (2009) confirmed this puzzling finding on the basis of data covering the whole period of telescopic observations (and renumbering cycles before 1790 in accordance with the lost cycle hypothesis); they also argue that cycle pair 22 – 23 does not deviate strongly from the even–odd correlation curve so it should not be considered a “real” exception to the even–odd rule.

The fact that shortly after its formulation by Gnevyshev and Ohl (1948), the (strong) even–odd rule was used by Kopecký (1950) to successfully predict the unusually strong cycle 19 made this method particularly popular for forecast purposes. However, forecasts based on the even–odd rule completely failed for cycle 23, overpredicting the amplitude by > 50% (see review by Li et al., 2001). Taken together with the implausibility of the suggested two-pack system, this shows that it is probably wiser to take the position that “extraordinary claims need extraordinary evidence” – which is yet to be provided in the case of the “strong” even–odd rule.

Finally, in the context of the even–odd rule, it is also worth mentioning the three-cycle regularity proposed by Ahluwalia (1998). Even though the evidence presented for the alleged triadic pattern is not overwhelming, this method resulted in one of the few successful predictions for the amplitude of cycle 23.

3.2.3 The Gleissberg cycle

Besides the changes in the length of the 11-year cycle related to the amplitude–cycle length correlation, even more significant are the variations in the period of the so-called Gleissberg cycle (Gleissberg, 1939). This “cycle”, corresponding to the 60 – 120 year “plateau” in Figure 8 was actually first noticed by Wolf, who placed it in the range 55 – 80 years (see Richard, 2004, for a discussion of the history of the studies of the Gleissberg cycle). Researchers in the middle of the 20th century characterized it as an 80 – 100 year variation. Figure 9 explains why so widely differing periods were found in different studies: the period has in fact shown a secular increase in the past 300 years, from about 50 years in the early 18th century, to a current value exceeding 140 years. This increased length of the Gleissberg cycle also agrees with the results of Forgács-Dajka and Borkovits (2007).

The detection of  100 year periods in a data set of 300 years is of course always questionable, especially if the period is even claimed to be varying. However, the very clear and, most importantly, nearly linear secular trend seen in Figure 9 argues convincingly for the reality of the period in question. This clear appearance of the period is due to the carefully optimized choice of the kernel function in the time–frequency analysis, a method resulting in a so-called pseudo-Wigner distribution (PWD). In addition, in their study Kolláth and Oláh (2009) present an extremely conscientious test of the reliability of their methods, effectively proving that the most salient features in their PWD are not artefacts. (The method was subsequently also applied to stellar activity, Oláh et al., 2009.) This is the most compelling evidence for the reality of the Gleissberg cycle yet presented.

3.2.4 Supersecular cycles

For the 210-year Suess cycle, McCracken and Beer (2008) present further evidence for the temporally intermittent nature of this marked peak in the spectrum of solar proxies. The Suess cycle seems to have a role in regulating the recurrence rate of grand minima. Grand minima, in turn, only seem to occur during < 1 kiloyear intervals (“Spörer events”) around the minimum of the  2400-year Hallstatt cycle.

For further discussion of long term variations in solar activity we refer the reader to the reviews by Beer et al. (2006) and Usoskin (2008).