text:
# 1 Introduction

## 1 Introduction

The widely used Fourier transform, although useful for stationary signals of simple dynamics consisted of a linear superposition of few independent, strong, non-evolving periodicities, has severe drawbacks for analyzing signals of the following two important categories, often encountered as output of physically interesting complex systems:

a. Signals that include transient or variable periodicities. The Fourier transform, since integrating over the whole time domain, does not allow monitoring of amplitude or frequency evolution in periodicities. A limited solution is to divide the signal in segments or else ’windows’ (a method called short-time, Gabor or windowed Fourier transform), however introducing the arbitrary, fixed window width that imposes a low frequency cut-off near which edge effects occur, thus artificially deforming the estimated spectrum.

b. Signals that significantly depart from stationarity, consisted of intermittent ’activity bursts’ (e.g. energy releases), or hierarchies of localized structures in space or time (e.g. in turbulent fluids reside eddies, discontinuities, filaments, sheets, shocks etc.). For each localized structure a wide frequency range of strongly phase coherent (i.e. ’synchronized’ in time) Fourier sine components are required to reproduce it, adding near and canceling away from the structure. Hence, each structure spreads and become undetectable over the whole Fourier spectrum. Moreover, such structures often exhibit ’self-similarity’, or else ’scale invariance’, i.e. consists of hierarchies of mutually similar structures, emerging through a cascade mechanism, and ’roughness’ over a wide scale range (well modeled by fractals and multi-fractals) that is not built in the sine waves of Fourier transform.

The projection of signals with continuous or discrete wavelet analyzing functions, that are by definition self-similar and well localized in both time (or space) and frequency (or wavenumber) domains, provides efficient descriptions, in terms of sparseness, of non-stationary, transient, or ’bursty’ records. It is proven a useful tool, used in an increasingly large variety of applications. Here we mention some of them indicatively according to the type of use: detection of periodic signals in noisy timeseries (e.g. Otazu et al., 2002), monitoring of period variations (e.g. Frick et al., 1997, Fligge et al., 1999), unevenly sampled period analysis (e.g. Foster, 1996), derivation of multifractal properties (e.g. Arneodo et al., 1988, Argoul et al., 1989, Muzy et al., 1991, 1993, Arneodo et al., 1995, 1998), random multifractal synthesis (e.g. Benzi etal, 1993), localized structure identification (e.g. Lucek & Balogh, 1997, Roux et al., 1999, Mouri et al., 1999), polarization analysis (e.g. Baumjohann et al., 1999), denoising (e.g. Fligge & Solanki, 1997, Komm et al., 1999), signal-to-noise ratio enhancement (e.g. Zhang & Paulson, 1997), information compression (e.g. Muhlmann & Hanslmeier, 1996), multi-resolution image decomposition (e.g. Mallat, 1989, Pantin & Starck, 1996), studies of random walks (Arneodo et al., 1996), flow structure analysis (Haynes & Norton, 1993), population classification (Bendjoya et al., 1991), clustering detection (e.g. Bijaoui et al., 1993, Girardi et al., 1997, Lima Neto et al., 1997), fast data query (e.g. Chakrabarti et al., 2001), fluid turbulence simulations (e.g. Schneider et al., 1997) etc.

The wavelet transform is now used extensively for the analysis of oscillations and periodicities in various solar structural and activity features, as the sunspot number (Fligge et al., 1999, Mordvidov & Kuklin, 1999, Sello, 2000, Rigozo et al., 2001, Nayar et al., 2002), sunspot groups (Ballester & Oliver, 1999), active regions (Ireland et al., 1999, O’Shea et al., 2001, 2002), coronal holes (Marsh et al., 2002), bursts (Schwarz et al., 1998, Meszarosova et al., 1999, Gimenez et al., 2001), flares (Aschwanden et al., 1998), polar plumes (Banerjee et al., 2000), coronal loops (De Moortel et al., 2000, 2002a, b), photospheric flows (Lawrence et al., 2001), chromosphere (Bocchialini & Baudin, 1995, Frick et al., 1997), transition region (Fludra, 2001), differential rotation (Hempelmann & Donahue, 1997, Soon et al., 1999, Hempelmann, 2002) etc.

In this article we demonstrate that two exactly defined, mutually independent, complementary wavelet skeleton spectra, extracted from the total continuous wavelet spectrum, allow to discriminate the components of (stable or transient) periodicities and of hierarchies of discontinuities in timeseries. An indicative application is made to sunspot number, since (a) is known to include both periodic and noise-like components, and (b) its study may contribute in understanding the features and mechanism of the solar activity. Section 2 includes a brief review of the continuous wavelet transform for the unfamiliar reader (for extended reviews see e.g. Hunt et al., 1993, Kaiser, 1994, Louis et al., 1997, Torrence & Compo 1998). The two wavelet skeleton spectra are introduced in Section 3. The application to the sunspot index is presented in Section 4. Conclusions and discussion are summarized in Section 5.

## 2 A brief review of the wavelet transform

### 2.1 Linear functional projection

One way to quantify the degree of ’similarity’ between two, generally complex function $`x(t),y(t)`$ of a variable $`t[0,T]`$ (e.g. time or distance) is by the amplitude (or else cross-correlation) of their functional projection over the domain of $`t`$:

$$a=\frac{1}{T}_0^Tx(t)y^{}(t)𝑑t$$
(1)
with $`y^{}`$ being the complex conjugate of $`y`$. The factor $`1/T`$ makes the amplitude $`a`$ independent of $`T`$ and finite in the limit $`T+\mathrm{}`$. The two functions are called ’orthogonal’ if $`a=0`$.

In analogy to the projection of a vector to a coordinate system in physical space, a given function $`x(t)`$ (e.g. the timeseries of a measured observable) thus can be analyzed to a family of functions $`y(t,p_i)`$, of $`n`$ characteristic parameters $`p_i,i=1\mathrm{}n`$:

$$a(p_i)=\frac{1}{T}_0^Tx(t)y^{}(t,p_i)𝑑t$$
(2)
The functional analysis $`(\text{2})`$ is linear in the sense that any $`x(t)`$ can be expanded to a linear superposition of a set of functions $`y(t,p_i)`$ that are mutually orthogonal and integrable.

Note that the continuous functional transform $`(\text{2})`$ is generally redundant since it spreads the information content of $`x(t)`$ from the one-dimensional time axis, to the n-dimensional space of the parameters $`p_i`$. The minimization of information spread in the parameter space $`p_i`$ by appropriate selection of the functional basis and the development of efficient redundancy reduction schemes is of central importance for spectral analysis, functional approximation and compression applications.

### 2.2 Appropriate selection of analyzing functions

