Wavelet analysis: a new significance test for signals dominated by intrinsic rednoise variability
Abstract
We develop a new statistical test for the wavelet power spectrum. We design it with purpose of testing signals which intrinsic variability displays in a Fourier domain a rednoise component described by a single, broken or doublybroken powerlaw model. We formulate our methodology as straightforwardly applicable to astronomical Xray light curves and aimed at judging the significance level for detected quasiperiodic oscillations (QPOs). Our test is based on a comparison of wavelet coefficients derived for the source signal with these obtained from the averaged wavelet decomposition of simulated signal which preserves the same broadband model of variability as displayed by Xray source. We perform a test for statistically significant QPO detection in XTE J1550–564 microquasar and active galaxy of RE J1034+396 confirming these results in the wavelet domain with our method. In addition, we argue on the usefulness of our new algorithm for general class of signals displaying type variability.
I Introduction
Fourier transform has become a very popular tool of timeseries analysis within last fifty years in almost all areas of research. Applied to test signal frequency content, it was quickly recognized as a helpful method for investigation of periodic variability [1]. It gained a particular interest in a wide field of astronomical research regarding different classes of objects with variable emission. In particular, timeseries provided by recent highenergy detectors onboard Xray orbiting NASA/ESA satellites [2]–[3] have been found to trace remarkably rapid variability coming from neutron star and blackhole systems [4]–[6]. Since the majority of Xray emission takes its origin in accretion process of matter and gas onto central object [7]–[8], therefore the observed signal variability is a direct manifestation of physical processes taking place in a very strong gravitational field [9].
The accretion onto blackholes is one of the most energetic process in the Universe [10]. Energy spectra for the corresponding objects uncover both thermal and nonthermal nature of emission (e.g. [11]). The analysis of Fourier power spectra (also referred to as Power Spectral Density; PSD) of Xray light curves has revealed that the overall shapes are very often well described, in general, by a powerlaw model, i.e. where denotes a slope ([5], [12]–[16]). In case of one claims about signals dominated by intrinsic variability of rednoise type. The rednoise character of source variability is generally rarely met in nature [17]. Nevertheless, it seems that blackhole systems emitting in Xrays enlarge the sample significantly.
To date the majority of efforts devoted to our curiosity of understanding the Xray variability of astronomical objects were based on the systematic study of PSD shapes and looking for occasionally observed quasiperiodic oscillations (QPOs; [5], [18]–[19]). Since Fourier PSD contains no information about time evolution of detected periodicities, the need of use of the timefrequency techniques emerged. The application of ShortTerm Fourier Transform (STFT; e.g. [20]) provided a good solution in this domain and has been found helpful in localizing QPOs (e.g. [21]–[22]). However, regardless the frequency of dectected QPOs, STFT keeps the same timefrequency resolution as dictated by a length of a sliding window. In most cases it prevents a detection of oscillations lasting shorter than the window span. On the contrary, a wavelet analysis occurred to be more successful where the timefrequency resolution was scale (frequency) dependant. Within recent fifteen years it has attracted attention of many researchers [23]–[24].
The wavelet power spectrum (scalogram) for discrete evenly sampled Xray signal (; ) given as countrate [cts s] can be defined as the normalized square of the modulus of the wavelet transform:
(1) 
where denotes a normalization factor and a discrete form of the continuous wavelet transform in a function of two parameters: scale and localized time index [24]–[26]:
(2) 
where the discrete Fourier transform of signal is given by:
(3) 
and denotes a frequency index,
(4) 
and a mean value of [26]. Following [27] one can assume normalization factor of to be which provides that integrated wavelet spectrum over time, also referred to as a global wavelet power spectrum [26],
(5) 
will hold the units of (rms/mean) Hz, i.e. it gives the variance, relative to the mean, within a given frequency range of integration [27]–[28]. For the computation of (1) the use of a dyadic grid of scales, i.e.
(6) 
can be implemented [26], where the smallest scale corresponds to the reversed Nyquist frequency, and
(7) 
A grid resolution is given by . Since the search for periodicity and quasiperiodicity in Xray light curves seems to be the main target, it is reasonable to assume the analyzing function to be Morlet wavelet,
(8) 
which oscillates due to a term , where parameter is set to 6.
In order to popularize wavelet analysis as a new tool in research over timeseries, an additional requirement in the form of significance test was strongly desired. Torrence & Compo [26] (TC98) were first who addressed that issue extensively providing the public with an excellent guide and software readytouse. In general, a map of wavelet power spectrum shows a distribution of spectral components in the time–scale (or time–frequency) plane. Single oscillation, of sufficiently high amplitude to be detected, marks itself in the map as a wavelet peak. To every peak a certain statistical significance can be assessed. We say that a wavelet peak is significant at assumed per cent of confidence when it is above a certain background spectrum. The latter can be defined by the mean Fourier power spectrum of analyzed timeseries.
TC98 in their work adopted the background power spectrum, , after [29] (see also Eq.(16) in [26]) which is valid for application as long as the analyzed signal is considered to be a realization of univariate lag1 autoregressive (AR1) process [30]–[32], , where is lag1 autocorrelation, denotes a random variable drawn from the Gaussian distribution (with zero mean and variance ) and . However, the PSD shape provided by [29] is not a good representation of broadband PSDs in case of considered here Xray sources. This is because of at least two reasons: (i) not every Xray light curve can be described by AR1 process [33] (although [34]), (ii) PSDs are dominated by the rednoise spectral components which shapes are mainly fitted with single, broken, doubly broken or even bending powerlaw models.
TC98 showed that for a timeseries modeled as the AR1 process its local wavelet spectrum follows the mean Fourier spectrum given by their Eq.(16) [26]. By the “local wavelet spectrum” they define a vertical slice (along the frequency axis) through a wavelet map. Empirical justification is done on the way of Monte Carlo simulations which were performed considering only one local wavelet power spectrum per simulation and, in addition, selected for a particular moment of time ( of 512 points). The corresponding distribution of the local wavelet spectrum (at each localized time index and scale ), ought to have distribution,
(9) 
where denotes the distribution with two degrees of freedom and an arrow “” means “is distributed as”. The above implication is performed based on the assumption that a random variable is normally distributed.
A noteworthy work in this domain has been recently done by Zhongfu Ge [35] who generalized TC98’s significance test. Ge left the assumption of signal to be the realization of AR1 process and proposed Gaussian White Noise process to serve as the null hypothesis.
What TC98 and Ge did not show is that the local wavelet power spectra, , at nearby times are correlated due to variable length of the wavelet function, and the correlation is stronger at low frequencies as shown by [36]. This introduces departures off distribution in any finite data set. Thus a more correct procedure is required to take this effect into account.
Maraun, Kurths & Holschneider [37] proposed another test for wavelet significance, namely, areawise test which aims at comparing the wavelet peak size to the expected timefrequency uncertainty. They noticed that for any process with unknown distribution a number of wavelet peaks forming in the wavelet map different shapes (circles, ovals, patches, etc.) become difficult in interpretation. The extracted picture is messy. Some peaks are due to rednoise process of the source intrinsic variability and the others are due to whitenoise. As shown by them, the areawise test may reduce a number of significant peaks even by 90 per cent.
To judge the reality of each wavelet peak it seems best to set the significance level highly enough in order to reject spurious results. However, not always this is a right way, for instance, in case of Xray sources of rather low countrate (low S/N ratio) and strong mixing between red and whitenoise components in highfrequency region. On the other hand, this also does not mean that the analysis does not contain any scientific results. For example, the same situation took place for COBE detections of microwave background peaks [38]. Most of them were spurious but the fluctuations were statistically different from expectations of pure noise [39].
None of this method can be applied as a fully trustful test working well for astronomical Xray signals as long as their assumptions rely on the analyses that use impropriete background spectra. In this paper we address that issue and propose a new significance test.
Ii Wavelet significance for rednoise signals
of powerlaw form
We design a method for determination of significance of wavelet spectral features (in short, wavelet peaks) calculated for any Xray light curve which intrinsic variability is a realization of a sum of red and whitenoise (Poisson) processes. The main idea standing for that is to compare the distribution of wavelet power between real data and simulated signal. Here, by the simulated light curve (also referred to as TK light curve) we will understand a timeseries which displays the same profile of rednoise variability as that one found for Xray data. Before providing the reader with readytouse algorithm (Sect. IIC), some essential information on the background assumptions need to be clarified. We provide them within the following two subsections (Sect. IIA and IIB).
Iia Simulations of process
Let’s assume that describes a part of stochastic continuous process generated by the physical system of unknown nature which produces the variability. Based solely on a reproduction of the overall information and properties of the system is very difficult. Recorded data can facilitate a determination of a model which would reproduce the system best. The observed light curve of Xray emission is the realization of the physical process(es) for which calculated Fourier power spectrum (PSD) reveals the variability of powerlaw form (see references in Section I). Therefore, we are able to approximate one of the system’s parameters (i.e. its variability) by a powerlaw model which, as we believe, constitutes the best PSD description of the underlying process.
Aperiodic broadband variability as displayed in Xray light curves can be thought, in the first approximation, as the realization of a linear stochastic process, weakly nonstationary. Timmer & Konig [40] proposed an algorithm (hereafter also referred to as TK) for generating linear aperiodic signals, , that exhibit powerlaw spectrum,
(10) 
where the phase is randomized . Since the amplitudes are taken as , [40] noticed that in order to create powerlaw signal it is essential to randomize also the amplitude.
A light curve obtained in this manner displays shape, however as argued by [33], it does not preserve fundamental statistical properties of Xray timeseries for accreting blackhole systems. The latter have been found to be rather nonlinear showing rms–flux relation at all timescales [33], [41] (rms is defined as a square root of variance calculated for a signal). According to [33], formally nonlinear type of variability for any linear input light curve can be obtained by taking simple exponential transformation of (10), i.e. . To make the picture complete every simulated light curve should include white noise. Following [25] it can be expressed as:
(11) 
where denotes simulated signal with Poisson noise, is a light curve resolution (bin time) and poisdev represents a subroutine generating a random number with Poisson distribution [42].
IiB Selection of best model for PSD of Xray source
A description of Fourier PSD of Xray object under investigation can be obtained in the process of numerical fitting of powerlaw model. As reviewed by [5], a small but substantial subset of Xray blackhole (and neutron star) systems apart from exhibiting like broadband variability, also show quasiperiodic oscillations (QPOs). Hardly ever met as tiny spikes, they are much often resolved as broadened features characterized by the quality factor of where denotes the peak’s full width at half maximum (FWHM) and QPO central frequency. The broader peak the higher thus the lower . So far, the majority of QPO peaks have been found to be surprisingly well described by Lorentzians, e.g. [12], [16]. Lorentzians are defined as with a centroid frequency and, if detected, they constitute an intergal component of the overall PSD model. In practice, fitting of powerlaw model as a background model may or may not include fitting of component simultaneously. In the latter case, PSD data points describing Lorentzian can be simply omitted and these two aforementioned approaches rely on the individual preferences and fitting experience.
In this work our wish will be to generate simulated Xray light curve which will own the same powerlaw background model of variability as fitted to the object’s PSD. This step requires that the overall rms and the average count rate for simulated timeseries must match the corresponding same values for Xray signal. For Fourier PSD normalized to unit of (rms/mean) Hz [27], the rootmeansquare (rms) can be obtained by integration of underlying PSD in a given band of frequencies [28] as follows:
(12) 
where . We argue that rms values, both obtained for simulated signal and Xray light curve ought to be calculated in the same range of frequencies what can guarantee the best match. The integration (12) should omit PSD data points corresponding to QPO component(s), if present.
IiC A new algorithm for wavelet significance
We define a new algorithm leading to determination of significance matrix for derived wavelet power spectrum by the following steps:

