the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 17 Jul 2026
Systematic deviation in wavelet-based covariance estimation and implications for eddy-covariance flux calculations
Abstract. Accurate estimation of variance and covariance is central to the analysis of geophysical time series and underpins key diagnostics such as turbulent fluxes in eddy-covariance (EC) measurements. Continuous wavelet transform (CWT) methods are widely used for this purpose due to their ability to represent non-stationary and multi-scale processes. However, the consistency of commonly used formulations with the theoretical properties of the CWT has not been fully clarified. In this work, we show that the widely used formulation of Torrence and Compo can lead to systematic deviations in covariance estimation, which can be traced to normalization choices that are not fully consistent with the CWT reconstruction framework. These deviations, of the order of 10–15 % under typical conditions, are observed across a range of signal types, including deterministic, white-noise, and long-range dependent processes. We derive an alternative formulation directly from the CWT reconstruction identity, ensuring consistency with wavelet energy conservation. This approach reduces bias and improves convergence properties. Application to EC data shows that the proposed formulation recovers classical covariance estimates under stationary conditions while providing a more consistent framework for the analysis of non-stationary signals. These results emphasize the importance of scale-consistent estimation of second-order moments for interpreting geophysical variability and turbulent fluxes, and provide a methodological basis for multi-scale analysis of environmental time series.
- Preprint
(984 KB) - Metadata XML
-
Supplement
(307 KB) - BibTeX
- EndNote
Status: final response (author comments only)
- RC1: 'Comment on egusphere-2026-2881', Yubin Li, 08 Aug 2026
-
RC2: 'Comment on egusphere-2026-2881', Anonymous Referee #2, 26 Aug 2026
The comment was uploaded in the form of a supplement: https://egusphere.copernicus.org/preprints/2026/egusphere-2026-2881/egusphere-2026-2881-RC2-supplement.pdf
-
RC3: 'Comment on egusphere-2026-2881', Anonymous Referee #3, 14 Sep 2026
General comments
The manuscript " Systematic deviation in wavelet-based covariance estimation and implications for eddy-covariance flux calculations" by Jonathan Bitton et al. (egusphere-2026-2881) presents a systematic deviation in commonly used formulation for wavelet-based covariance showing implications for the eddy covariance field and proposing a better formulation.
The manuscript is very clear. I praise the authors for their successful effort in making what could be a mathematically heavy subject into a digestible, direct and clear paper.
It is in Atmospheric Measurement Techniques’ scope and it brings clear contributions to the flux community interested in frontier methods.
For the manuscript sections, introduction and results are clear and sound. They all show good scientific effort.
On the methods some points are raised:
- Should you mention the actual wavelet implementation, for example ψ(t)=2π−3/4exp(−t2/2)exp(iω0t) ?
- It’s striking the resemblance between the classical formulation and the real-valued wavelets formulation. Could we summarise the methods around using an empirical (𝐶𝛿) or theoretical (𝐶𝜓/ln(2)) reconstruction factor? For example in TC’s 𝐶𝛿 for Morlet (6) is 0.776, while 𝐶𝜓 is 1.11953135173, so the ratio between TC and CWT should be (1/𝐶𝛿) / (2*ln(2)/𝐶𝜓), so (1.28865979381443298) / (1.238280963706522869254) = 1.040684490502956708?
- If this is true why the gap in Figure 1 and Figure 3 between CWT and TC seems wider? In both figures it seems to be bigger than 1.1 (and not 1.04 as the previous logic would sugest).
- Also fed from the previous comment, the 𝐶𝜓 reconstruction factor in section 2.1.2 seems to be a key piece and deserves more detail. In here or in appendix A, you should develop a bit more in detail how do you get C𝜓 ≈ 1.1988 for Morlet (6) and 𝐶𝜓 ≈ 2.1349 for DOG (2). You may want to include Morlet’s mother wavelet, for example ψ(t)=2π−3/4exp(−t2/2)exp(iω0t).
Some changes in the discussion section may help highlight the paper’s ideas. Section 4.3 is too long. The most direct discussion to the flux community you are targeting is the last but 4.3’s length gets in the way. I propose a few suggestions below. Section 4.4 could also be reworded for readability.
Minor comments
- In Figure 3, it's unclear the non-stationarity (b, d) role here. The interesting comparison is within EC plausible range, one would not use EC under non-stationarity. All is explained by the authors in the text (l. 240-247), but this paragraph adds little to the paper’s message.
- 265, can you be more concrete “long-range dependent geophysical processes”? Possibly add “(eg. ...)”.
- 272 where does this come from: “Under non-stationary conditions, covariance contributions become intermittent and strongly scale-dependent.”?
- 272-276. Do you mean non-stationary conditions and pink-noise are somewhat similar?
- 288-289. Replace “with the one-to-one relationship.” By “(Fig. 3)”
- 338-340. A too long concluding remark. Suggestion to split it “These results underline the importance of theoretically consistent formulations in multi-scale analysis. They indicate that wavelet-based estimates of geophysical variability and turbulent fluxes may benefit from a careful reassessment within this framework.”
- 329-331. The phrasing can be improved using something more down to earth like “help identify the temporal scales driving observed surface–atmosphere exchanges.” or by given a real-world example before the final sentence.
Suggestions on the form
I took the liberty to give a few suggestions to the authors on the rewording of some passages. They come from goodwill to strengthen the readability and so impact of the manuscript, they don't come to correct misconceptions.
Section 4.3
- 291-296. What do you want to say here? In the beginning it feels like you are going for “theoretically Morlet is better for non-stationarity”, but you can’t truly base that on the results. So you conclude “and broader analytic wavelets for improved sensitivity to scale-dependent and non-stationary processes” which sound a bit too vague.
Suggestion; “In contrast, Morlet detects departures from bulk covariance earlier because its broader support captures low-frequency and transitional variability. This is not a loss of consistency, but a difference in what each wavelet resolves. Therefore wavelet choice should match the objective: localized wavelets when reproducing bulk estimates and broader wavelets when resolving scale-dependent non-stationarity.” - 280-286. “The Morlet wavelet shows faster convergence and lower dispersion across signal types, consistent with its improved frequency localization and more uniform energy distribution across scales. Its complex coefficients may further reduce estimator variance through their two-degree-of-freedom X² distribution, compared with the one degree of freedom of real-valued wavelets (Ge, 2007; Torrence and Compo, 1998). This effect remains, however, secondary to spectral structure and sample size.”
- 297-306. The same as above. “A persistent DOG2 bias was observed in the white-noise experiments, likely due to its broad, overlapping frequency response and imperfect spectral partitioning (Percival and Walden, 2000). This causes uneven weighting of frequency bands and systematic bias in wavelet variance estimates. The more localized Morlet response provides a cleaner spectral partition and reduced bias. For pink-noise processes, this filter-induced bias becomes secondary to the limited degrees of freedom at large scales.”
Section 4.4
- L.308-311. “The proposed formulation provides a physically consistent framework for estimating variance and covariance across scales, where they are particularly suited to long, correlated, and multiscale geophysical time series (Aubinet et al., 2012; Seo et al., 2023; Stull, 1988).”
Citation: https://doi.org/10.5194/egusphere-2026-2881-RC3
Model code and software
Wavelet analysis python/matlab codes Jonathan Bitton https://doi.org/10.5281/zenodo.20187583
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 210 | 83 | 34 | 327 | 52 | 32 | 30 |
- HTML: 210
- PDF: 83
- XML: 34
- Total: 327
- Supplement: 52
- BibTeX: 32
- EndNote: 30
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
This study revisits wavelet-based variance/covariance estimation for eddy-covariance flux calculations. The authors identify that the Torrence and Compo formulation introduces 10–15% systematic bias due to normalisation issues, and propose a corrected CWT-derived estimator. The theoretical derivation is rigorous, synthetic validation is comprehensive, and the EC data application demonstrates practical relevance. The work addresses an important methodological gap in atmospheric flux measurements. Overall, the paper is well structured and clearly written. However, several minor issues require attention before publication.
1. The authors claim that the TC formulation exhibits a systematic bias of approximately 10–15%, whereas the proposed CWT-based estimator is essentially unbiased. For synthetic signals, this claim is well supported. However, in the EC application, the authors treat the classical time-averaged EC covariance as the reference benchmark (page 12: "Under stationary conditions... covariance reconstructed with the CWT-based formulation closely matches EC estimates... TC covariance exhibits larger deviations..." and "Under non-stationary conditions... the agreement decreases substantially"). The concern is that the classical EC estimator itself relies on Reynolds decomposition and strict stationarity assumptions; it is not an absolute reference, particularly under non-stationary conditions. Using a compromised reference to evaluate wavelet methods in the non-stationary regime risks circularity. The observed discrepancy may reflect the breakdown of the bulk EC average's representativeness rather than a limitation of the wavelet approach. A clearer distinction should be drawn between (i) the theoretical truth (synthetic signals) and (ii) the EC reference value (a conventional, operationally defined estimate, not a gold standard). The non-stationary results would be better framed as a diagnostic of multi-scale covariance structure and the limitations of ergodicity assumptions, rather than as a "deviation from the true flux."
2. Synthetic experiments show Morlet outperforms DOG2 in convergence and frequency localisation (page 10, Figure 2), suggesting theoretical superiority. Yet in the EC application under stationary conditions (Figures 3a, 3c), DOG2 shows closer 1:1 agreement than Morlet—an apparent reversal not adequately explained by invoking DOG2's "stronger temporal localisation." Possible hypotheses: (i) DOG2's broad frequency response may compensate for the 2/3 inertial subrange slope; (ii) Morlet's wider temporal support may increase sensitivity to residual non-stationarities in "stationary" EC periods. Further investigation (e.g., residual spectral decomposition) would clarify whether this advantage is physically meaningful.
3. In the TC formulation (page 3), the reconstruction factor is denoted as Cδ, while in the proposed formulation (page 4), the admissibility constant is denoted as Cψ. Although the authors indicate that these are conceptually distinct—the former being an empirical discretisation factor and the latter the theoretical admissibility constant—their numerical values differ substantially for a given wavelet. This notational proximity is likely to cause confusion. An explicit statement is needed, where Cδ is first introduced, clarifying that it is not equivalent to Cψ, but rather an empirical normalisation scalar arising from the discretisation scheme of Torrence and Compo (1998). A small table listing the numerical values of both constants for the two wavelets used would further aid clarity.
4. On page 6, the authors state: "Stationary and non-stationary conditions were identified using the combined steady-state and integral turbulence characteristics test (SSITC), with SSITC=0 indicating stationary conditions and SSITC=2 indicating non-stationary conditions (Foken et al., 2004)." However, the threshold values or criteria used to assign these classes are not described. It would be helpful to specify, either here or in a supplementary section, the exact numerical thresholds employed for the steady-state test and the integral turbulence characteristics test, as different implementations of the SSITC use different tolerance levels.
5. The reduction from 44 (DOG2) to 42 (Morlet) non-stationary periods (page 6) is modest given Morlet requires 5× longer extensions. It should be clarified whether this reflects data availability or other filtering criteria, and the edge-effect exclusion criteria for each wavelet should be briefly explained.
6. The term "pseudo-ensemble" (page 6) is non-standard. An ensemble typically means multiple realisations from the same parameter setting; here it appears to mean discrete parameter sampling (frequency, amplitude, variance) instead. For stochastic signals, a single realisation itself contains sampling variability, so the "pseudo-ensemble" mixes parametric and stochastic variability. Clarification is needed on whether each parameter combination corresponds to one or multiple realisations, and justification of the statistical representativeness.
7. In Figure 3, a legend is provided only in subpanel (d); subpanels (a)–(c) lack one. This hinders interpretation when comparing across conditions or wavelet types. A consistent legend should be added to all subpanels, or a single comprehensive legend placed outside the panel array.