Although the linear functional projection $`(\text{2})`$ is mathematically valid with any set of functions, an appropriate projection should aim in selecting such functions to match as much as possible (i.e. being nearly orthogonal to most of the functions so that the representation is sparse) the analyzed function, which often is or consists of the following:

(a) Analog signal, as the independent variables of a system or model of relatively simple dynamics. The usually smoothly varying variables of systems of low dimensionality (i.e. small number of independent variables) exhibit oscillatory behavior, thus can often be sparsely expanded to a linear superposition of a few periodic functions, or, more generally, transient, quasi-periodic components. Hence, it is desirable that the analyzing functions have well-defined, single-peaked, narrow frequency spectral content.

(b) Noise, as this due to observational ’random’ fluctuations (often additive) of high-frequency content, e.g. of flat spectrum in uncorrelated (white) or more generally of wide range power-law decaying spectrum in auto-correlated (’colored’) noise. Noise-like contributions may also caused by fast-varying components of a deterministic complex system of large number of independent variables (such as macroscopic systems) or the interaction of a low-dimensional open system with an environment of high dimensionality. The observables of high-dimensional systems, often involving strongly turbulent fluids, are characterized by intermittency or ’singularities’ due to coherent structures or intermittent energy bursts and, exact or weaker, self-similarity between different scales (often termed scale-invariance), well modeled by fractals and multi-fractals. Although such systems are characterized by strongly non-linear physical mechanisms, for an appropriate linear projection of such observables, the analyzing functions should at least posses the element of self-similarity.

### 2.3 Fourier transform

The widely used Fourier transform can be interpreted as mapping, or complete ’rotation’ (in analogy to the vector rotation in physical space) from the time domain to the frequency domain, by means of the complex periodic plain wave functions:

$$y(t,f)=e^{2\pi ift}$$
(3)
The (generally complex) spectral amplitude (see equ. $`(\text{2})`$):

$$a(f)=\frac{1}{T}_0^Tx(t)y^{}(t,f)𝑑t$$
(4)
thus represents an estimate of the amplitude of periodicities of frequency $`f`$, that occurred during the time interval $`T`$. The, often studied, Fourier power spectral density is defined as follows:

$$p(f)=a^2(f)$$
(5)
A useful property of the Fourier transforms is the power theorem: any linear functional projection $`(\text{2})`$ can be calculated via Fourier transforms in the frequency domain (e.g. Hunt et al., 1993, p.6):

$$a(p_i)=_0^Tx(t)y^{}(t,p_i)𝑑t=_{\mathrm{}}^+\mathrm{}\widehat{x}(f)\widehat{y}^{}(f,p_i)𝑑f$$
(6)
where $`\widehat{x},\widehat{y}`$ the Fourier transforms of $`x,y`$ respectively (see $`(\text{4})`$). Hence, the functional projection $`(\text{2})`$ can be treated as a linear filter acting on $`\widehat{x}`$.

### 2.4 Time-frequency analysis

Note that each Fourier basis function $`y`$ is totally localized in frequency but nowhere in time, since having infinite duration (i.e. $`y`$ in $`(\text{3})`$ do not vanish for $`t\pm \mathrm{}`$). As a consequence, the Fourier transform does not allow monitoring of transient (periodic or aperiodic) components in the signal, since they are spread over the whole Fourier spectrum.

In the framework of the transform $`(\text{2})`$ a solution to this problem is to introduce analysing functions of well-defined frequency but are also localized in time i.e. of significant amplitude for a finite time duration.

Perfect localization in both time and frequency in however inhibited by the uncertainty principle (e.g. Kaiser, 1994, p. 50): $`(\text{2})`$ keeps the information of $`x(t)`$ confined in cells, within which there is total uncertainty about the information of $`x(t)`$, and the size of these cells cannot be infinitesimally small, having a lower limit depending on the properties (shape, rate of decay etc.) of the selected analysing functions.

Any time-frequency transform represents a mapping of $`x(t)`$ on the $`tf`$ plane. Given a set of (generally complex-valued) functions $`y(t,t^{},f)`$, of characteristic frequency $`f`$, localized in the vicinity of time $`t^{}`$, the linear functional projection $`(\text{2})`$ is:

$$a(t,f)=\frac{1}{T}_0^Tx(t^{})y^{}(t,t^{},f)𝑑t^{}$$
(7)
The graph of the power $`a^2(t,f)`$ on the $`tf`$ plane is oftenly called scalogram. The total spectrum:

$$p(f)=_0^Ta^2(t,f)𝑑t$$
(8)
is expected to generally exhibit similar features and shape (e.g. spectral peaks, slope) as the corresponding Fourier spectrum (equ. $`(\text{5})`$).

### 2.5 Windowed Fourier transform

In the windowed Fourier transform (Gabor, 1946) the Fourier plain wave functions are modulated by a ’support’ $`\varphi `$ (often called ’window’ function) of significant value only in the vicinity $`\sigma `$ of time $`t^{}`$:

$$y(t,t^{},f)=e^{2\pi ift}\varphi (\frac{tt^{}}{\sigma })$$
(9)
Its Fourier transform has spectral width of the order of $`1/\sigma `$. Therefore, the $`tf`$ plane is resolved in cells of fixed shape and area at least $`\sigma `$ by $`1/\sigma `$.

A commonly used support function is a Gaussian of width $`\sigma `$:

$$\varphi (\frac{tt^{}}{\sigma })=e^{\frac{(tt^{})^2}{2\sigma ^2}}$$
(10)
The resulting analysing function $`(\text{9})`$ has Fourier transform:

$$\widehat{y}(f,f^{})=\sigma e^{\frac{(ff^{})^2\sigma ^2}{2}}$$
(11)
The Gaussian support can be proven to give optimal localization in the sense that among all functions it covers smallest area permitted by the uncertainty principle in the time-frequency plane (Louis et al., 1997). Other window functions are also used, aimed at reducing the leakage of power of spectral peaks to nearby frequencies (e.g. Press et al., 1992).

However the windowed Fourier transform always introduces the arbitrary duration of the support function ($`\sigma `$ in equ. $`(\text{10})`$) which itself does not depend on $`f`$ and must be selected appropriately before the analysis. Also, the analyzing functions $`y`$ are not self-similar since include increasingly many oscillations within their support with increasing frequency.

### 2.6 Continuous wavelet transform

The continuous wavelet transform resolves the information on the $`tf`$ plane in cells of variable shape depending on frequency, with size of the order of $`1/f`$ by $`f`$. With this partition of the $`tf`$ plane the arbitrary selection of the windowed Fourier transform is avoided and the analyzing functions (wavelets) are exactly self-similar (or else scale-invariant) since each includes a constant number of oscillations independent of frequency. The wavelet transform at lower frequencies provides better resolution in frequency but worse in time (since the analysing functions are wider) and at high frequencies better time localization (since the analysing functions are narrower) but more uncertainty in frequency, in all acting as a ’mathematical microscope’.