Twodimensional matrix of wavelet power spectrum (1) for a timeseries is a set of wavelet coefficients where corresponds to (7) and denotes a number of data points. Let us denote the results of computation of (1) for Xray source as:
(13) Please notice that in practice instead of time–scale plane one often refers to wavelet spectrum displayed in time–frequency plane. If required, for Morlet wavelet, a transformation from scales to Fourier frequencies can be applied according to relation [26].

For a given discrete signal calculate its Fourier power spectrum and find the best model which describes it (Section IIB). If the model for considered data set is known in advance (e.g. from the existing literature), use the best fitted parameters as the reference ones.

Generate one long simulated light curve (with the same ), best at least two or three orders of magnitude longer than your signal under analysis. Use Timmer & König’s method which allows to obtain simulated signals from Fourier power spectra of assumed shape (Section IIA). Take care to keep the rms variability at the same level as in the data (Section IIB).

From the simulated TK light curve, select randomly signals of the same time duration as in original data set, i.e. , and compute for each of them the wavelet power spectrum according to (1). Shall us denote the resulting matrices as:
(14) In calculation of (14) for normalization factor (see Sect. I) assume every time a mean count rate of the local TK signal .

Find the maximum value between all matrices:
(15) 
For every scale of every calculate the histogram of wavelet power distribution from 0 to with a histogram resolution defined as . In consequence, the following matrices shall be created:
(16) consisting of histograms.

