Multivariate Hadamard self-similarity: testing fractal connectivity
Herwig Wendt, Gustavo Didier, S\'ebastien Combrexelle and, Patrice Abry
TL;DR
This paper introduces Hadamard fractional Brownian motion to model multivariate signals with scale invariance but varying inter-component fractal connectivity, proposing wavelet-based estimators and tests validated through simulations.
Contribution
It presents a new class of multivariate Gaussian processes and develops wavelet-based methods for estimating and testing fractal connectivity in multivariate signals.
Findings
Wavelet-based estimator for Hurst matrix is asymptotically normal.
Finite sample procedures for confidence intervals and testing are effective.
Simulation validates the estimation and testing methods across parameters.
Abstract
While scale invariance is commonly observed in each component of real world multivariate signals, it is also often the case that the inter-component correlation structure is not fractally connected, i.e., its scaling behavior is not determined by that of the individual components. To model this situation in a versatile manner, we introduce a class of multivariate Gaussian stochastic processes called Hadamard fractional Brownian motion (HfBm). Its theoretical study sheds light on the issues raised by the joint requirement of entry-wise scaling and departures from fractal connectivity. An asymptotically normal wavelet-based estimator for its scaling parameter, called the Hurst matrix, is proposed, as well as asymptotically valid confidence intervals. The latter are accompanied by original finite sample procedures for computing confidence intervals and testing fractal connectivity from one…
| 2 |
| var | 65.33 | 9.51 | 1.87 | 0.42 | 69.27 | 10.05 | 1.97 | 0.44 | 18.76 | 2.73 | 0.54 | 0.12 | 12.42 | 1.81 | 0.35 | 0.08 |
| term 1 | 0.74 | 0.68 | 0.66 | 0.66 | 0.69 | 0.65 | 0.63 | 0.62 | 0.72 | 0.67 | 0.65 | 0.64 | 0.71 | 0.66 | 0.65 | 0.64 |
| term 2 | 0.15 | 0.14 | 0.14 | 0.14 | 0.20 | 0.19 | 0.18 | 0.18 | 0.17 | 0.16 | 0.16 | 0.16 | 0.18 | 0.17 | 0.16 | 0.16 |
| term 3 | 0.11 | 0.17 | 0.20 | 0.20 | 0.11 | 0.17 | 0.19 | 0.19 | 0.11 | 0.17 | 0.19 | 0.20 | 0.11 | 0.17 | 0.19 | 0.2 |
| co | 54.47 | 7.92 | 1.55 | 0.35 | 33.12 | 4.82 | 0.95 | 0.21 | 34.1 | 4.95 | 0.97 | 0.22 | ||||
| term 1 | 0.72 | 0.66 | 0.65 | 0.64 | 0.73 | 0.67 | 0.66 | 0.65 | 0.71 | 0.66 | 0.64 | 0.63 | ||||
| term 2 | 0.18 | 0.17 | 0.16 | 0.16 | 0.16 | 0.15 | 0.15 | 0.15 | 0.19 | 0.18 | 0.17 | 0.17 | ||||
| term 3 | 0.11 | 0.17 | 0.19 | 0.20 | 0.11 | 0.17 | 0.19 | 0.20 | 0.11 | 0.17 | 0.19 | 0.2 | ||||
| Ratio of | theo/MC | est/MC | |||||||
|---|---|---|---|---|---|---|---|---|---|
| Var | |||||||||
| Var | |||||||||
| Cov | |||||||||
| Var | |||||||||
| Var | |||||||||
| Var | |||||||||
| Cov | |||||||||
| Var | |||||||||
| Var | |||||||||
| Var | |||||||||
| Cov | |||||||||
| Var | |||||||||
| HFBM – , | ||||||
| WCF – , | ||||||
| HFBM – , | ||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| WCF – , | ||||||||||||
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Multivariate Hadamard self-similarity:
testing fractal connectivity
Herwig Wendt1, Gustavo Didier2, Sébastien Combrexelle1, Patrice Abry3
1 IRIT-ENSEEIHT, CNRS (UMR 5505), Université de Toulouse, France, [email protected]
2 Mathematics Department, Tulane University, New Orleans, USA, [email protected]
3 Université de Lyon, ENS de Lyon, Université Claude Bernard, CNRS, Laboratoire de Physique, F-69342 Lyon, France, [email protected]
Abstract.
While scale invariance is commonly observed in each component of real world multivariate signals, it is also often the case that the inter-component correlation structure is not fractally connected, i.e., its scaling behavior is not determined by that of the individual components. To model this situation in a versatile manner, we introduce a class of multivariate Gaussian stochastic processes called Hadamard fractional Brownian motion (HfBm). Its theoretical study sheds light on the issues raised by the joint requirement of entry-wise scaling and departures from fractal connectivity. An asymptotically normal wavelet-based estimator for its scaling parameter, called the Hurst matrix, is proposed, as well as asymptotically valid confidence intervals. The latter are accompanied by original finite sample procedures for computing confidence intervals and testing fractal connectivity from one single and finite size observation. Monte Carlo simulation studies are used to assess the estimation performance as a function of the (finite) sample size, and to quantify the impact of omitting wavelet cross-correlation terms. The simulation studies are shown to validate the use of approximate confidence intervals, together with the significance level and power of the fractal connectivity test. The test performance and properties are further studied as functions of the HfBm parameters. Keywords. Multivariate self-similarity, Hadamard self-similarity, fractal connectivity, wavelet analysis, confidence intervals, hypothesis testing
P.A. and H.W. are supported in part by the French ANR grant MultiFracs. G.D. was supported in part by the ARO grant W911NF-14-1-0475. G.D. gratefully acknowledges the support of ENS de Lyon for his visits, and thanks Alexandre Belloni for the insightful mathematical discussions.
1. Introduction
1.1. Scale Invariance
The relevance of the paradigm of scale invariance is evidenced by its successful use, over the last few decades, in the analysis of the dynamics in data obtained from a rather diverse spectrum of real world applications. The latter range from natural phenomena – physics (hydrodynamic turbulence (?), out-of-equilibrium physics), geophysics (rainfalls), biology (body rhythms (?), heart rate (?, ?), neurosciences and genomics (?, ?, ?, ?)) – to human activity – Internet traffic (?, ?), finance (?), urban growth and art investigation (?, ?, ?).
In essence, scale invariance – also called scaling, or scale-free dynamics – implies that the phenomenical or phenomenological dynamics are driven by a large continuum of equally important time scales, rather than by a small number of characteristic scales. Thus, the investigation’s focus is on identifying a relation amongst relevant scales rather than picking out characteristic scales.
Historically, self-similarity was one of the first proposed mathematical frameworks for the modeling of scale invariance (e.g., (?)). A random system is called self-similar when dilated copies of a single signal are statistically indistinguishable, namely,
[TABLE]
where stands for the equality of finite-dimensional distributions. An example of a stochastic process that satisfies the property (1.1) is fractional Brownian motion (fBm). Indeed, the latter is the only self-similar, Gaussian, stationary increment process, and it is the most widely used scaling model for real-world signals (?).
Starting from (1.1), the key parameter for quantifying scale-free dynamics is the scaling, or Hurst, exponent . The estimation of is the central task in scaling analysis, and it has received considerable effort and attention in the last three decades (see (?) for a review). The present contribution is about wavelet-based estimation (?, ?). It relies on the key scaling property
[TABLE]
where is the wavelet coefficient of an underlying self-similar stochastic process and is the number of available coefficients. In other words, the sample wavelet variance of the stochastic process behaves like a power law with respect to the scale .
1.2. Multivariate scaling
In many modern fields of application such as Internet traffic and neurology, data is collected in the form of multivariate time series. Univariate-like analysis in the spirit of (1.2) – i.e., independently on each component – does not account for the information stemming from correlations across components. The classical fBm parametric family, for example, provides at best a model for component-wise scaling, and thus cannot be used as the foundation for a multivariate modeling paradigm.
To model self-similarity in a multivariate setting, a natural extension of fBm, called Operator fractional Brownian motion (OfBm), was recently defined and studied (see (?, ?, ?)). An OfBm satisfies the -variate self-similarity relation
[TABLE]
where the scaling exponent is a matrix , and stands for the matrix exponential . Likewise, the wavelet spectrum of each individual component is not a single power law as in (1.2); instead, it behaves like a mixture of distinct univariate power laws. In its most general form, OfBm remains scarcely used in applications; recent efforts have tackled many difficulties that arise in the identification of its parameters (?, ?).
1.3. Entry-wise multivariate scaling
We call an OfBm entry-wise scaling when the Hurst parameter is simply a diagonal matrix . This instance of OfBm has been used in many applications (e.g., (?, ?)) and its estimation is thoroughly studied in (?). Since is diagonal, the relation (1.3) takes the form
[TABLE]
which is reminiscent of the univariate case. This implies that the extension of (1.2) to all auto- and cross-components of -variate data can be written as
[TABLE]
where is as in (1.2).
1.4. Fractal connectivity
Yet, entry-wise scaling OfBm is a restrictive model since the cross-scaling exponents , , are determined by the auto-scaling exponents and , i.e.,
[TABLE]
In this situation, called fractal connectivity (?, ?, ?), no additional scaling information can be extracted from the analysis of cross-components. However, in real world applications, cross-components are expected to contain information on the dynamics underlying the data, e.g., cross-correlation functions (?, ?). As an example, recent investigation of multivariate brain dynamics in (?) produced evidence of departures from fractal connectivity, notably for subjects carrying out prescribed tasks. This means that there is a clear need for more versatile models than entry-wise scaling OfBm (see also Remark 2.2 below). The covariance structure of the new model should satisfy the following two requirements:
- (1)
all auto- and cross-components are (approximately) self-similar; 2. (2)
departures from fractal connectivity are allowed, i.e., the exponents of the cross-components are not necessarily determined by the exponents of the corresponding auto-components.
Hereinafter, a departure from fractal connectivity (1.6) on a given covariance structure entry will be quantified by means of the parameter
[TABLE]
where nonnegativeness is a consequence of the Cauchy-Schwarz inequality (see (2.9) below). It is clear that when .
1.5. Goals, contributions and outline
Our contribution comprises four main components. First, we propose a new class of multivariate Gaussian stochastic processes, called Hadamard fractional Brownian motion (HfBm), that combines scale-free dynamics and potential departures from fractal connectivity. Moreover, we provide a precise discussion of the issues entailed by the presence of these two properties (Section 2). Second, we study the multivariate discrete wavelet transform (DWT) of HfBm, define wavelet-based estimators for the scaling exponents and the fractal connectivity parameter , mathematically establish their asymptotic performance (i.e., asymptotic normality and covariance structure), and computationally quantify finite sample size effects (Section 3). Third, starting from a single sample path, we construct approximate confidence intervals for the proposed estimators. The procedure is backed up by the mathematical identification of the approximation orders as a function of the sample path size and the limiting coarse scale. This is further investigated by means of Monte Carlo simulations, as well as by means of a study of the ubiquitous issue of the impact of (partially) neglecting the correlation amongst wavelet coefficients (Section 4). Beyond being of interest in multivariate modeling, the study sheds more light on the same issue for the univariate case. Fourth, we devise an efficient test for the presence of fractal connectivity from a single sample path. In addition, we assess the finite sample performance of the test in terms of targeted risk by means of Monte Carlo simulations (Section 5). Finally, routines for the synthesis of HfBm, as well as for estimation, computation of confidence intervals and testing will be made publicly available at time of publication. All proofs can be found in the Appendix.
2. Hadamard fractional Brownian motion
For Hadamard fractional Brownian motion, defined next, the fractal connectivity relation (1.6) does not necessarily hold.
Definition 2.1**.**
A Hadamard fractional Brownian motion (HfBm) is a proper, Gaussian (stationary increment) process whose second moments can be written as
[TABLE]
For , the matrix exponent H=\Big{(}h_{q_{1}q_{2}}\Big{)}_{q_{1},q_{2}=1,\ldots,m} satisfies the conditions
[TABLE]
The matrix-valued function is a spectral density of the form
[TABLE]
i.e., the Hadamard scaling parameters are given by
[TABLE]
where , . The real-valued functions satisfy
[TABLE]
[TABLE]
for constants , where
[TABLE]
Example 2.1**.**
An HfBm with parameters (i.e., fractally connected), , , is an entry-wise scaling OfBm with diagonal Hurst matrix (see (?, ?, ?)).
By the known properties of spectral densities, (symmetric positive semidefinite) a.e. and satisfies
[TABLE]
((?), p.436). The relation (2.8) further implies that
[TABLE]
Whenever convenient we also write
[TABLE]
The name “Hadamard” comes from Hadamard (entry-wise) matrix products. If one rewrites HfBm componentwise as , then the conditions (2.1), (2.2) and (2.6) yield the asymptotic equivalence
[TABLE]
, for some . In other words, over large scales, the covariance between each pair of entries of an HfBm approaches that of a univariate fBm, up to a change of sign (see also Proposition 3.1, ()). In this sense, for large , an HfBm behaves like its ideal counterpart , defined as a generally non-existent stochastic process satisfying the also ideal Hadamard (entry-wise) self-similarity relation
[TABLE]
where denotes the Hadamard (entry-wise) matrix product and c^{\circ H}:=\Big{(}c^{h_{q_{1}q_{2}}}\Big{)}_{q_{1},q_{2}=1,\ldots,m}. The process can be viewed as a heuristic tool for developing intuition on multivariate self-similarity. Mathematically, though, it can only exist in fractally connected instances, the reason being that distinct (spectral) power laws cross over. Indeed, we must have for close to 0 and for large , whence . This shows that HfBm is a perturbation of its virtual counterpart, where the regularization functions in (2.3) introduce high-frequency corrections.
Example 2.2**.**
An illustrative subclass of HfBm is obtained by setting , , and , . Note that satisfies (2.7) with . In this case, the expression for the main diagonal spectral entry of an HfBm is identical to that of an ideal-HfBM, and the difference lies on the off-diagonal entries:
[TABLE]
In this case, each individual entry , , of HfBm is by itself a fBm with Hurst parameter . In particular,
[TABLE]
However, it is generally not true that
[TABLE]
Otherwise, would necessarily be fractally connected.
Remark 2.1**.**
When simulating HfBm via Circulant Matrix Embedding, one verifies that regularization is rarely necessary for ensuring the positive definiteness of the covariance matrix for finite sample sizes. In other words, ideal-HfBM is also a useful approximation in practice.
Remark 2.2**.**
Let , where is an HfBm. Then, is a (discrete time) stationary process with spectral density
[TABLE]
where is the HfBm spectral density (2.1). Under (2.5), (2.6) and (2.7), we obtain the entry-wise limiting behavior
[TABLE]
In particular, does not fall within the scope of the usual definitions of multivariate scaling behavior or long range dependence (?, ?, ?, ?), which are restricted to fractally connected processes.
3. Wavelet-based analysis of HfBm
In this section, we construct the wavelet analysis and estimation for HfBm. Due to the mathematical convenience of the notion of Hadamard (approximate) self-similarity, most of the properties of wavelet-based constructs resemble their univariate analogues.
3.1. Multivariate discrete wavelet transform
Throughout the rest of the paper, we will make the following assumptions on the underlying wavelet basis. Such assumptions will be omitted in the statements.
Assumption : is a wavelet function, namely,
[TABLE]
Assumption (): the functions (a bounded scaling function) and correspond to
[TABLE]
Assumption : for some ,
[TABLE]
and
[TABLE]
Under (3.1)–(3.3), is continuous, is everywhere differentiable and
[TABLE]
(see (?), Theorem 6.1 and the proof of Theorem 7.4).
Definition 3.1**.**
Let be an HfBm. For a scale parameter and a shift parameter , its (-normalized) wavelet transform is defined by
[TABLE]
Under (3.1)-(3.3) and the continuity of the covariance function (2.1), the wavelet transform (3.5) is well-defined in the mean-square sense and , , (see (?), p. 86).
3.2. Multivariate wavelet spectrum
Fix , . Because of the approximate nature of Hadamard self-similarity, analysis and estimation must be considered in the coarse scale limit, by means of a sequence of dyadic numbers satisfying
[TABLE]
Example 3.1**.**
An example of a scaling sequence satisfying (3.6) for large enough is
[TABLE]
In other words, a wide parameter range or a low regularity parameter value implies that must grow slowly by comparison to .
Remark 3.1**.**
For a fixed octave range associated with an initial scaling factor value and sample size , define the scale range for a general sample size (where is assumed dyadic). Then, under (3.6), for every the range of useful octaves is constant and given by , where the new octaves are , .
Definition 3.2**.**
Let be an HfBm. Let be its wavelet coefficients, where
[TABLE]
is the number of wavelet coefficients at scale , of a total of (wavelet) data points. The wavelet variance and the sample wavelet variance, respectively, at scale are defined by
[TABLE]
Let
[TABLE]
The standardized counterparts of (3.7) are
[TABLE]
, where \varrho(a(n)2^{j})=\Big{(}\varrho^{(q_{1}q_{2})}(a(n)2^{j})\Big{)}_{q_{1},q_{2}=1,\ldots,m}. Entry-wise, for ,
[TABLE]
Proposition 3.1, stated next, provides basic results on the moments and asymptotic distribution, in the coarse scale limit , of the wavelet transform (3.5) and variance (3.9). In particular, in regard to the limits in distribution, the vector of random matrices can be intuitively interpreted as an asymptotically unbiased and Gaussian estimator of its population counterpart . The celebrated decorrelation property of the wavelet transform (e.g., (?), Proposition II.2) lies behind the Gaussian limits in the proposition, as well as of the fact that the random matrices have weak inter-component, intra-scale and inter-scale correlations. Also, each entry , , displays asymptotic power law scaling. In Section 3.3 below, these properties are used to define estimators of the scaling exponents and to analytically establish their asymptotic performance.
In the statement of Proposition 3.1, we make use of the operator
[TABLE]
In other words, vectorizes the upper triangular entries of the symmetric matrix .
Proposition 3.1**.**
Let be an HfBm. Consider and . Then,
- ()
for fixed , , , , we can write
[TABLE]
for some matrix-valued function that depends on . In particular, for fixed scales , , is a stationary stochastic process;
- ()
for ,
[TABLE]
i.e., (3.9) is well-defined;
- ()
for as in (3.12), , and , as ,
[TABLE]
In particular, the wavelet spectrum (see (3.7)) can be approximated by that of an ideal-HfBM (see (2.10)) in the sense that
[TABLE]
- ()
for as in (3.12) and , as ,
[TABLE]
[TABLE]
where is given by (3.14);
- (v)
for , , as ,
[TABLE]
for as in (1.7), where
[TABLE]
[TABLE]
and
[TABLE]
- (vi)
as ,
[TABLE]
where
[TABLE]
and is a vector of ones. Each entry of the asymptotic covariance matrix in (3.19) is given by the terms in (3.17) for appropriate values of , , , , and .
Remark 3.2**.**
Note that, up to a change of sign, each entry on the right-hand side of (3.14) corresponds to the covariance between the wavelet coefficients at octaves and of a fBm with parameter .
Remark 3.3**.**
In Proposition 3.1, , the asymptotic covariance between same-entry wavelet variance terms is always nontrivial, irrespective of the values of fractal connectivity parameters. For example, when , it is given by
[TABLE]
This is not true for terms associated with different pairs of indices. For example, if , then
[TABLE]
In other words, the phenomenon of the asymptotic decorrelation of wavelet variance terms is only observed for instances involving departures from fractal connectivity.
3.3. Estimation of the scaling exponents
3.3.1. Definition of the estimators
As in the univariate case, the fact that the wavelet variance (3.7) satisfies the Hadamard scaling relation in the coarse scale limit (see (3.15)) points to the development of a log-regression regression method based on the sample wavelet variance across scales.
Estimators can be defined in a standard way by means of the log-regression relations
[TABLE]
[TABLE]
which are well-defined with probability 1, where , , are linear regression weights satisfying
[TABLE]
(e.g., (?, ?)). The derived estimator
[TABLE]
for the parameter (1.7) will be used in Section 5 in the construction of a hypothesis test for fractal connectivity (see also (?, ?, ?)). With the standard Gaussian asymptotics for the sample wavelet variances established in Proposition 3.1, the system of equations (3.21) and (3.23) is expected to yield efficient estimators. In fact, this is proved in the next section.
3.3.2. Asymptotic distribution
In the next result, Theorem 3.1, we build upon Proposition 3.1 to show that the vector is an asymptotically unbiased and Gaussian estimator of the scaling exponents (whence the analogous statement holds for the estimator (3.23)). As in the aforementioned proposition, the decorrelation property of the wavelet transform contributes to the Gaussianity of the asymptotic distribution.
Theorem 3.1**.**
Let be the weight terms in (3.21) and let be as in (3.19). Suppose
[TABLE]
(see (2.3)). Then, as ,
[TABLE]
Remark 3.4**.**
Condition (3.24) is needed to ensure the consistency of the estimator. Furthermore, it is clear that for , since otherwise the process is not proper or its spectral density is not positive semidefinite a.e.
Remark 3.5**.**
The convergence rate in (3.25) depends on the unknown fractal connectivity parameters . In practice, the latter can be replaced by their estimates or, in some cases, ignored, since they exponentiate the slow growth term .
3.3.3. Numerical simulation setting
We conducted Monte Carlo experiments over 1000 independent realizations of bivariate HfBm. Though several parameter settings were tested, results are reported only for two representative cases: fractally connected ideal-HfBm ; and non-fractally connected HfBm (see Example 2.2) with . HfBm copies were synthesized using the multivariate process synthesis toolbox described in (?, ?) and available at www.hermir.org. For the wavelet analysis, least asymmetric Daubechies wavelet with vanishing moments were used (?). Estimation by means of weighted linear regression was performed on the octave range . This choice of regression range (with fine scale fixed and ) is used here in order to, first, obtain results that are consistent and comparable with those reported in the literature for the estimation of univariate scaling exponents, cf., e.g., (?) and references therein, and, second, avoid the issue of fixing , and in (3.6).
Monte Carlo parameter estimates are reported in Figures 1 and 2 in terms of bias, standard deviation, skewness and kurtosis, as functions of (the of) sample sizes .
3.3.4. Finite sample performance
Figures 1 and 2 show that for all four estimates , the bias becomes negligible as the sample size grows. It can also be seen that the bias is slightly larger for non-fractally connected data, notably for the parameter .
The figures further show that standard deviations for all estimates decrease as (the latter trend being plotted as superimposed dashed lines), with no significant difference between fractally and non-fractally connected instances. For , there is no detectable dependence of the outcomes on the actual values of the parameters themselves (in accordance with theoretical calculations reported in Table 3). Accordingly, one observes that the variances for and are identical and depend on the sample size, but not on any parameter of the stochastic process. For , the variance does not depend on , but results not shown demonstrate that it does vary with , as confirmed by approximate calculations (see (4.9)). For , further simulations not displayed indicate that Var depends on roughly proportionally to . Departures from fractal connectivity tend to imply a decrease in the variance of ; this is a counterintuitive phenomenon that has yet to be fully explained.
In regard to convergence in distribution, Figures 1 and 2 show that the Monte Carlo skewness and kurtosis estimates for are very close to [math] even at very small sample sizes, both with and without fractal connectivity. This characterizes a fast convergence to limiting normal distributions, and lies in agreement with the asymptotics presented in Theorem 3.1. However, the weak convergence is observed to be much slower in practice for the cross-exponents , yet with no noticeable difference between fractally and non-fractally connected instances, i.e., different values of the parameter .
4. Confidence intervals
To complement the asymptotic and finite sample results in Sections 3.3.2 and 3.3.4, we now construct confidence intervals for the estimators (3.21) and (3.23). In particular, we investigate the effect of omitting the covariance among sample wavelet variance terms, as well as the dependence of the asymptotic variance in (3.25) and confidence intervals on the unknown parameter values. This will allow us to efficiently compute approximate confidence intervals in practice based on a single sample path.
4.1. Asymptotic covariance
In Theorem 4.1 below, we provide a finite sample approximation for the covariance between the logarithms of sample wavelet variances, coupled with the convergence rate of the covariance approximation to its asymptotic expression in terms of ideal-HfBm covariances (3.14). As in Theorem 3.1, in part due to the decorrelation property of the wavelet transform, the estimators’ second order properties are comparable to those observed in the i.i.d. case (e.g., (?), p. 430). In particular, the finite sums provided in Theorem 4.1 mathematically justify heuristic expressions obtained from Taylor expansions (c.f. (4.5) below), and the theorem also accurately quantifies the approximation error as a function of the sample size and scale.
In the theorem, we use indicator functions to regularize around the origin. The truncation converges fast to zero and is of no impact in practice. Indeed, one can redefine the estimator (3.21) with the indicators and obtain the same asymptotic distribution (3.25).
Theorem 4.1**.**
Let be generic indices, and suppose condition (3.24) holds. Consider , as in (3.10) and as in (3.14). For notational simplicity, denote
[TABLE]
and also recall the definition of in (2.7). Consider a sequence satisfying
[TABLE]
for any and any . Then,
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
Remark 4.1**.**
Note that the finite summations on the right-hand side of (4.3) converge as . For example,
[TABLE]
where is defined by (3.18) (this can be shown by the same argument for establishing (A.16) and (A.25) in the proof of Proposition 3.1, ).
Remark 4.2**.**
In (4.3), the main residual term satisfies
[TABLE]
under the condition (3.6). In practice, the modeling of instances involving extreme deviations from fractal connectivity ( close to 1) may require greater regularity from the functions (see (2.3)) or the wavelet basis. Otherwise, for conservative choices of the regularity parameters (e.g., and for small ), the expression (4.4) may converge slower to zero than the main terms in (4.3).
Remark 4.3**.**
In (3.5), we assume that an HfBm continuous time sample path is available. However, it is well known that, under mild assumptions, the availability of discrete time observations does not generally alter the nature of the asymptotic distribution for wavelet estimators (see, for instance, (?), Section III; (?), Section 3.2; (?), Section C). Establishing this for the case of HfBm is a topic for future work.
Remark 4.4**.**
When , for the univariate case the expansion (4.3) appears implicitly in (?), expression (86). Also, when computing the moments of , the truncation is unnecessary (see Remark C.2 below).
4.2. Variances and covariances of and
4.2.1. Closed-form approximations
Turning back to expression (4.3), by ignoring the scaling factor (), the truncation and setting all the convergence order terms to zero we obtain the approximation
[TABLE]
[TABLE]
[TABLE]
where
[TABLE]
Note that the approximation (4.5) is a finite sample one, and ignores the potential asymptotic decorrelation effect stemming from the shifting scaling factor (see Proposition 3.1, , and (4.3)). Consider the normalization
[TABLE]
Expressions (3.21), (4.5) and (4.6) yield the covariance approximation
[TABLE]
In particular, (4.7) further allows us to compute
[TABLE]
The variances (4.8) will be used in Section 5 in the construction of a test for fractal connectivity, i.e., for the hypothesis . Table 1 summarizes the closed-form approximations for , , and established in (4.7).
4.2.2. Impact of inter- and intra-scale correlations
The expression of can further be split into three terms, namely,
[TABLE]
where the first term in the sum over reflects the variance only of wavelet coefficients, the second term the covariance of wavelet coefficients at a given scale, and the third term, the covariance of wavelet coefficients at different scales. In other words, if wavelet coefficients were independent, the second and third terms would equal zero.
The relative contributions of the three terms to the final variances are quantified by means of Monte Carlo simulations conducted following the same protocol and settings as those described in Section 3.3.3. Table 2, reporting the relative contributions of each of the three terms for various sample sizes under fractal connectivity (i.e., ), clearly shows that the second and third terms (intra- and inter-scale covariances) cannot be neglected, namely, one cannot use only the first term (variance) in the construction of confidence intervals for fBm, as proposed in (?). Identical conclusions, not shown here, are drawn under departures from fractal connectivity (i.e., ).
4.2.3. First order approximations for the variances and covariances of
It is of interest to further examine the leading order approximations for the variances and covariances of and , corresponding to neglecting all intra- and inter-scale correlations amongst wavelet coefficients. The first order approximations neglecting all inter- and intra- scale correlations amongst wavelet coefficients of and of are summarized in Table 3.
The results show that and , , do not depend on the actual values of the scaling exponents , which corroborates the numerical performance reported in Section 3.3.4.
Note that for an ideal-HfBm with fractal connectivity, it is straightforward to show that and when . While does not depend on correlations , as expected, does vary with according to , showing that when . This can be interpreted as the fact that when , the scaling exponent , , loses its meaning. Furthermore, depends on as , not surprisingly indicating that when (no correlation amongst components), (no correlation amongst estimates).
Moreover, the first order approximation of is observed not to depend on the actual value of , while clearly depends on based on the numerical simulations reported in Section 3.3.4. This can be interpreted as the fact that, for , the first order approximation (neglecting all intra- and inter-scale correlations amongst wavelet coefficients) is not sufficient to approximate well , as opposed to what is observed for and .
To finish with, varies with as . This shows again that when , parameter becomes irrelevant. Moreover, it also shows that when , . This can be understood as the fact that when , a departure from fractal connectivity is no longer permitted, as indicated by (2.11). Thus, implies , which is then no longer a random variable.
4.3. Practical computation of the variances and covariances of and
4.3.1. Computation of and
Evaluation of (4.7) requires knowledge of the covariance between wavelet coefficients (4.6), which will be developed here explicitly for the discrete and the dyadic wavelet transforms. Let and , , be the coefficients of the high pass and low pass filters of the discrete wavelet transform, respectively, and let and be the dyadic upsampling and decimation operators. The wavelet transform of the discrete-time process component yields, at each scale , sequences of approximation coefficients and detail coefficients , . The corresponding dyadic coefficients are given by and , respectively. Pick the initialization (see (?) for a discussion of the initialization of the discrete wavelet transform). At scale , a_{q}(1,\cdot)=h\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}a_{q}(0,\cdot) and d_{q}(1,\cdot)=g\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}a_{q}(0,\cdot), at scale , a_{q}(2,\cdot)=\uparrow_{2}[h]\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}a_{q}(1,\cdot)=\uparrow_{2}[h]\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}h\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}a_{q}(0,\cdot) and d_{q}(2,\cdot)=\uparrow_{2}[g]\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}a_{q}(1,\cdot)=\uparrow_{2}[g]\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}h\mathop{\scalebox{1.2}{\raisebox{-0.08612pt}{\ast}}}a_{q}(0,\cdot). By iteration, we obtain the sequences of detail coefficients at each scale , i.e.,
[TABLE]
where
[TABLE]
Now let , , be the covariances of univariate fBm of indices ,
[TABLE]
Then (cf. (?); note that a change in time scale would only result in a multiplicative constant, which we assume to be absorbed in and which cancels out in the final expression for )
[TABLE]
where
[TABLE]
Consequently,
[TABLE]
and, for the dyadic wavelet transform,
[TABLE]
For given values of and , these expressions can be easily evaluated numerically.
4.3.2. Practical estimation of (co)variances of and of
Evaluating (4.7) and (4.8) for HfBm in practice requires the unknown parameter values and in (4.11), and hence we replace them by their estimates and . The former are defined in (3.21), and estimates for for can be readily obtained as the cross-correlation coefficients of the difference processes (fGn) of the process components and However, note that the expressions for the (co)variances of in the previous sections are derived assuming knowledge of the true parameter values and can only be expected to be approximations when these are replaced by estimates. This will be studied numerically in the next section.
4.3.3. Assessment of the estimated (co)variances of and of by means of Monte Carlo experiments
Monte Carlo studies were conducted following the protocol and settings described in Section 3.3.3, aiming to evaluate the quality of the estimated approximations (4.7) and (4.8) for the (co)variances of and of . The simulations involved independent realizations of each of two general instances of HfBm with components, one using the true values of the parameters and , and the other, their estimates and . Four different sample sizes, and three different values are investigated for the set of exponents .
Table 4 summarizes the square roots of the ratios of the averages over realizations of (co)variance estimates and of the Monte Carlo (co)variances. The first four columns, labeled “theo/MC”, report results obtained when using theoretical parameter values and yield the following conclusions. First, even for small sample size and weak correlation , the quality of the approximations (4.7) and (4.8) is very good for the variances of exponents , , and satisfactory for the “cross” exponent , the covariance parameters and the connectivity parameter . Second, when the sample size and correlation level increase, the approximation of variances and covariances becomes excellent, with maximum errors of the order of for and strong correlation . Finally, the last four columns of Table 4, labeled “est/MC”, report results obtained when using estimates and . They indicate that replacing the true parameter values and with estimates has very little impact on the quality of approximations (4.7) and (4.8). Indeed, the average values of the (co)variance estimates are essentially equal to those obtained when using true parameter values.
5. Statistical test for fractal connectivity
5.1. Procedure
The mathematical and computational results in Sections 3 and 4 enable us to construct component-wise fractal connectivity tests, i.e., for the hypotheses
[TABLE]
We assume for . As a consequence of Theorem 3.1, the distribution of under can be approximated over finite samples by
[TABLE]
where, in turn, can be approximated by (4.8). Therefore, a simple two-sided test with significance level can be defined as
[TABLE]
where is the inverse cumulative distribution function of the standard Normal distribution. In addition, the -value of the test statistic (i.e., the probability of observing an absolute value at least as large as for the test statistic under ) is given by
[TABLE]
This test can be performed by evaluating (5.1) with an estimate for obtained from the procedure detailed in Section 4.1.
5.2. Monte Carlo assessment of the test performance
We assess the performance of the test by applying it to independent realizations of HfBm with exponent values and exponent values detailed below for sample sizes and correlation levels . For simplicity of illustration and without loss of generality, we consider again HfBm with components. For each realization, the test decision (5.1) and the -value (5.2) are evaluated using (4.8) with approximations (4.7) to obtain an estimate of the . Estimates of the expected values of the test decisions and -values, denoted by and , are then obtained as the averages over realizations of test decisions and -values (5.1) and (5.2).
We now compare the performance of the proposed test, denoted hereinafter HFBM (not to be confused with the stochastic process HfBm), to that of the test put forward in (?). The latter relies on the intuition that the wavelet coherence function of two components of a multivariate Gaussian scale invariant random process approximately behaves as
[TABLE]
The test itself, denoted WCF (for wavelet coherence function), is formulated without the rigorous statistical framework developed above. Rather, it is built on the observation that is the Pearson product-moment correlation coefficient of the time series and , and hence that the Fisher’s statistics of , , are approximately Gaussian, with known variances and, in the case of fractal connectivity, with equal means across scales. The test for fractal connectivity is then formulated as the UMPI test for the equality of means of Gaussian random variables, cf. (?) for details.
5.2.1. Performance under
We first consider the case that is true, i.e., . The significance level is set to , results are reported in Table 5 for the proposed test (top) and for the test in (?) (bottom). Note that, under , averages of test decisions should equal the preset significance level , and averages of -values should equal . HFBM rejects with slightly larger probability than the prescribed value , yet the differences between empirical significance levels and never exceed ; similarly, average -values are slightly below . For large sample size and large , average test decisions and -values are very close to the theoretical values and . These remarks are consistent with the results reported in Table 4, where a small but systematic underestimation of for small sample sizes and is observed.
In contrast, the empirical significances of WCF strongly differ from the preset value by values of up to , and this difference is especially pronounced for large sample sizes for which one would expect the test to perform better. One reason for this poor performance may lie in the fact that the test in (?) was designed for fGn, rather than fBm. Note that the asymptotic calculations developed above can be adapted to the easier case of fGn (and, in principle, any other Gaussian process with stationary increments) without difficulty by simply changing the covariance function in the calculations leading to the expressions (4.12) and (4.13).
5.2.2. Test power
We assess the power of the test under the alternative hypotheses with . Yet, a direct comparison of the power of HFBM and WCF is only meaningful for identical rejection probabilities under . Indeed, a test for which under is expected to display an artificially large power. Since the performance of HFBM and WCF differ, as discussed above (cf. Table 5), we therefore adjust, for each sample size and correlation level, the prescribed significance to the value for which the average rejection rate under equals . Using this adjusted level of significance , the power of the test is then estimated as the average of the test decisions when is true. Results are reported in Table 6 and yield the following conclusions. First, the power of each test systematically increases with the magnitudes of the deviation from , of the correlation level and of the sample size , as expected. Second, HFBM is systematically and significantly more powerful. Indeed, it enables us to detect a non-zero value for up to two times as often as WCF. For instance, for the small sample size of and the low correlation level of , it permits the detection of a deviation of from the null value with probability , as compared to a probability of for the test in (?).
Overall, these results confirm that the proposed developments can be relevantly applied in the assessment of scaling and fractal connectivity in multivariate time series.
6. Conclusion
The present contribution introduces a versatile class of multivariate stochastic processes, named Hadamard fractional Brownian motion (HfBm), for scale invariance modeling. HfBm provides a stochastic framework within which cross-component scaling laws are not directly controlled by the scaling along the main diagonal, namely, HfBm is not necessarily fractally connected.
Interestingly, the theoretical study of HfBm reveals that exact entry-wise scaling on both auto- and cross-components and departures from fractal connectivity are mathematically incompatible. In other words, there is a dichotomy in multivariate scaling modeling: either there is exact entry-wise scaling in every component combined with fractal connectivity, or departures from fractal connectivity are allowed at the price of approximate (i.e., asymptotic) scaling on the cross-components.
Our main mathematical results consist of an asymptotically normal, wavelet-based linear regression estimator for the scaling exponents, as well as asymptotically valid confidence intervals with convenient mathematical expressions. Furthermore, the Taylor expansions used in the development of the asymptotic confidence intervals lead to the construction of practical procedures for the numerical calculation of the variance of the estimates. These approximate calculations enable the study of the ubiquitous issue of the impact of neglecting intra- or inter-scale correlations amongst wavelet coefficients in the computations of variances and covariances for the estimates. We also devised an asymptotically normal hypothesis test for fractal connectivity. Again, a major feature of the designed test procedure is the fact that it can be applied to a single observed HfBm data path.
For both fractally and non-fractally connected instances, simulations demonstrate the satisfactory performance of the estimators of the scaling and fractal connectivity parameters, even for small sample size data. The estimation bias is shown to be negligible, and the variance decreases according to the inverse of the sample size. In addition, the practical computations of approximated variances and covariances of the estimates are shown to be of excellent quality, irrespective of sample size, and the Monte Carlo significance levels and powers are very close to their theoretical counterparts.
The tools developed in the present contribution pave the way for novel analysis and modeling perspectives on multivariate scaling in real-world data, in the spirit of the accomplishments in (?). Routines for the synthesis of HfBm, as well as for estimation, computation of confidence intervals and fractal connectivity testing will be made publicly available at time of publication.
Appendix A Section 3.2
Proof of Proposition 3.1: For simplicity, we will write instead of throughout the proof.
To show (), the change of variable in (3.5) and the harmonizable representation of HfBm yield
[TABLE]
where f(x):=\Big{(}|x|^{-2(h_{q_{1}q_{2}}+1/2)}g_{q_{1}q_{2}}(x)\Big{)}_{q_{1},q_{2}=1,\ldots,m}. Statement () is a consequence of the formula (A.1) with
[TABLE]
Statement follows from Lemma C.1, , below under condition (3.6).
Turning to , establishing (3.16) is equivalent to showing that
[TABLE]
[TABLE]
.
Consider the entry , of the matrix-valued function as in (A.1). Based on a change of variable , recast the first term of the sum on the left-hand side of (A.3) as
[TABLE]
Fix and denote , . By (3.1), (3.4) and a Taylor expansion with Lagrange residual of the real and imaginary parts of , there exist functions , on such that
[TABLE]
Therefore, and extending this reasoning to , ,
[TABLE]
For as in (A.4), we now show that
[TABLE]
We will only develop expressions for , since analogous developments hold for . For mathematical convenience, rewrite as in (A.4) as
[TABLE]
Hence,
[TABLE]
[TABLE]
[TABLE]
Note that, by conditions (3.1) and (2.7),
[TABLE]
Then, is also smooth around zero and . Moreover, by (A.5) and (A.10) with , for fixed and ,
[TABLE]
as , where the constant does not depend on . Similarly, by (2.6), (A.5) and (A.10) with ,
[TABLE]
[TABLE]
This proves (A.6), as desired. Next, note that
[TABLE]
and , by the continuity of and condition (3.2). Therefore, by the dominated convergence theorem, , for appropriate constants . Consequently, by (3.2),
[TABLE]
So, fix . By (2.5) and (A.13),
[TABLE]
Thus, in view of (A.4), (A.6), (A.8), (A.9), (A.10) (with ) and (A.14), by integrating by parts twice we obtain
[TABLE]
[TABLE]
where the last inequality is a consequence of (2.5), (A.5) and (A.13).
Now consider the first summation term in (A.3). We proceed as in the proof of Proposition 3.3, , in (?) to establish that
[TABLE]
We outline the main steps for the reader’s convenience. By Theorem 1.8 in Jones and Jones (?), p. 10, the set of values to the equation , , is given by . Therefore, a pair is a solution to
[TABLE]
for some if and only if it is a solution to
[TABLE]
for the same . Therefore, we can replace with in Lemmas B.2 and B.3, (?), and reexpress the first summation term on the left-hand side of (A.16) as
[TABLE]
In (A.19), is the range for such that the pairs satisfying (A.18) lie in the region
[TABLE]
and is the number of such solution pairs given some . Moreover, for any sufficiently large , let be the smallest number such that solves (A.17) (for some ), where
[TABLE]
From the proof of Lemma B.2, (?), the set of such solutions to (A.18) has the form
[TABLE]
In light of (A.22), define the function , . In particular, is a solution pair for (A.18). Let . Then, by (A.22), is the rightmost solution for (A.17) within the first-entry range . Moreover, given , the number of solution pairs in the region (A.20) is , where
[TABLE]
In addition, by (A.15), and (3.14),
[TABLE]
By expression (A.19), the dominated convergence theorem and (A.23), the limit (A.16) holds.
Next recall that, up to a change of sign, corresponds to the cross moment of the wavelet transform of a fBm. Thus, by an analogous procedure, we also obtain
[TABLE]
This establishes (A.3).
To show statement (), first note that
[TABLE]
which follows from (3.6) and the fact that
[TABLE]
Thus, by expression (C.42) for (established in the proof of Lemma C.5 below),
[TABLE]
[TABLE]
[TABLE]
[TABLE]
By (A.26) and (A.27), statement holds.
The argument for showing () is an adaptation of the proof of Theorem 3.1 in (?) (see also (?), pp.510-513; (?), p.997; and (?), Lemma 2). For notational simplicity, we only write the proof for , where the entries are indexed 1 and 2, and denote generic entries. The proof is by means of the Cramér-Wold device. Under (.2), form the vector of wavelet coefficients
[TABLE]
[TABLE]
where
[TABLE]
Let
[TABLE]
where and , . Now form the block-diagonal matrix defined by
[TABLE]
where
[TABLE]
In (A.32), it can be understood that the slow growth factors for the main diagonal terms are . We would like to establish the limiting distribution of the statistic
[TABLE]
[TABLE]
[TABLE]
where it suffices to consider in (A.30) such that
[TABLE]
(see (?), pp.211 and 214). The matrix in (A.33) is obtained from (3.17) and can be written in block form as , corresponding to block entries of the vector . Let , and consider the spectral decomposition , where is diagonal with real, and not necessarily positive, eigenvalues
[TABLE]
and is an orthogonal matrix. Now let . Then,
[TABLE]
Assume for the moment that
[TABLE]
[TABLE]
[TABLE]
Therefore, there exists a constant such that, for large enough , . In view of condition (A.35),
[TABLE]
The claim (3.19) is now a consequence of Lemma B.4 in (?).
So, we need to show (A.35). The first step is to establish the bound
[TABLE]
Let and let . We can break up the vector v into two-dimensional subvectors to reexpress . Based on the block-diagonal structure of expressed in (A.31),
[TABLE]
[TABLE]
where the constant comes from a change of matrix norms and only depends on the fixed dimension . By taking on both sides of (A.37), we arrive at (A.36).
The second step towards showing (A.35) consists of analyzing the asymptotic behavior of the right-hand side of (A.36), as . For this, we will assume the result of Lemma C.1 below. So, note that
[TABLE]
Moreover, by (C.14) below, the maximum eigenvalue of is bounded by , where is the matrix Euclidean norm. Therefore, in view of (A.38) and (C.14), the right-hand side of (A.36) is bounded by . In turn, the latter expression divided by is bounded by C\Big{(}\frac{a^{4(h_{\max}-h_{\min})+1}}{n}\Big{)}^{1/2}\frac{1}{a^{\delta_{12}}}. By condition (3.6), this implies (A.35), and as a result, also (3.19).
Remark A.1**.**
As a consequence of expression (A.27) and the fact that under (3.6),
[TABLE]
Appendix B Section 3.3
Proof of Theorem 3.1: Fix . Based on (3.22), rewrite
[TABLE]
[TABLE]
[TABLE]
By Lemma C.1, (below), and an application of the mean value theorem, for some between \Big{|}\mathbb{E}\!\left[S^{(q_{1}q_{2})}_{n}(a2^{j})\right]\Big{|}/a^{\alpha_{q_{1}q_{2}}} and the second term in the sum (B.1) can be bounded in absolute value by
[TABLE]
[TABLE]
as , where the limit is a consequence of condition (3.6). Also note that, after a change of variable in the expression for (see (3.14)),
[TABLE]
Therefore, by (3.22), the third term in the sum (B.1) can be written as
[TABLE]
So, in regard to the first term in the sum (B.1), consider the weight matrix defined by
[TABLE]
where is a identity matrix and is given by (3.20). We would like to show that
[TABLE]
where \log\circ|A|:=\Big{(}\log|A_{q_{1}q_{2}}|\Big{)}_{q_{1},q_{2}=1,\ldots,m} for any real matrix , and the term post-multiplying the matrix in (B.2) is a -dimensional random vector.
For any , fix a pair , and an octave , which specifies one of the entries of the random vector on the left-hand side of (B.2). Define the set . Under (.2), Proposition 3.1, , implies that as . Thus, for large enough , in the set the mean value theorem gives the almost sure expression
[TABLE]
for a random variable between and 1. Since , then . By considering (B.3) for all and , relation (B.2) is now a consequence of Proposition 3.1, (), and Slutsky’s theorem.
Appendix C Section 4.1
C.1. General results
For and , consider the jointly Gaussian vector
[TABLE]
[TABLE]
(see (A.29) for the definition of ). Let
[TABLE]
be the associated covariance matrix and its matrix spectral decomposition. The following lemma provides the finite-sample and asymptotic properties of the covariance structure of wavelet coefficients, both from and . Note that such covariance structure does not in general correspond to a multivariate stationary stochastic process when multiple octaves are considered.
Lemma C.1**.**
For , , and , the following statements hold.
Consider (see (3.12)). For and , there is a continuous spectral density such that we can write
[TABLE]
Moreover,
[TABLE]
for a constant that does not depend on .
Let be as in (3.14), , and fix any . Then,
[TABLE]
where and are defined in expressions (2.6) and (3.3), respectively.
Let be as in (C.1), and let , , be the eigenvalues of the covariance matrix (see (A.29) and (C.2)). Then, for some ,
[TABLE]
Proof.
We will use the shorthand and throughout the proof.
We first show . By making the change of variable in (A.1), we obtain (C.3) with
[TABLE]
Moreover, by (3.1)–(3.3), is continuous, and hence, by (2.5), (3.3), (3.4) and the dominated convergence theorem, the periodic function in (C.7) is also continuous on , and hence, bounded. This shows (C.3) and (C.4).
We now turn to . It suffices to show that, for any ,
[TABLE]
since for small enough , conditions (2.7) and (3.3) imply that
[TABLE]
In fact, by (A.1), (3.14) and a change of variable ,
[TABLE]
[TABLE]
[TABLE]
for any . For large enough , by (2.6), (2.7) and (3.1)–(3.3), the absolute value of the first integral in the sum (C.9) can be bounded by
[TABLE]
On the other hand, by (3.2) the absolute value of the second term in the sum (C.9) can be bounded by
[TABLE]
Expressions (C.10) and (C.11) yield (C.5).
We turn to . For notational simplicity, consider only two octaves , whence . Then, from (A.29),
[TABLE]
Fix . For notational simplicity, divide the summation range into the subranges
[TABLE]
[TABLE]
Define the octave and index functions
[TABLE]
respectively, which reflects the order of appearance of different octave and index values in the vector (C.1). By (A.1),
[TABLE]
[TABLE]
By expanding the double summation in (C.12) into each pair in the Cartesian product , we end up with 16 double summation terms under the sign of the integral, i.e., 8 pairs of conjugates. To the cross terms, namely, those involving distinct summation ranges in the index , we can apply the elementary inequality . One such pair, associated with the summation regions and , is
[TABLE]
[TABLE]
[TABLE]
since , and analogous bounds hold for the remaining terms. Therefore, (C.12) is bounded by
[TABLE]
[TABLE]
[TABLE]
[TABLE]
By the change of variable in each integral, breaking them up (in ) into subregions of length and using the periodicity of Fourier sums, we obtain
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
for some . The inequality (C.13) is a consequence of (C.3) (from Lemma C.1, ) and Proposition 3.1, , as applied to . Thus, the claim (C.6) holds.
∎
Remark C.1**.**
For and , let and be as in (A.28) and (C.1), respectively. Then,
[TABLE]
[TABLE]
[TABLE]
for some permutation matrix . Moreover, since is a real symmetric matrix, by Lemma C.1, and , we obtain the bound
[TABLE]
for the maximum eigenvalue of , where is the matrix Euclidean norm.
For any , consider the Gaussian vector as in (C.1) but with only one octave . It will be convenient to reexpress the sum as in (3.10) (with and ) based on a quadratic form. Define the permutation matrix
[TABLE]
Consider the spectral decomposition (C.2). Then,
[TABLE]
[TABLE]
where . Let
[TABLE]
which is a real symmetric matrix. Write out its spectral decomposition
[TABLE]
Expression (C.16) becomes
[TABLE]
[TABLE]
where and are the indices for which is nonnegative or negative, respectively. In particular, note that
[TABLE]
where
[TABLE]
and , , are as in expression (C.19).
The following lemma is an exponential inequality for -like distributions, and corresponds to Lemma 1 in (?). It will be used in the ensuing Lemma C.3 to establish a bound on the rate of convergence to zero of the probability , where (under (.2)).
Lemma C.2**.**
(?) Let and , not all zero. Let and be the Euclidean square and sup norms of the vector . Also, define the random variable . Then, for every ,
[TABLE]
[TABLE]
Lemma C.3**.**
Let be as in (C.19), and fix . Then, for any ,
[TABLE]
for large enough .
Proof.
Expression (C.19) implies that
[TABLE]
For notational simplicity, let and , which are zero mean random variables. Now fix . By the independence of and ,
[TABLE]
where is the density function of . The first integral in (C.24) is bounded from above by . Let
[TABLE]
Then, by (C.42) below (under (C.27)),
[TABLE]
for a constant . Moreover, since , then
[TABLE]
[TABLE]
[TABLE]
In (C.26), does not depend on and the last two inequalities are a consequence of (C.22) and (C.25) below, respectively, where
[TABLE]
stems from (2.9). On the other hand, for any and large enough the second integral in (C.24) is bounded from above by
[TABLE]
Suppose, without loss of generality, that
[TABLE]
(otherwise, for such that ). Therefore, by (C.25),
[TABLE]
for any , as . Consequently, for large enough , (C.28) is bounded by
[TABLE]
[TABLE]
where the last inequality follows from (C.21). From (C.26) and (C.29), since does not depend on , we obtain (C.23).
∎
C.2. Asymptotic covariances
We will establish Theorem 4.1 at the end of this section, after proving Lemmas C.4–C.7. The latter establish the asymptotic behavior of the first four moments of the -valued random variables , . So, consider the function
[TABLE]
where the truncation works as a regularization of the log function around the origin. In the event , we interpret that
[TABLE]
Throughout this section, we will make use of the Isserlis theorem (e.g., (?)). For a zero mean, Gaussian random vector , the theorem states that
[TABLE]
The notation stands for adding over all possible -fold products of pairs , where the indices partition the set .
The following lemma shows that the high order centered cross moments of , are negligible with respect to their low order counterparts.
Lemma C.4**.**
Let , , and fix . Then, as ,
[TABLE]
[TABLE]
Proof.
We first show that
[TABLE]
, . We will only establish (C.33) for and , since the remaining cases can be tackled by a similar argument.
The left-hand side of (C.33) can be rewritten as
[TABLE]
[TABLE]
Starting from the assumption (.2), for notational simplicity denote the generic terms under the summation sign in (C.34) as , , , , , . Then,
[TABLE]
[TABLE]
By applying the Isserlis relation (C.31) to the six-fold and four-fold products in (C.35), the latter expression becomes
[TABLE]
[TABLE]
[TABLE]
[TABLE]
The asymptotic behavior of the summation of each term in the seven-fold sum (C.36) can be established in the same way, so we only study the first one. Up to a constant, the latter is asymptotically equivalent to
[TABLE]
[TABLE]
[TABLE]
[TABLE]
However, the summation in in expression (C.37) is bounded by
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
where does not depend on . To justify the last inequality, we only look at the second term in (C.38), since the remaining terms can be analyzed in a similar way. It suffices to proceed as in (?), in particular around expression (B.31) in the latter reference. Indeed, starting from (3.2), suppose without loss of generality that
[TABLE]
For , decompose
[TABLE]
[TABLE]
[TABLE]
The first sum term in (C.39) only contains a finite number of terms, where such number does not depend on or . Moreover,
[TABLE]
where does not depend , or . The latter statement follows from the fact that is, up to a change of sign, the wavelet variance of a fBm and from (A.24). Moreover, by Proposition II.2 in (?) and again by (A.24), the second sum term in (C.39) is bounded by . This shows that (C.39) is bounded by a constant not depending on , or . Hence, (C.38) holds.
Consequently, up to a constant the absolute value of (C.37) is bounded from above by
[TABLE]
where we used the fact that . This establishes (C.33).
So, rewrite
[TABLE]
[TABLE]
[TABLE]
[TABLE]
The asymptotic behavior of every term on the right-hand side of (C.40) can be established in the same way, so we only look at the first one. For any , the Cauchy-Schwarz inequality, Lemma C.3 and expression (C.33) yield
[TABLE]
[TABLE]
[TABLE]
where does not depend on . Moreover, under (3.6),
[TABLE]
Therefore, (C.32) holds.
∎
The following lemma expresses, up to a residual term, the cross-covariance (first cross-moment) between sample wavelet variances in terms of the functions (3.14).
Lemma C.5**.**
Let and be as in (3.14) and (4.1), respectively. For ,
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
as .
Proof.
As in the proof of Lemma C.4, we first drop the indicator functions on the left-hand side of (C.41) and investigate the limit. We will show that
[TABLE]
[TABLE]
[TABLE]
[TABLE]
In fact, the left-hand side of (C.46) can be written as
[TABLE]
[TABLE]
By the Isserlis relation (C.31), the first term in the argument of the sum (C.43) can be reexpressed as
[TABLE]
[TABLE]
[TABLE]
By Lemma C.1, , and (C.44), we can rewrite (C.43) as
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
By adding and subtracting the counterparts for each term, up to the factor the expression (C.45) can be written as
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
It remains to justify the order of the error term in (C.46). So, by adapting the proof of (A.16), we obtain
[TABLE]
[TABLE]
The same bound holds for
[TABLE]
and
[TABLE]
[TABLE]
By extending this analysis to the remaining terms of (C.46), we obtain analogous bounds and the error term o\Big{(}\frac{a^{2\max\{h_{13}+h_{24},h_{14}+h_{23}\}-2(h_{12}+h_{34})}}{n_{*}}\Big{)}, as claimed.
In view of (C.45) and (C.46), it suffices to show that the indicator functions on the left-hand side of (C.41) do not affect the approximation order. In fact,
[TABLE]
[TABLE]
[TABLE]
[TABLE]
Define
[TABLE]
For , by the Cauchy-Schwarz inequality, expression (C.33) (from Lemma C.4) and Lemma C.3, the first term on the right-hand side of (C.47) is bounded by
[TABLE]
[TABLE]
[TABLE]
[TABLE]
where the constant does not depend on and the last equality is a consequence of (3.6). Similar bounds hold for the remaining terms on the right-hand side of (C.47). Therefore, the expression (C.41) follows.
∎
The following lemma describes the decay rate of the first individual truncated moment of the wavelet variance.
Lemma C.6**.**
For any and ,
[TABLE]
[TABLE]
Proof.
Notice that
[TABLE]
Hence,
[TABLE]
However, by the Cauchy-Schwarz inequality,
[TABLE]
By expressions (C.50), (C.42) and Lemma C.3, the expression (C.49) follows.
∎
The following lemma establishes the decay of the covariances between truncated terms (C.30) and indicators involving wavelet variances, or between indicators only.
Lemma C.7**.**
For any and ,
[TABLE]
[TABLE]
[TABLE]
[TABLE]
and
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
In (C.52), does not depend on .
Proof.
We first show (C.51). By the Cauchy-Schwarz inequality, the left-hand side of (C.51) is bounded from above by
[TABLE]
Moreover, for and the octave , Lemma C.3 implies that
[TABLE]
[TABLE]
The same bound holds for the octave in (C.55). Thus, (C.51) holds.
To prove (C.52), note that is bounded by
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
This establishes (C.52).
To prove (C.53), note that is bounded by
[TABLE]
[TABLE]
[TABLE]
However, for some ,
[TABLE]
where, by (A.39), and
[TABLE]
Therefore, by the dominated convergence theorem for convergence in probability,
[TABLE]
This establishes (C.53).
To show (C.54), again by applying the Cauchy-Schwarz inequality, the left-hand side of (C.54) is bounded from above by
[TABLE]
Hence, the bound follows from (C.51) and (C.53).
∎
We are now in a position to establish Theorem 4.1.
Proof of Theorem 4.1: Fix and recall that . Then,
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
The last two equalities in (C.58) are a consequence of Lemma C.1, , as applied to and , and of the bound (C.54) (from Lemma C.7, ) under the condition (4.2).
Therefore, it suffices to show that the main term on the right-hand side of (C.58) is equal to the main term on the right-hand side of (4.3). By accounting for absolute values, the covariance term in the former can be broken up into a sum of four terms, namely,
[TABLE]
plus the remainder
[TABLE]
[TABLE]
[TABLE]
By the Cauchy-Schwarz inequality, the bounds (C.52) and (C.53) (from Lemma C.7, ) and condition (4.2), the absolute value of the second term in the sum (C.60) is bounded by
[TABLE]
[TABLE]
By a similar argument, the same bound holds for the remaining terms in the sum (C.60). Thus, it suffices to focus on (C.59). In the following developments, expressions involving individual sample wavelet variance terms will be expressed in terms of , but analogous expressions hold when substituting for .
Fix . For any given , write out the almost sure Taylor expansion
[TABLE]
[TABLE]
where . Then,
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
For , recast
[TABLE]
[TABLE]
[TABLE]
Therefore, we can rewrite the fourth term in (C.62) as
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
By Lemma C.4, the fourth term term in the sum is bounded by
[TABLE]
By the Cauchy-Schwarz inequality, Lemma C.4 and condition (4.2), the first term in the sum (C.64) is bounded by
[TABLE]
[TABLE]
[TABLE]
[TABLE]
[TABLE]
since . The second term in the sum (C.64) is bounded by
[TABLE]
[TABLE]
[TABLE]
[TABLE]
An analogous bound holds for the third term in the sum (C.64). Therefore, by (C.62), (C.65), (C.66), (C.67) and Lemma C.5, we conclude that (C.58) is equal to the right-hand side of (4.3), as claimed.
Remark C.2**.**
For , a.s. (see (C.19)). Then, the existence of the moment , , can be directly established by applying relation (96) in (?). Moreover, the analysis of moments in this section can be extended without the truncation based on the sequence (4.2).
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1[1] P. Abry, R. Baraniuk, P. Flandrin, R. Riedi, and D. Veitch. Multiscale nature of network traffic. IEEE Signal Processing Magazine , 19(3):28–46, 2002.
- 2[2] P. Abry and G. Didier. Wavelet estimation for operator fractional Brownian motion: supplementary material. To appear in Bernoulli , pages 1–15, 2015.
- 3[3] P. Abry and G. Didier. Wavelet estimation for operator fractional Brownian motion. To appear in Bernoulli , pages 1–30, 2016.
- 4[4] P. Abry and P. Flandrin. On the initialization of the discrete wavelet transform. IEEE Signal Processing Letters , 1(2):32–34, 1994.
- 5[5] P. Abry, P. Flandrin, M. Taqqu, and D. Veitch. Wavelets for the analysis, estimation and synthesis of scaling data. In Self-similar Network Traffic and Performance Evaluation , pages 39–88. Wiley, spring 2000.
- 6[6] P. Abry, S. G. Roux, H. Wendt, P. Messier, A. G. Klein, N. Tremblay, P. Borgnat, S. Jaffard, B. Vedel, J. Coddington, and L. Daffner. Multiscale anisotropic texture analysis and classification of photographic prints: art scholarship meets image processing algorithms. IEEE Signal Processing Magazine , 32(4):18–27, July 2015.
- 7[7] P. Abry and D. Veitch. Wavelet analysis of long-range dependent traffic. IEEE Transactions on Information Theory , 44(1):2–15, 1998.
- 8[8] P. Abry, H. Wendt, and S. Jaffard. When Van Gogh meets Mandelbrot: multifractal classification of painting’s texture. Signal Processing , 93(3):554–572, 2013.