While many wavelet families are so far proposed (also depending on the nature of application), of the most commonly used in physical applications is the Morlet wavelet (Morlet et al., 1982), being the generalization of the windowed Fourier transform (see equ. $`(\text{9})`$):

$$y(t,t^{},f)=e^{2\pi ift}\varphi (f(tt^{}))$$
(12)
with Gaussian support (in analogy with $`(\text{10})`$) of width $`1/f`$:

$$\varphi (f(tt^{}))=e^{\frac{f^2(tt^{})^2}{2}}$$
(13)
Its Fourier transform is also Gaussian of width $`f`$:

$$\widehat{y}(f,f^{})=\frac{1}{f}e^{\frac{(ff^{})^2}{2f^2}}$$
(14)
Among all possible wavelet analysing functions, the Morlet wavelet due to its Gaussian support inherits optimality as regarding the uncertainty principle (Louis et al., 1997).

## 3 The two wavelet skeleton spectra

### 3.1 Signal-noise decomposition:

Successful decomposition of the timeseries into ’signal’ and noise-like component is of fundamental importance in various applications, for the separate study of the signal and noise properties (e.g. locating spectral peaks and deriving scaling laws), denoising, compression etc. The wavelet transform is promising for that decomposition, since allowing monitoring of the evolution of time-localized wavepackets of specific frequency, hence also of stable or transient periodicities and noise-like spikes.

The overall wavelet scalogram generally contains a mixture of the signal and noise characteristics. Imposing some thresholding logic (e.g. hard or soft level, trend etc) to the wavelet components (i.e. a type of wavelet band-pass filter) is expected to have limited success in the case of strong, auto-correlated noise and/or weak, transient periodicities. The deeper reason that such filters have limited success in those cases is that the selection of the wavelet components that are attributed to the signal or the noise is to some extend arbitrary, not adaptively taking into account the information content in the original timeseries.

As commented in section $`(\text{2.1})`$ for any continuous functional projection, the continuous wavelet transform is always redundant in both time and frequency, since mapping the one dimensional information in $`x(t)`$ on the (two dimensional) $`tf`$ plane. It contains correlations between the wavelet components that do not exist in the timeseries but is significant at those (small) regions and scales where the wavelet functions cross-correlation is also significant. As a consequence, the higher the achieved resolution (i.e. sampling in frequency), the higher the smoothness of the resulting spectrogram due to the increased cross correlations between wavelet amplitudes at nearby frequencies. The presence of noise-like components, e.g. due to observational noise, contaminates (’blurs’) the scalogram, by making more uncorrelated wavelet amplitudes of nearby times and introducing ridges across wide ranges of frequencies.

In contrast to the continuous, the discrete wavelet transform (the review of which is outside the scope of this article) leads to uncorrelated wavelet amplitudes if the discrete wavelet analysing functions are orthogonal. However the lack or redundancy is not an advantage for the signal-noise decomposition, since there is no profound criterion for attributing the wavelet components to signal or to noise. Some thresholding criterion can be imposed, as commonly used in signal compression applications, again having a degree of arbitrariness.

From the above discussion emerges a plan for an efficient signal-noise decomposition using the continuous wavelet transform, consisted of the following steps:

(a) Continuous wavelet transform of the timeseries with sufficiently high frequency resolution (sampling), thus producing a highly redundant wavelet scalogram.

(b) A non-parametric, information dependent selection criterion that disposes effectively the redundancy, in a way that selects separately those wavelet components that can be with high probability attributed to signal and to those to noise.

### 3.2 Scale maximal wavelet skeleton spectrum

A transient structure or burst of finite, characteristic duration $`\tau `$, occurring (i.e. having center of mass) in the timeseries at time $`t`$ is expected to cause a local increase of the amplitudes $`|a(t,f)|`$ of the timeseries wavelet transform in the vicinity of time $`t`$, at wavelet period of the order of $`\tau `$ (or else frequency $`1/\tau `$), while can extent across different timescales, e.g. a discontinuity ($`\tau =0`$) is wavelet-transformed to a superposition of wavelets with the same center of mass at time $`t`$, at all timescales (and frequencies).

For the study of such structures the scale maximal wavelet skeleton spectrum is extracted from the overall wavelet spectrum, keeping only those wavelet components of which the complex amplitude is locally maximum across time at any given time-scale, i.e. those for which:

$$\frac{|a(t,f)|}{t}=0,\frac{^2|a(t,f)|}{t^2}<0$$
(15)
It is already recognized that the so defined spectrum (otherwise called wavelet transform modulus maxima) contains all the important information about the existence of hierarchies of discontinuities, usually appearing as continuous lines across scales pointing at the time of the discontinuity, hence can be used for detection of singularities (e.g. Mallat, 1999, p. 176), or even approximate reconstruction of timeseries (e.g. Mallat, 1999, p. 185 and references therein). It is also used for the study scaling (e.g. possible multifractal) properties in timeseries: the lines of each singularity is often, at least at high frequencies, consisted of wavelet amplitudes that follow power-law dependence on frequency across many orders of frequency range (i.e. for these ranges there is no characteristic scale, or else there is some scale-invariance in the timeseries), well modeled by multifractals: If the timeseries is a realization of a multifractal then the wavelet amplitudes of the locally scale-invariant skeleton spectrum, at least at high frequencies, $`f+\mathrm{}`$, follow power-law of the form (for a review see e.g. Arneodo et al., 1995):

$$|a(t,f)|f^{h(t)}$$
(16)
where $`h(t)`$ is the H$`\ddot{o}`$lder (or else local Hurst) exponent, being a local measure of ’burstiness’ in the timeseries at time $`t`$ (for details see e.g. Muzy et al., 1993).

Power-law spectral power dependance over frequency, at least over several orders of frequency magnitude and especially tails at high frequencies (or small scales) are fairly common in all physical observables. In the special case of uni-fractals $`h(t)`$ has the same value for all $`t`$. Of specific interest is the case of $`1/f`$ noise (i.e. with $`h=0`$), constituting the hallmark of self-organized criticality, (SOC, Bak et al., 1987) i.e. the meta-stable state close to a critical point of a universal class of systems characterized by intermittency, well modeling aspects of solar flares, earthquakes, geological formation, biological evolution, economy etc. (for a review see Bak, 1996).

### 3.3 Instantly maximal wavelet skeleton spectrum

A transient or stable periodic component of instant period $`T`$ at time $`t`$ in the timeseries is expected to cause increase of the wavelet amplitudes $`|a(t,f)|`$ in the vicinity of wavelet period of the order of $`T`$ (or else frequency $`1/T`$), while may extend across time as far as its period $`T`$ remains constant in time.