Calculate a matrix containing averaged histograms for every scale :
(17) This step is introduced in order to smooth out different distribution of wavelet power appearing in every for each scale. Please note the higher the more smoothed averaged wavelet power histograms one obtains.

For every of histograms stored in , determine the th quantile . Namely, let be a continuous random variable, then for , the th quantile is the such that [44]. For instance, in deriving of the significance of wavelet peaks at the 95 per cent confidence level one aims in determination of 0.95th quantile in a function of scale (frequency). In practice, a numerical integration of normalized^{1}^{1}1A term normalized means that an integral of must equal 1. must be performed until an integral value exceeds 0.95. This determines quantile with an error of for each histogram:
(18) 
Expand matrix to have dimensions by calculating a linear algebraic product of the matrices and :
(19) where is a matrix of ones of dimension. Denote the result as a background wavelet matrix .

Finally, a matrix of significance, , will identify significant wavelet peaks appearing in computed wavelet power spectrum when the following condition will be met:
(20)
In proposed procedure a distribution of wavelet power at a given wavelet scale (frequency) is directly compared to the distribution of power obtained for simulated light curves. Since the Fourier PSD model and rms variability in simulated timeseries are assumed to match the same quantities as for the Xray source data, they serve as a good reference level to test any excess of wavelet power above it. In addition, by testing a wavelet power scale by scale, one takes into account an influence of correlation of wavelet power over a range of frequency and time due to varying length of a wavelet function (Sect. I, [36]).
Iii Application to Xray observations
In this section we provide the reader with two examples of application of our new algorithm to real observations of blackhole systems emitting in Xrays: microquasar XTE J1550–564 [45] and active galactic nuclei of RE J1034+396. In case of the former object, the source is known because of displaying a strong Xray aperiodic variability with a variety of lowfrequency (0.08–18 Hz) and highfrequency QPOs (100–285 Hz) [46]–[48]. Here, we examine its activity peaked at the frequency of 4 Hz and derive corresponding QPO significance in the wavelet domain. As for the latter source, we select the active galaxy of RE J1034+396 which has been recently found to uncover the very first ever significant signature (over ) of quasiperiodicity among all active galaxies [49]. Inspired by Fourier results of [49], we verify their findings by determination of wavelet significance for detected Hz variability.
Iiia Xte J1550–564
In order to examine XTE J1550–564 light curve variability which has been dominated by intrinsic rednoise component in – Hz frequency band, after [50] we select RXTE/PCA observation of 40099012401 (19980929) from the public HEASARC archive. We reduce the data with LHEASOFT package ver. 6.6.1 applying the standard PCA selection criteria for above data set (see http://heasarc.gsfc.nasa.gov for details). As the final product we extract a 2–13 keV light curve with a bin size of s for the first part of PCA observation (i.e. s).
Fig. 1 presents Fourier power spectrum density calculated from 78 averaged spectra based on 1024 point data segments. The overall PSD shape can be split into two main components: powerlaw underlying spectrum and broad line profile representing rednoise source variability and QPO component, respectively. In the fitting process we found that a model composed of threesegment broken powerlaw (i.e. with two break frequencies) and Lorentzian profile line [51] fit the PSD best yielding reduced of 10.9 at 59 degrees of freedom. The corresponding powerlaw slope values and break frequencies have been determined as follows: Hz and , whereas for Lorentzian the central frequency of Hz and its FWHM Hz were fitted providing QPO quality factor of . We marked both model components in Fig. 1 by red solid line and dotted line, respectively, where the overall best model has been indicated by black solid line.
Based on the above powerlaw model we generated a point long simulated signal (64 segments of 1024 points) with TK method. In simulation we assumed the sampling time of s and ensured that resulting rms value (here, derived in 0.5–2 Hz and 8–12 Hz frequency band; see Sect. IIB for details) was nearly the same as compared to corresponding rms values for Xray source. We calculated the averaged Fourier PSD for TK timeseries and plotted it in Fig. 1 by blue circle markers. One can clearly notice that TK PSD reproduces well the assumed double broken powerlaw model component.
Following the recipe given in Section IIC we aim at determination of wavelet confidence at 0.9995 level. From the entire XTE J1550–564 light curve, we select a 3072 point long test signal (96 s) for which we compute wavelet significance matrix (Eq. 20). Upper plot of Fig. 2 presents a dependence of 0.9995th wavelet quantile (Eq. 18) in a function of Fourier frequency in 1–10 Hz band, both for Xray data segment (black solid line) and TK simulation (red dashed line). It is straightforward to note that at the significance of the quasiperiodic variability around Hz is high above the mean (simulated) rednoise level. Finally, the bottom map of Fig. 2 presents the wavelet power spectrum (Eq. 1) for considered data segment with overlaid contours of derived significance.
IiiB Re J1034+396
A Seyfert 1 galaxy of RE J1034+396 has been found and confirmed as the first active galaxy showing highly significant quasiperiodic oscillation (QPO) at the frequency Hz [49]. We reproduce their PSD in Fig. 3 as computed for signal segment of XMMNewton EPIC observation of RE J1034+396 conducted on 20070531. Solid black line crossing PSD denotes bestfitted powerlaw spectrum described as:
(21) 
(M. Gierlinski, private communication). Upper two dashed lines mark computed confidence levels for detected periodicity at 99.73% () and 99.99% confidence levels. In order to verify aforementioned significance of QPO detection, we downloaded from XMMNewton Science Operations Centre the same EPIC pn+MOS1/MOS2 observation data set and reduced it in the same manner as described in [49]. We limit original Xray light curve (0.3–10 keV) to 601 point long timeseries (corresponding to segment 2 of [49]) where the sampling time of s has been used.
We compute wavelet power spectrum (1) in – Hz band and make use of our new algorithm in order to determine two contour plots marking 99.73% and 99.99% confidence for wavelet results. For the purpose of wavelet background simulation, we generate point TK signal based on assumed single powerlaw model (Eq. 21). Since the underlying PSD model is known in advance, we perform an integration of it between Hz and Hz to find the corresponding rms value in this frequency range (rms%). Within the simulation process of TK timeseries we check to make sure that in the same frequency band the rms of TK signal stays in agreement with derived value. Based on that, we determine two significance matrices (Eq. 20) and plot them together with wavelet power spectrum. The final form of all calculations presents Fig. 4 where an outter and inner contour denote 99.73% and 99.99% confidence level, respectively.
Remarkably, a detected quasiperiodic oscillation at Hz appears to be very coherent in time though showing small but noticeable period drift. We checked that the entire periodicity remains significant at the confidence level of 99% throughout the signal segment duration. In contrast, as marked in Fig. 4, only a small fraction of QPO activity exceeds top level of significance whereas about 70% of detected QPO signal stays significant at level.
Iv Conclusions
In this paper we reviewed a recent status of knowledge devoted to wavelet significance tests and proposed a new method designed especially for signals which intrinsic variability displays in Fourier power spectral density a powerlaw form. Since a large number of astronomical Xray sources has been found to be characterized by the similar PSD shapes, we focused our attention to that class of objects and their timeseries analysis.
We presented a flexible algorithm helpful in determination of the significance of quasiperiodic oscillations in the wavelet domain. Our method has been based on a comparison of wavelet coefficients calculated for Xray timeseries with those computed for simulated signal, where the latter displays the same PSD broadband model of variability as Xray object.
We argue that the usefulness of our new algorithm can be easily extended to nonastronomical timeseries which PSD can be described by powerlaw model as well. If required, the omission of exponential transformation of linear signal generated with the use of TK method and inclusion of Poisson noise should be avoided (Section IIA).
Acknowledgment
The author would like to thank Prof. Didier Barret for inspiration which led to the accomplishment of this paper and Dr. Marek Gierliński for a kind providing with an original version of Fig. 3.
References
 [1] R. N. Bracewell, ”The Fourier Transform and Its Applications”, McGrawHill, 1965, 2nd ed. 1978, revised 1986
 [2] J. H. Swank, “The XTE Mission”, Bulletin of the American Astronomical Society, vol. 26, pp. 1420, 1994
 [3] C. Gabriel, M. Guainazzi and L. Metcalfe, “XMMNewton: Passing Five Years of Successful Science Operations”, Astronomical Data Analysis Software and Systems XIV, vol. 347, pp. 425, 2005
 [4] P. Uttley, I. M. McHardy and I. E. Papadakis, “Measuring the broadband power spectra of active galactic nuclei with RXTE”, Monthly Notices of the Royal Astronomical Society, vol. 332, issue 1, pp. 231250, 2002
 [5] W. H. G. Lewin and M. van der Klis, “Compact stellar Xray sources”, Cambridge Astrophysics Series, No. 39. Cambridge, UK: Cambridge University Press, 2006
 [6] R. Sunyaev, and M. Revnivtsev, “Fourier power spectra at high frequencies: a way to distinguish a neutron star from a black hole”, Astronomy and Astrophysics, vol. 358, pp. 617623, 2000
 [7] N. I. Shakura and R.A. Sunyaev, “Black holes in binary systems. Observational appearance”, Astronomy and Astrophysics, vol. 24, pp. 337, 1973
 [8] C. M. Urry and P. Padovani, “Unified Schemes for RadioLoud Active Galactic Nuclei”, Publications of the Astronomical Society of the Pacific, vol. 107, p. 803
 [9] M. van der Klis, “Neutron Star QPOs as Probes of Strong Gravity and Dense Matter”, Xray Timing 2003: Rossie and Beyond. AIP Conference Proceedings, vol. 714, pp. 371378, 2004
 [10] J. Frank, A. King A. and D. J. Raine, “Accretion Power in Astrophysics”, Cambridge University Press, 2002
 [11] C. Done and M. Gierliński, “Observing the effects of the event horizon in black holes”, Monthly Notice of the Royal Astronomical Society, vol. 342, issue 4, pp. 10411055, 2003
 [12] K. Pottschmidt, J. Wilms, M. A. Nowak, G. G. Pooley, T. Gleissner, W. A. Heindl, D. M. Smith, R. Remillard and R. Staubert, “Long term variability of Cygnus X1. I. Xray spectraltemporal correlations in the hard state”, Astronomy and Astrophysics, vol. 407, pp. 10391058, 2003
 [13] A. Lawrence, M. G. Watson, K. A. Pounds,and M. Elvis, “Lowfrequency divergent Xray variability in the Seyfert galaxy NGC4051”, Nature, vol. 325, p. 694, 1987
 [14] I. M. McHardy, I. E. Papadakis, P. Uttley, M. J. Page and K. O. Mason, “Combined long and short timescale Xray variability of NGC 4051 with RXTE and XMMNewton”, Monthly Notices of the Royal Astronomical Society, vol. 348, issue 3, pp. 783801, 2004
 [15] T. Belloni and G. Hasinger, “Variability in the noise properties of Cygnus X1”, Astronomy and Astrophysics, vol. 227, no. 2, L33L36, 1990
 [16] M. A. Nowak, B. A. Vaughan, J. Wilms, J. B. Dove and M. C. Begelman, “Rossi XRay Timing Explorer Observation of Cygnus X1. II. Timing Analysis”, The Astrophysical Journal, vol. 510, issue 2, pp. 874891, 1999
 [17] W. Press, ”Flicker noises in astronomy and elsewhere”, Comments on Astrophysics, vol. 7, p. 103, 1978
 [18] M. van der Klis, “The QPO phenomenon”, Astronomische Nachrichten, vol. 326, issue 9, p.798803, 2005
 [19] M. van der Klis, “Comparing Black Hole and Neutron Star Variability”, Astrophysics and Space Science, vol. 300, issue 13, pp. 149157, 2005
 [20] L. Cohen, “Timefrequency analysis”, Englewood Cliffs, NJ: PrenticeHall, 1995
 [21] J. Wilms, M. A. Nowak, K. Pottschmidt, W. A. Heindl, J. B. Dove and M. C. Begelman, “Discovery of recurring softtohard state transitions in LMC X3”, Monthly Notices of the Royal Astronomical Society, vol. 320, issue 3, pp. 327340, 2001
 [22] D. Barret, W. Kluźniak, J. F. Olive, S. Paltani and G. K. Skinner, “On the high coherence of kHz quasiperiodic oscillations”, Monthly Notices of the Royal Astronomical Society, vol. 357, issue 4, pp. 12881294, 2005
 [23] Y. Meyer and S. Roques, “Progress in wavelet analysis and applications”, 1993, Proceedings of the International Conference ”Wavelets and Applications”, Toulouse, France, June 1992, GifsurYvette: Editions Frontieres, ed. by Meyer, Yves; Roques, Sylvie, 1993
 [24] P. S. Addison, “Illustrated Wavelet Transform Handbook”, IoP, 2002
 [25] P. Lachowicz and B. Czerny, “Wavelet analysis of millisecond variability of Cygnus X1 during its failed state transition”, Monthly Notices of the Royal Astronomical Society, vol. 361, issue 2, pp. 645658, 2005
 [26] C. Torrence and G. P. Compo, “A Practical Guide to Wavelet Analysis”, Bulletin of the American Meteorological Society, vol. 79, issue 1, pp. 6178, 1998
 [27] S. Vaughan, R. Edelson, R. S. Warwick and P. Uttley, “On characterizing the variability properties of Xray light curves from active galaxies”, Monthly Notices of the Royal Astronomical Society, vol. 345, issue 4, pp. 12711284, 2003
 [28] M. van der Klis, “Quasiperiodic oscillations and noise in lowmass Xray binaries”, Annual review of astronomy and astrophysics. Volume 27. Palo Alto, CA, Annual Reviews, Inc., pp. 517553, 1989
 [29] D. L. Gilman, F. J. Fuglister and J. M. Mitchell Jr., “On the Power Spectrum of Red Noise”, Journal of Atmospheric Sciences, vol. 20, issue 2, pp. 182184, 1963
 [30] S. L. Marple, Jr., “Digital Spectral Analysis with Applications”, Prentice Hall, Englewood Cliffs, 1987, Chapter 8
 [31] C. Chatfield, “The Analysis of Time Series. An Introduction.”, Chapman & Hall, A CRC Press Company, Sixth Ed., 2004
 [32] S. M. Kay, “Modern Spectral Estimation: Theory and Application”. Englewood Cliffs, NJ: PrenticeHall, 1988
 [33] P. Uttley, I. M. McHardy and S. Vaughan, “Nonlinear Xray variability in Xray binaries and active galaxies”, Monthly Notices of the Royal Astronomical Society, vol. 359, p. 345, 2005
 [34] K. Pottschmidt, M. Koenig, J. Wilms and R. Staubert, “Analyzing shortterm Xray variability of Cygnus X1 with Linear State Space Models”, Astronomy and Astrophysics, vol. 334, p. 201, 1998
 [35] Z. Ge, “Significance test for wavelet power and the wavelet power spectrum”, Ann. Geophys. , vol. 25, p. 2259, 2007
 [36] D. Maraun and J. Kurths, “Cross Wavelet Analysis. Significance Testing and Pitfalls”, Nonlin. Proc. Geoph., vol. 11(4), pp. 505514, 2004
 [37] D. Maraun, J. Kurths and M. Holschneider, “Nonstationary Gaussian processes in wavelet domain: Synthesis, estimation, and significance testing”, Physical Review E, vol. 75, issue 1, id. 016707, 2007
 [38] G. F. Smoot, et al., “Structure in the COBE differential microwave radiometer firstyear maps”, Astrophysical Journal, vol. 396, L1 , 1992
 [39] R. Fabbri and S. Torres, “Peak statistics on COBE maps”, Astronomy and Astrophysics, vol. 307, pp. 703707, 1996
 [40] J. Timmer and M. Koenig, “On generating power law noise”, Astronomy and Astrophysics, vol. 300, p. 707, 1995
 [41] P. Uttley and I. M. McHardy, “The fluxdependent amplitude of broadband noise variability in Xray binaries and active galaxies”, Monthly Notices of the Royal Astronomical Society, vol. 323, L26, 2001
 [42] W. H. Press, S. A. Teukolsky, W. T. Vetterling W.T. and B. P. Flannery B.P., “Numerical recipes in FORTRAN. The art of scientific computing”, Cambridge: University Press, 2nd ed., 1992
 [43] M. van der Klis, “Quantifying Rapid Variability in Accreting Compact Objects”, Statistical Challenges in Modern Astronomy II, 321, 1997
 [44] Ch. Clapham, The Concice Oxford Dictonary of Mathematics, Oxford University Press, 1996
 [45] D. A. Smith, “XTE J1550564”, IAU Circ., 7008, 1, 1998
 [46] J. A. Orosz, et al. “Dynamical Evidence for a Black Hole in the Microquasar XTE J1550564”, The Astrophysical Journal, vol. 568, p. 845, 2002
 [47] W. Cui, S. N. Zhang, W. Chen and E. H. Morgan, “Strong aperiodic Xray variability and quasiperiodic oscillationin Xray nova XTE J1550564, The Astrophysical Journal Letters, vol. 512, L43, 1999
 [48] R. A. Remillard, J. E. McClintock, G. J. Sobczak, C. D. Bailyn, J. A. Orosz, E. H. Morgan and A. M. Levine, “XRay Nova XTE J1550564: Discovery of a Quasiperiodic Oscillation near 185 Hz”, The Astrophysical Journal Letters, vol. 517, L127, 1999
 [49] M. Gierliński, M. Middleton, M. Ward, and C. Done, “A periodicity of 1hour in Xray emission from the active galaxy RE J1034+396”, Nature, vol. 455, issue 7211, pp. 369371, 2008
 [50] A. A. Zdziarski and M. Gierliński, “Radiative Processes, Spectral States and Variability of BlackHole Binaries”, Progress of Theoretical Physics Supplement, No. 155, pp. 99119, 2004
 [51] K. A. Arnaud, K.A., Astronomical Data Analysis Software and Systems V, eds. Jacoby G. and Barnes J., p. 17, ASP Conf. Series vol. 101, 1996