For the detection of periodicities in timeseries here we introduce the instantly maximal wavelet skeleton spectrum, extracted from the overall wavelet spectrum, if keeping only those wavelet components that are locally of maximum amplitude at a given time, i.e. only those for which:

$$\frac{|a(t,f)|}{f}=0,\frac{^2|a(t,f)|}{f^2}<0$$
(17)
From the above discussion becomes clear that the so defined spectrum can be used for detecting transient periodicities of varying period, and of course those of stable period in the timeseries. Note that the instantly maximal skeleton spectrum excludes the wavelet components that instantly extend across time scales, thus it is less contaminated by the noise-like scale-invariant component of the timeseries. On the other hand, the scale maximal spectrum includes only those wavelet components that extend across time scales, thus the scaling properties derived from it are not distorted by possible transient periodicities in the timeseries due to the presence of components of analog signal type. In that sense the two skeleton spectra are complementary, permitting the separate study of the noise-like and signal-like timeseries component properties separately. Noise-like components of an hierarchy of discontinuities also affects strongly the complex phases, $`arg(a(t,f))`$, of the wavelet spectrum since the wavelet expansion requires ’synchronization’ of wavelets (i.e. to be placed nearby in time). For this reason only the wavelet amplitudes are considered in $`(\text{17})`$. However, auto-correlated random noise has significant slowly varying components, i.e. non-periodic trends at various scales (and in that sense the two skeleton spectra do not completely separate signal from noise), nevertheless introducing a power-law background without concentration of power in spectral peaks. In all, the instantly maximal skeleton spectrum is expected to exhibit more clearly the periodicities of a signal-like component of timeseries, practically excluding uncorrelated hierarchies of discontinuities due to noise.

## 4 Application to the sunspot index

### 4.1 A review on solar activity and sunspots

The solar activity (e.g. sunspots, flares, coronal mass ejections etc) is rich, involving complex, not well understood magnetofluid processes, extending in a wide range of spatial scales and durations. The sunspot number is widely regarded as an easily collectable quantitive solar activity measure. Sporadic, naked-eye observations exist in Chinese dynastic histories since 28 BC, while systematic sunspot observations are kept since the discovery of telescope in $`1610`$, while not uniformly reliable (Waldmeier 1961, Eddy, 1977, Sonett 1983).

The sunspot number is dominated by the 11-year periodicity (well-known Schwabe in 1843) connected to the 22-year cycle of the overall solar magnetic field. The 11-year cycle is stable, exhibiting only small variations of each cycle duration, as also shown using wavelet analysis of sunspot number (Fligge et al., 1999). However, the maximum number of sunspots in each cycle varies significantly and quite unpredictably, even nearly, but not totally, vanishing for many cycles. In all, the sunspot record includes enigmatic epochs of suppression, as the Maunder Minimum from 1645 to 1715 (e.g. Ribes & Nesme-Ribes, 1993) and Sp$`\ddot{o}`$rer Minimum in the $`15^{th}`$ century, as well as epochs of enhanced activity, as the nowadays Modern Maximum and Medieval Maximum in the $`12^{th}`$ century. The reason for the large variations of the cycle’s amplitudes and the epochs of suppressed activity of Maunder minimum type is unclear, e.g. proposed being due to chaotic behaviour of the nonlinear dynamo equations (e.g. Ruzmaikin, 1983, K$`\ddot{u}`$ker et al., 1999) or due to stochastic instabilities forcing the solar dynamo leading to on-off intermittency (e.g. Hoyng, 1988, 1993, Ossendrijver & Hoyng, 1996, Schmitt et al., 1996).

The average shape of each sunspot cycle is systematically asymmetric, taking less time rising to maximum than reaching the next minimum, implying that the solar cycle is intrinsically non-linear (e.g. Veselovski & Tarsina, 2002). Probably the most significant relation between the cycle shape characteristics is the Waldmeier effect: cycles with larger amplitudes are more asymmetric, taking less time to reach maximum (e.g. Hathaway, et al., 1994, Li, 1999, Veselovski & Tarsina, 2002). That the 11-year oscillation is not harmonic has an important implication: it introduces infinite number of spectral peaks (harmonics) with periods $`11/n`$ years, $`n=2,3,\mathrm{}`$ (e.g. Polygiannakis et al., 1996, Mursula et al., 1997). All those general features were successfully described by a simplified mono-parametric Van der Pol non-linear RLC oscillator model, further shown to accurately reproduce the observed sunspot record and extended epochs of suppressed activity, given the maximum of each cycle (Polygiannakis et al., 1996).

Systematic search for periodicities other than the basic 11-year in sunspot number, solar activity indices and, consequently, in the interplanetary, geomagnetic, earth rotation fluctuations and climate parameters lead to numerous claims, also depending on the nature of the analyzed parameter (as regarding its sensitivity to the solar activity periodicities), the periodicity detection method and the finite size of the records:

(a) Periodicities longer than 11 years, such as 180 years, 90 years (’Gleiseberg cycle’), 45 and 22 years (e.g. Cohen & Lintz, 1974, Sonnet 1982), perhaps subharmonics (i.e. integer multiples) of the basic 11-year periodicity, and several others (Norderman & Trivedi, 1992), intended to describe the long-term modulation of the sunspot cycle amplitudes as periodic, in contrast to the chaotic or stochastic theoretical scenarios discussed above.

(b) 11-year periodicity harmonics, namely of 5.5 (FYO, five-year oscillation) and 3.7 years in sunspot number, other solar activity and solar-terrestrial connection parameters as the geomagnetic field and earth rotation disturbances (e.g. Currie, 1976, Courtillot & Le Mouel, 1976, Suguira, 1980, Carta et al., 1982, Kondor, 1993, Djurovic & Paquet, 1996, Mursula et al., 1997, Nayar et al., 2002).

(c) Numerous shorter periodicities of transient character, detected over one or few 11-year cycles, with most frequently appearing in bibliography the 2-year (QBO, quasi-biennial) oscillation in sunspot number (Shapiro & Ward, 1962, Apostolov, 1985), solar radio flux (Hughes & Kesteven, 1981), solar magnetic helicity (Bao & Zhang, 1998), other solar activity, solar-terrestrial connection and atmospheric parameters and earth rotation rate (e.g. Chao, 1989, Djurovic & Paquet, 1993, Kane, 1997). Oftenly discussed are also periodicities of 154, 129, 103, 77, 51 days in solar flare occurrence and sunspot area (Riegel et al., 1984, Bogard & Bai, 1985, Dennis, 1985, Kile & Cliver, 1991, Bai & Surrock, 1991, Carbonel & Ballester, 1992).

The 11-year sunspot (equivalently, 22-year solar magnetic) cycle phenomenon, is now believed to be due to a magneto-hydrodynamic dynamo action, periodically regenerating the solar magnetic field at the basis of the convection zone. The cause and importance of periodicities other than the 11-year remain so far unclear. The periodicities of 154, 129, 103, 77, 51 days were proposed to be subharmonics of a 25.5 days fundamental periodicity of some ’clock mechanism’ in the sun, modeled by an oblique rotator, (Bai & Surrock, 1991, Sturrock & Bai, 1992). Alternatively it was proposed that the 154-day periodicity and others can be related to the normal modes of the solar interior oscillations (Wolff, 1983, 1992). The source of the QBO may be connected with periodic poleward streams of magnetic flux, first discussed by Howard and LaBonde (1981). Based on magnetic helicity observations (Bao & Zhang, 1998), it was proposed a double-cycle solar magnetic dynamo model, with one dynamo operating at the basis of the convection zone (11-year cycle) due to radial shear and the other at the top of the convection zone (with 2-year periodicity) due to latitudal shear (Benevolenskaya, 1998, 2000). Also supporting this model, further analysis of the solar magnetic field data showed that the active region magnetic field exhibits an intense, short-scale component varying with 2-year periodicity (Erofeev, 2001). It was also proposed that the sunspot number can be well fitted by a superposition of the usual 11-year cycles and wave trains with periodicity continuously varying from 38 months at solar maximum, to 21 months towards solar minimum (Krasotkin, 2001).

The sunspot number and all other solar activity parameters also exhibit non-periodic intense fluctuations that should be attributed to the episodic character of short-timescale solar phenomena. The solar magnetic field on the photosphere exhibits complex organization and sudden, intermittent energy emissions as the solar flares, in all well modeled as if the solar corona is in a self-organized critical state (SOC). The corresponding cellular automata models are in good agreement with the observed power-law distributions of flare released energy and occurrence rate (e.g. Dennis, 1985, Lu & Hamilton, 1991, Lu et al., 1993, Georgoulis & Vlahos, 1998, Isliker et al., 2000, 2001, Anastasiadis, 2002). The emergence of a flare drastically reorganizes the magnetic fields in solar active regions, hence is expected to affect the underlying sunspot and active region lifetimes, thus leaving an imprint of a, characteristic of SOC, $`1/f`$ tail at high frequencies in all the solar-activity related parameters, as the here considered sunspot index.

From the timeseries point of view the sunspot number record provides a good example of mixture of, both stable and transient, periodicities of significantly varying amplitudes and of high-frequency, noise-like components (e.g. Watari, 1996). From the physical point of view it promises to provide useful information about the, not completely understood, nature of the solar activity, as connected to the overall solar magnetic field regeneration mechanism and, at smaller timescales, the generation of active, transient solar phenomena, such as flares and coronal mass ejections.

### 4.2 Wavelet spectra of sunspot index

For this application we the record of 2048 recent monthly sunspot number timeseries extending from June 1831 to February 2002, compiled from the SIDC RWC Belgium World Data Center for the Sunspot Index catalogues. The number of points was selected being a power of 2 in order to avoid edge effects in the fast Fourier transforms used to derive the wavelet scalograms, thus also omitting the record before 1830 which are of questionable reliability (e.g. Sonett 1983). The use of monthly-averaged numbers excludes effects related to the solar differential rotation, with period of 27 days at the solar equator, and the time-life of each sunspot, persisting for one or for a few solar rotations.

The wavelet transform was calculated using the Morlet wavelet (equ. $`(\text{13})`$), filtering the timeseries for each wavelet scale in the frequency domain for $`N=100`$ total wavelet periods, sampled logarithmically between the minimum sampling time (1 month) and the total timeseries duration (2048 months) and then inverse Fourier transforming it in the time domain (power theorem, see $`(\text{6})`$). Note that the total number of period samples, $`N`$, used to derive the wavelet scalogram is arbitrary, depending on the information content in the record (defining the complexity of the scalogram). Too large $`N`$ may lead to mode splitting effects and too small to poorly sampled scalogram. Generally, selecting $`N`$ being of the order of $`1/10`$ of the total number of points in the timeseries (here 2048) gives satisfactory results. For a review on continuous wavelet transform and the details of the related numerical techniques used to derive the wavelet scalogram, see Torrence & Compo (1998). Once the wavelet scalogram is calculated it is parsed twice, once across time at each scale and once across scales for each time instance, testing successive triads of amplitudes and keeping only the locally maximal, thus obtaining the two skeleton scalograms.

In figure 1(a) we present the monthly-averaged sunspot index. The well-known 11-year periodicity and the asymmetry of each cycle are prominent. Also, the fluctuations of the index are larger near each cycle’s maximum.

In figure 1(b) the corresponding wavelet scalogram is shown. Notice the stable 11-year periodicity, having nearly constant period with time. The spectrum is rich in higher-frequency components too, and the intensity of higher frequencies is larger at cycle’s maxima, since the fluctuations are more intense there. Note that there is some arbitrarily in the image of the scalogram (as this applies to the representation of any wavelet scalogram and more generally time-frequency distribution), in the sense that appropriate color or grey scales must be attributed to the range of wavelet amplitudes, resulting in visually richer or poorer images.

The two skeleton scalograms, extracted from the overall wavelet one, are shown in figures 1(c) and 1(d). In the locally scale-invariant spectrum it can be noticed that the lines of maxima proceed towards larger periods at the cycle’s maxima, since they represent the center of each cycle as an overall structure. The instantly maximal skeleton spectrum clearly shows the 11-year periodicity as a nearly horizontal line and others, also quite stable. Periodicities of about 22-years are also prominent, probably being sub-harmonics of the 11-year one. Near the period of 1 year the periodicities appear having significant evolution with time, lasting for many decades (several cycles). The periodicities of smaller periods seem becoming increasingly unstable and intermittent (i.e. the lines have more and more gaps), probably due to their small amplitude, and the presence of self-similar noise-like fluctuations.

The total wavelet power of the overall and the two skeleton spectra is shown in figure 2. The main spectral peaks are indicated with arrows, noting the corresponding period in years. On the first column of table 1 we give the period of each of these peaks. The corresponding period error estimates represent the width between successive wavelet periods, as sampled to derive the wavelet transform (equ. $`(\text{13})`$).

The overall wavelet spectrum is poor, exhibiting only few spectral peaks, namely at the 11, 5.5, 2 and 1 year, and an $`1/f`$ tail at periods smaller than 1-year. However, the instantly maximal skeleton spectrum exhibits several spectral peaks, and this verifies its value: apart from the dominant 11-year periodicity, 3 of its harmonics, expected due to the 11-year cycle rise-fall asymmetry, with periods $`11/n,n=2,3,4`$ (i.e. $`5.5,3.67`$ and $`2.75`$ years) and rapidly decreasing amplitudes at shorter periods are also prominent. The weaker peaks of $`16,13`$ and $`8.8`$ years seem related with the main 11-year periodicity, probably being a mode-splitting effect. The long, $`28,24`$ and $`19`$ year periodicities are probably due to mode-splitting of a 22-year periodicity, being the main 11-year sub-harmonic, or can be interpreted as modulation effects of the cycle’s amplitudes (see the discussion for periods longer than 11-years in the previous section).

At shorter periods, a set of 10 spectral peaks with less rapidly decreasing amplitudes (nearly as $`1/f`$ with frequency) exists, with periods given in Table 1. Note that the QBO and the periodicities of 154, 129, 103, 77 days, i.e. 0.42, 0.35, 0.28, 0.21 years (previously reported in various solar activity parameters, as discussed in the previous section), are recovered.

As discussed in the previous section, according to the ’clock mechanism’ scenario (Bai & Surrock, 1991, Sturrock & Bai, 1992) the 154, 129, 103, 77 day periodicities constitute the $`6^{th},5^{th},4^{th}`$ and $`3^{th}`$ subharmonics of a 25.5 days fundamental periodicity. However the remaining 6 peaks are not easily explained in a similar way as even higher subharmonics. The QBO should be the $`29^{th}`$ subharmonic of the 25.5 day periodicity, however the subharmonics in between are not observed in the spectrum.

Alternatively, here we propose that all these periodicities are actually harmonics of the QBO periodicity, which then itself should be non-harmonic. In table $`1`$ we present the periodicities found and also those predicted by the hypothesis of only two non-harmonic fundamental oscillations, of 11 and 2 years. Note that $`12`$ out of $`15`$ detected spectral peaks, with periods of 11-year and shorter can be explained solely with this hypothesis. However:

(a) Two of the predicted QBO harmonics, of 0.29 and 0.22 years are not prominent in the spectrum, possibly being merged, due to noise, with the expected 0.33 and 0.2 year peaks respectively.

(b) The expected single spectral peak at 1-year period is wide, extending from about 0.9 to 1.6 years, within which 3 merged peaks of $`0.97,1.07,1.4`$ period are detectable. Possibly this peak widening is due to evolution of ephemeral, evolving periodicities (as seen in fig. 1(d)), indicative of transient underlying solar activity, or, alternatively, the $`0.97`$ and $`1.4`$ year peaks are additional periodicities (independent of the QBO) of unknown origin.

The scale maximal skeleton spectrum is affected only by the main 11-year periodicity, appearing as a dominant spectral peak, while at timescales of 1-3 years, where the main QBO peak resides, the scale maximal spectral power is reduced. At periods smaller than 1 year the scale maximal spectrum exhibits a plateau of nearly zero-slope. This power law-behavior at high frequencies is expected from $`(\text{16})`$ for Hurst exponent $`h=0`$. The noise-like component in the sunspot index thus exhibits power-law (scale-invariant) properties characterizing the $`1/f`$ noise, supporting the hypothesis of self-organized criticality (SOC) in the solar activity, previously proposed for the solar flare occurrence, as discussed in the previous section.

The observation that all the spectral peaks at periods of 2 years and shorter, here proposed to be harmonics of the QBO, have amplitudes decreasing nearly as $`1/f`$ with frequency (clearly differently from the 11-year cycle, the harmonics of which decay much faster), also emerging in the same frequency range as the $`1/f`$ noise-like component (apart from the 2-year peak itself) indicates a connection of these oscillations to surface (photospheric) solar activity phenomena, such as the solar flares. From theoretical point of view this observation again supports the model of double solar dynamo (Benevolenskaya, 1998, 2000) in which the 2-year oscillation is due to latitudal shear at the top of the convection zone, hence also the photosphere. The nature and properties of the surface dynamo producing both the QBO and intermittent bursts of $`1/f`$-type (SOC), and the coupling with the main 11-year dynamo at the basis of the convection zone is not within the scope of the present article, yet are to our opinion of great interest for solar physics and deserve further investigation.

## 5 Conclusions

We proposed the parallel use of two complementary skeleton spectra, extracted from the overall continuous wavelet spectrum, that are useful for studying the existence and evolution of stable or transient periodicities and noise-like components in timeseries. As an indicative physical application, we analyzed the monthly averaged sunspot index. The results showed:

(a) Clearer spectral peaks in the instantly maximal, compared to the overall, wavelet spectrum, that can be further attributed to the action of only two non-harmonic oscillators, namely the well-known 11-year and the 2-year (QBO), supporting the double-cycle solar dynamo model (Benevolenskaya, 1998, 2000).

(b) The existence of a $`1/f`$ noise background at timescales of 1 year and shorter, possibly being an imprint of the previously proposed self-organized critical state (SOC) of the solar magnetic fields and the intermittent energy releases, such as solar flares (e.g. Lu & Hamilton, 1991, Lu et al., 1993).

We hope that the two wavelet skeleton spectra will be further useful in various signal processing applications, as data compression, filtering, data fitting and modeling.

## Acknowledgments

This research of the author J.Polygiannakis is funded by the Greek National Foundation for Scholarships (IKY). The author J.Polygiannakis wishes to thank Dr. V.Agapitou of Queen Mary College of London, School of Mathematical Sciences, for helpful discussions. We express our thanks to the University of Athens Research Foundation. We are also grateful to RWC Belgium World Data Center for making the Sunspot Index catalogues available to the public.

## References

Anastasiadis, A.: 2002, J. Atm. Solar-Terr. Phys., 64, 481-488.

Apostolov, E.M.: 1985, Bull. Astron. Inst. Czechosl., 36, 97-102.

Argoul, F., Arneodo, A., Grasseau, G., Gagne, Y., Hopfinger, E.J., Frisch, U.: 1989, Nature, 338, 51.

Arneodo, A., Grasseau, G., Holschneider, M.: 1988, Phys. Rev. Lett., 20, 2281.

Arneodo, A., Bacry, E., Muzy, J.F.: 1995, Physica A, 213, 231-75.

Arneodo, A., d’Aubenton-Carafa, Y., Bacry, E., Graves, P.V., Muzy, J.F., Thermes, C.: 1996, Physica D, 96, 291.

Arneodo, A., Manneville, S., Muzy, J.F.: 1998, Eur. Phys. J. B, 1, 129.

Aschwanden, M.J., Kliem, B., Schwarz, U., Kurths, J., Dennis, B.R., Schwarz, R.A.: 1998, Astroph. J., 505, 941-956.

Bak, P., Tang, C., Wiesenfeld, K.: 1987, Phys. Rev. Lett., 59, 381-384.

Bak, P.: 1996, How nature works, the science of self organized criticality, Copernicus, Springer-Vergag, New York.

Bai, T., Surrock, P.A.: 1991, Nature, 350, 141.

Bao, S., Zhang, H.: 1998, Astroph. J., 496, L43.

Ballester, J.L., Oliver, R.: 1999, Astroph. J., 522, L153-L156.

Banerjee, D., O’Shea, E., Doyle, J.G.: 2000, Solar Phys., 196, 63.

Baumjohann, W., Treumann, R.A., Georgescu, E., Haerendel, G., Fornacon, K.-H., Auster, U.: 1999, Ann. Geophys., 17, 1528.

Bendjoya, Ph., Slezak, E., Froeschle: 1991, Astron. Astrophys., 251, 312.

Benevolenskaya, E.: 1998, Astrophys. J., 509, 49.

Benevolenskaya, E.: 2000, Solar Phys., 191, 247-255.

Benzi, R., Biferale, L., Crisanti, A., Paladin, G., Vergassola, M., Vulpiani, A.: 1993, Physica D, 65, 352.

Beylkin, G., Coifman, R., Rokhlin, V.: 1991, Comm. on Pure and Appl. Math., 44, 141.

Bijaoui, A., Slezak, E., Mars, G.: 1993, Wavelets, Fractals and Fourier Transforms, Oxford University Press, Oxford., p. 213.

Bocchialini, K., Baudin, F.: 1995, Astron. Astroph. 299, 893.

Bogard, R.S., Bai, T. : 1985, Astroph. J., 299, L51.

Carbonel, M., Ballester, J.L.: 1992, Astron. Astroph., 255, 350-362.

Carta, F., Chlistovski, F., Mannara, A., Mazzolenni, F: 1982, Astron. Astroph., 114, 388.

Castaing, B., Gagne, Y., Hopfinger, E.J.: 1990, Physica D, 46, 117.

Chakrabarti, K., Garofalakis, M., Rastogi, R., Shim, K: 2001, VLDB Journ., 10, 199-223.

Chao, B.F.: 1989, Science, 243, 923-925.

Cohen, T.J., Lintz, P.R.: 1974, Nature, 250, 398.

Cohen, A., Kovacevic, J.: 1996, Proc. IEEE, 84, 4, 514.

Courtillot, V., Le Mouel, J.-L.: 1976, J. Geophys. Res., 81 (17), 2941.

Currie, R.G.: 1976, Astroph. Space Sci., 39, 251.

De Moortel, Ireland, J., Walsh, R.W.: 2000, Astron. Astroph., 355, L23.

De Moortel, Hood, A.W., Ireland, J. : 2002a, Astron. Astroph., 381, 311-323.

De Moortel, Ireland, J., Hood, A.W., Walsh, R.W. : 2002b, Astron. Astroph., 387, L13-L16.

Dennis, B.R.: 1985, Solar Phys., 100, 65.

Djurovic, D., Paquet, P.: 1993, Astron. Astroph., 277, 669-676.

Djurovic, D., Paquet, P.: 1996, Solar Phys., 167, 427.

Eddy, J.A. : 1977, The Ancient Sun, ed. Pepin,R.A., Eddy, J.A., Merrill, R.B., Pergamon, New York.

Erofeev, D.V.: 2001, Solar Phys., 198, 31-50.

Fligge, M., Solanki, S.K.: 1997, Astron. Astroph. Supl., 124, 579.

Fligge, M., Solanki, S.K., Beer, J.: 1999, Astron. Astroph., 346, 313.

Fludra, A. : 2001, Astron. Astroph., 368, 639-651.

Foster, G.: 1996, Astron. J., 112 (4), 1709.

Frick, P., Galyagin, D., Hoyt, D.V., Nesme-Ribes, E., Schatten, K.H., Sokoloff, D., Zakharov, V.: 1997, Astron. Astroph., 328, 670.

Frick, P., Baliunas, S.L., Galyagin, D., Sokoloff, D., Soon, W.: 1997, Astroph. J. 483, 426.

Gabor, D.: 1946, J. Inst. Electr. Engin., 93, 429.

Georgoulis, M.K., Vlahos, L: 1998, Astron. Astroph., 336, 721-735.

Girardi, M., Fadda, D., Escalera, E.; Giuricin, G., Mardirossian, F.; Mezzetti, M.: 1997, Astroph. J., 490, 56.

Gimenez de Castro, C.G., Raulin, J.P., Mandrini, C.H., Kaufmann, P., Magun, A.: 2001, Astron. Astroph., 366, 317-325.

Hathaway, D.H., Wislon, R.M., Reichmann, E.J.: 1994, Solar Phys., 151, 177.

Haynes, P.H., Norton, W.A.: 1993, Wavelets, Fractals and Fourier Transforms, Oxford University Press, Oxford., p. 229.

Hempelmann, A., Donahue, R.A.: 1997, Astron. Astroph., 322, 835-840.

Hempelmann, A. : 2002, Astron. Astroph., 388, 540-545.

Howard, R., LaBonde, B.J.: 1981, Solar Phys., 74, 131.

Hoyng, P.: 1988, Astroph. J., 332, 857.

Hoyng, P.: 1993, Astron. Astroph., 272, 321-339.

Hughes, V.A., Kesteven, M.J.L.: 1981, Solar Phys., 71, 259.

Norderman, D.J.R., Trivedi, N.B.: 1992, Solar Phys., 142, 411-414.

Hunt, J.C.R., Kevlahan, N.K.R., Vassilicos, J.C., Farge, M.: 1993, Wavelets, Fractals and Fourier Transforms, Oxford University Press, Oxford.

Isliker, H., Anastasiadis, A., Vlahos, L.: 2000, Astron. Astroph., 363, 1134-1144.

Isliker, H., Anastasiadis, A., Vlahos, L.: 2001, Astron. Astroph., 377, 1068-1080.

Ireland, J., Walsh, R.W., Harrison, R.A., Priest, E.R.: 1999, Astron. Astroph., 347, 355-365.

Kaiser, G.: 1994, A friendly guide to wavelets, Birkhauser, Boston.

Kane, R.P.: 1997, Ann. Geophys., 15, 1581-1594.

Kile, J.N, Cliver, E.W.: 1991, Asroph. J, 370, 442-448.

K$`\ddot{u}`$ker, M., Arlt, R., R$`\ddot{u}`$diger, G.: 1999, Astron. Astroph., 343, 977-982.

Komm, R.W., Gu, Y., Hill, F., Stark, P.B., Fodor, I.K.: 1999, Astroph. J., 519, 407.

Kondor, N.N.: 1993, Adv. Space Res., 13 (9), 417.

Krasotkin, S.A.: 2001, SOLSPA 2001 Euroconf.: Solar and Space weather, Proc. ESA SP-477.

Lawrence, J.K., Cadavid, A.C., Ruzmaikin, A.: 2001, Solar Phys., 202, 27-39.

Li, K.: 1999, Astron. Astroph., 345, 1006-1010.

Lima Neto, G. B., Pislar, V., Durret, F., Gerbal, D.; Slezak, E.: 1997, Astron. Astroph., 327, 81-89.

Louis, A.K., Maa$`\beta `$, P., Rieder, A.: 1997, Wavelets, Theory and Applicatons, Wiley & Sons, Chichester, UK., p. 16.

Lu, E.T., Hamilton, J.R.: 1991, Astroph. J., 380, L89.

Lu, E.T., Hamilton, J.R., Mc Tiernan, J.M., Bromund, K.R.: 1993, Astroph. J., 412, 841.

Lucek, E.A., Balogh, A.: 1997, Geophys. Res. Lett., 24, 18, 2387.

Mallat, S. G.: 1989, IEEE Trans. on pattern analysis and machine intelligence, 11, 7, 674.

Mallat, S.G.: 1999, A wavelet tour of signal processing, Second Edition, Academic Press, CA.

Marsh, M. S., Walsh, R. W., Bromage, B. J. I.: 2002, Astron. Astroph., 393, 649-659.

Meszarosova, H., Jiricka, K., Karlicky, M.: 1999, Astron. Astroph., 348, 1005.

Muhlmann, W., Hanslmeier, A.: 1996, Solar Phys., 166, 445.

Monin, A.S, Yaglom, A.M.: 1975, Statistical fluid mechanics: mechanics of turbulence, vol. 2, Lumley, J.L. (ed.), MIT Press, Cambridge, Mass.

Mursula, K.; Usoskin, I.; Zieger, B., 1997, Solar Phys., 176, 1, 201.

Mordvinov, A.V., Kuklin, G.V.: 1999, Solar Phys., 187, 223-226.

Mordvinov, A.V., Plyusnina, L.A.: 2000, Solar Phys., 197, 1.

Morlet, J.G., Arehs, I, Fourgeau, Giard, D.: 1982, Geophysics, 47, 203.

Mouri, H., Takaoka, M., Kubotani, H.: 1999, Phys. Lett. A, 461, 82.

Muzy, J.F., Bacry, E., Arneodo, A.: 1991, Phys. Rev. Lett., 67, 25, 3515.

Muzy, J.F., Bacry, E., Arneodo, A.: 1993, Phys. Rev. E, 47, 2, 875.

Nayar Prabhakaran, S.R., Radhika, V.N., Revathy, K., Ramadas, V.: 2002, Solar Phys., 208, 359-373.

O’Shea, E., Banerjee, D., Doyle, J.G., Fleck, B., Murtagh, F.: 2001, Astron. Astroph., 368, 1065.

O’Shea, E., Murtagh, F., Fleck, B.: 2002, Astron. Astroph., 387, 642-664.

Ossendrijver, A.J.H., Hoyng, P.: 1996, Astron. Astroph., 313, 959-970.

Otazu, X., Ribo, M., Peracaula, M., Paredes, J.M., Nunez, J.: 2002, Mon. Not. R. Astron. Soc., 333, 365-372.

Pantin, E., Starck, J.L.: 1996, Astron. Astroph. Suppl., 118, 575.

Polygiannakis, J.M., Moussas, X.: 1996, Solar Phys., 163, 193-203.

Press, W.H., Flannery, B.P., Teukolski, S.A., Vetterling, W.T.: 1992, Numerical Recipes, Cambridge University Press, Cambridge, p. 420.

Ribes, J.C., Nesme-Ribes, E.: 1993, Astron. Astroph., 276, 549-563.

Riegel, E., Share, G.H., Forrest, D.J., Kanbach, G., Reppin, C., Chupp, E.L.: 1984, Nature, 312, 623-625.

Rigozo, N.R., Echer, E., Vieira, L.E.A., Nordemann, D.J.R. : 2001, Solar Phys., 203, 179-191.

Roux, S., Muzy, J.F., Arneodo, A.: 1999, Eur. Phys. J. B, 8, 301.

Rusmaikin, A.A.: 1983, Stellar and Planetaty magnetism, Soward, A.M., (ed.), Gordon & Breach, New York, p. 151.

Ruzmaikin, A.A., Lyannaya, I.P., Styashkin, V.A., Yeroshenko, E.: 1993, J. Geophys. Res., 98, 13303.

Saito, N., Beylkin, G.: 1993, IEEE Trans. on Signal Proc., 41, 12, 3585.

Schmitt, D., Schussler, M., Ferriz-Mas: 1996, Astron. Astroph., 311, L1-L4.

Schwarz, U., Kurths, J., Kliem, B., Kruger, A., Urpo, S.: 1998, Astron. Astroph. Suppl., 127, 309.

Sello, S.: 2000, Astron. Astroph., 363, 311-315.

Schneider, K., Kevlahan, N.K.-R., Farge, M.: 1997, Theor. Comp. Fluid Dyn., 9. 191-206.

Shapiro, R., Ward, F.: 1962, J. Atmos. Sci., 19, 127.

Sonett, C.P. : 1983, J. Geoph. Res., 88, 3225.

Soon, W., Frick, P., Baliunas, S.: 1999, Astroph. J., 510, L135-L138.

Sturrock, P., Bai, T.: 1992, Astroph. J., 397, 337.

Torrence, C., Compo, G.P.: 1998, Bull. Am. Met. Soc., 79, 61.

Veselovski, I.S., Tarsina, M.V.: 2002, Adv. Space Res., 29, 3, 417-420.

Waldmeier, M. : 1961, The sunspot activity in the years 1610-1960, Zurich Schulthess and Company AG.

Watari, S.: 1996, Solar Phys., 168, 413-422.

Wolff, C.L.: 1983, Astroph. J., 264, 667.

Wolff, C.L.: 1992, Solar Phys., 142, 187.

Zhang, Y., Paulson, K.V.: 1997, Pure Appl. Geophys., 149, 405.

## Figure captions

Figure 1. The monthly-averaged sunspot index. (b) The corresponding wavelet scalogram in grey color scale (whiter color corresponds to larger wavelet amplitudes). (c) The scale maximal and (d) the instantly maximal wavelet skeleton spectra.

figure 2. The total wavelet power of the overall (upper curve), the instantly maximal (middle curve, the period in years of the spectral peaks is also noted) and scale maximal wavelet skeleton spectra in common, arbitrary units.