Evidence for a Possible Weak 1–2 Year Quasi-Periodic Variability Band in Enceladus Plume Activities
Updated: Sep 1
Jicheng Gu*1 ； Yi Gu*1； Jeff Yu Gu*2
*1 GU Institute of Earthquake Prediction, Boston, U.S.A.
*2 Dept. of Physics, University of Alberta, Alberta, Canada
Preprint: This is a preprint (pre-submission version) and has not yet undergone peer review. The official version published in the journal will prevail if it is formally accepted and published.
Abstract
Enceladus exhibits persistent cryovolcanic plume activity driven by tidal forcing and internal
energy dissipation. Previous studies have primarily focused on short-timescale variability,
including orbital tidal modulation, plume brightness fluctuations, and local jet dynamics.
However, the possibility of longer-timescale organized variability remains insufficiently
explored.
In this exploratory note, we construct a composite plume activity index (X(t)) from normalized
proxies of plume mass flux, thermal output, and optical brightness. Using Fourier analysis,
Lomb–Scargle periodograms, bootstrap resampling, permutation testing, AR(1) red-noise
surrogate analysis, weight-robustness scans, and sliding-window spectral analysis, we investigate whether a persistent low-frequency structure exists within the available observational record.
Across multiple analysis methods and index formulations, we consistently observe enhanced
spectral power within a broad 1–2 year band, with local maxima near approximately 1.3–1.8
years and weaker secondary structure near ~2.1 and ~3.6 years. However, the peak location is
unstable under bootstrap resampling, and the observed spectral power does not significantly
exceed simple red-noise surrogates.
We therefore do not claim detection of a statistically significant coherent periodicity. Instead, the present results motivate a weaker but potentially important hypothesis: Enceladus plume activity may contain a weak, drifting, quasi-periodic low-frequency component superposed on persistent emission and red-noise-like variability.
1. Introduction
Enceladus, an icy moon of Saturn, exhibits continuous cryovolcanic plume activity sourced from fractures near its south polar terrain. Since the Cassini mission, plume variability has been extensively studied through measurements of thermal emission, particle flux, optical brightness, and in situ compositional sampling.Previous studies have demonstrated that plume activity varies on multiple timescales. These include orbital tidal modulation associated with Enceladus’ eccentric orbit, transient changes in individual jets, and possible multi-year variability related to internal structural or tidal processes.
However, most investigations have focused on single observables or local plume behavior.
The present work asks a different question: whether the plume system, viewed as a coupled
macroscopic energy-release process, exhibits any organized low-frequency variability detectable through a composite state variable.
To explore this possibility, we construct a first-pass composite activity index (X(t)) inspired by
system-level proxy approaches used in geophysical and climate sciences. The purpose of this
study is not to claim the discovery of a strict periodicity, but rather to examine whether the
available data contain evidence for a weak and potentially nonstationary quasi-periodic structure.
Existing studies do not describe Enceladus’ plume activity as strictly featureless or purely
random. Instead, they indicate persistent emission with variability on diurnal, monthly, multi-
year, and possibly decadal timescales. The present note therefore investigates whether a weak 1–2 year variability band may also be present.
2. Composite Activity Index
2.1 Motivation
Single plume observables are often strongly affected by local turbulence, viewing geometry,
instrumental effects, and short-timescale fluctuations. A composite index may better capture the large-scale dynamical state of the plume system.
Following this idea, we define a normalized composite activity index:
[ X(t)=0.5_x(t)+0.4_t(t)+0.1_t(t) ]
where:
• (_x(t)): normalized plume mass-flux proxy;
• (_t(t)): normalized thermal-output proxy;
• (_t(t)): normalized optical-brightness proxy.
All quantities are scaled to the interval ([0,1]).
The weighting scheme is heuristic and intended only as a first-pass exploratory construction. In particular, the brightness term is assigned a relatively small weight because optical brightness is influenced by observational geometry, scattering effects, and instrumental conditions in addition to intrinsic plume activity.
2.2 Brightness Term and Observational Effects
The observed plume brightness depends on several factors, including:• phase angle and illumination geometry;
• line-of-sight integration effects;
• instrumental response and calibration;
• particle-size distribution and scattering properties;
• background scattering and stray light.
Therefore, optical brightness is not treated as a direct measure of plume intensity, but rather as a weak auxiliary observational proxy.
Ideally, the brightness component would be corrected through phase and geometric calibration:
[ _{}(t)= ]
where:
• (P()): phase function;
• (G(t)): geometric correction factor;
• (B_{}): background signal.
Because fully synchronized corrected brightness data are not yet available, the present work uses a simplified first-pass approximation.
3. Data and Methods
Figure 1. First-pass Enceladus activity index and FFT comparison. (a) Time series of normalized plume activity proxies and composite index components. (b) FFT spectrum on period axis. (c) FFT spectrum on frequency axis. (d) Dominant spectral peaks
identified in the composite activity index.
3.1 Data
The present analysis uses publicly available plume-related observations derived from Cassini-era datasets and associated plume activity proxies.
The observational record is sparse and unevenly sampled, spanning approximately one decade.
Because the effective number of independent cycles within a 1–2 year band is limited, all
conclusions should be regarded as preliminary.
3.2 Fourier Analysis
As an initial exploratory step, Fourier transforms were applied to interpolated versions of the
activity index in order to identify candidate low-frequency spectral structure.
3.3 Lomb–Scargle Periodograms
Because the data are irregularly sampled, Lomb–Scargle periodograms were used as the primary spectral method. This approach avoids artifacts associated with direct FFT analysis of unevenly sampled data.
3.4 Bootstrap and Permutation Tests
Bootstrap resampling was used to evaluate the stability of peak locations within the candidate 1–2 year band.
Permutation testing was then applied by randomizing temporal ordering while preserving the
value distribution, providing a nonparametric null model for comparison.
3.5 AR(1) Red-Noise Surrogates
To evaluate whether the observed low-frequency structure exceeds expectations from correlated stochastic variability, AR(1) red-noise surrogate series were generated and analyzed using the same Lomb–Scargle procedure.
3.6 Weight Robustness
The robustness of the candidate low-frequency structure was tested by scanning large numbers of alternative weighting combinations:
[ X_w(t)=w_1_x+w_2_t+w_3_t ]
subject to:
[ w_1+w_2+w_3=1 ]
3.7 Sliding-Window Spectral Analysis
Sliding-window Lomb–Scargle analysis was performed to evaluate possible temporal drift of
local peak frequencies.
4. Results
4.1 Spectral Structure
Both FFT and Lomb–Scargle analyses reveal enhanced spectral power within a broad 1–2 yearband. Local maxima are typically found near:
• ~1.3 years;
• ~1.5 years;
• ~1.7–1.8 years.
Additional weaker components occasionally appear near:
• ~2.1 years;
• ~3.5–3.6 years.
Importantly, these structures appear across multiple formulations of the activity index and across different spectral methods（Figure 2）.
Figure 2. Lomb–Scargle spectrum and bootstrap stability analysis of the candidate 1–2
year variability band. (a) Observed Lomb–Scargle spectrum. (b) Bootstrap distribution of
peak periods. (c) Bootstrap distribution of peak powers. (d) Summary statistics and
confidence intervals.
4.2 Bootstrap Stability
Bootstrap resampling demonstrates that the precise peak location within the 1–2 year band is
unstable. Rather than converging to a narrow dominant period, the bootstrap distribution spreads broadly across approximately 1.1–1.9 years.
This behavior is more consistent with a broad quasi-periodic band than with a sharply defined
coherent oscillator(Figure 3).
Figure 3. Cross-method consistency checks for the candidate 1–2 year variability band. (a) Sparse unevenly sampled activity index data. (b) FFT after interpolation. (c) Lomb–Scargle spectrum computed directly from uneven data. (d) Permutation-test null distribution for the 1–2 year band.
4.3 Permutation Testing
Permutation testing indicates that randomized datasets frequently produce spectral peaks of
comparable strength within the same low-frequency range.
The observed 1–2 year peak therefore cannot presently be distinguished from stochastic low-
frequency variability at conventional significance thresholds.
4.4 AR(1) Red-Noise Analysis
AR(1) surrogate testing similarly shows that correlated red-noise processes can reproduce
comparable low-frequency spectral enhancement.
The observed peak power does not significantly exceed the red-noise null distribution.
4.5 Weight Robustness
Despite the statistical limitations above, the presence of enhanced low-frequency power within the 1–2 year band remains relatively robust under large variations of the weighting scheme.
This suggests that the observed low-frequency structure is not solely an artifact of the original
heuristic weights (Figure 4).
Figure 4. Extended robustness tests for the candidate low-frequency variability band. (a) Weight-sensitivity scan across randomized composite-index coefficients. (b) AR (1) red-noise surrogate distribution compared with the observed Lomb–Scargle peak power. (c) Sliding-window Lomb–Scargle peak drift. (d) Observed Lomb–Scargle spectrum withhighlighted 1–2 year band.
4.6 Sliding-Window Drift
Sliding-window Lomb–Scargle analysis reveals gradual drift of local peak positions over time.
Rather than remaining fixed at a single frequency, the candidate peak shifts between
approximately 1.3 and 1.8 years depending on the temporal window analyzed.
This behavior is consistent with weak nonstationary or quasi-periodic variability.
5. Discussion
5. 1 Theoretical Framework: Generalization to Fluid-Dissipative
Systems via Unified Metric Dynamics (UMD)
To resolve the physical origin of the weak, nonstationary 1–2 year quasi-periodic variability
band observed in Enceladus’ composite plume activity index X(t), we transcend purely
descriptive statistical methodologies and reformulate the problem through the lens of Nonlinear Statistical Physics and Continuum Mechanics. We introduce a generalized fluid-dissipative formulation of the Unified Metric Dynamics (UMD) framework—originally established for deterministic tidal phase-locking in terrestrial and lunar seismic networks (Gu, 2026).
While celestial seismicity governed by stress intensity factors and subcritical fatigue crack
propagation manifests as a second-order inertial-fracture manifold, the internal transport of mass and thermodynamic energy within Enceladus’ sub-surface hydrothermal environment represents a first-order fluid-dissipative system. We conceptualize the active cryovolcanic chamber beneath the South Polar Terrain (SPT) as an open, non-equilibrium dissipative structure. The temporal evolution of the collective sub-surface excess pressure is captured by a first-order non-linear time-delayed differential equation governing the dimensionless state variable 𝛹(𝜏) :
Here, rather than serving as mere scalar fitting parameters, the matrices (or tensors) 𝐀,𝐁 and 𝐂 act as dynamic coupling metric operators that project different statistical fields onto the phase-space manifold of the system.
• The Order Metric (𝐀) dictates the geometric transfer efficiency of directional
macroscopic energy injected by long-range gravitational fields into the localized
structural fluid chamber.
• The Intrinsic Vitality Metric (𝐁) defines the non-linear internal scaling constraint that
bounds the local degrees of freedom, determining the deterministic chaotic dissipation
rate of the fluid convection network.
• The Environmental Noise Metric (𝐂) represents the stochastic fluctuations of the local
ambient stress-temperature field, filtering outer random perturbations into the system’s
metric space.
Together, this triad (𝐀, 𝐁, 𝐂) defines the mathematical metric tensor that governs the
irreproducible trajectory of the non-equilibrium dissipative state variable. 𝑇orbit
where the intrinsic clock of the system is rigidly metricized by the orbital period of Enceladus
𝑡 around Saturn, normalizing the time variable as 𝜏 =(𝑇orbit = 1.37 days), thereby enforcing the base driving frequency. 𝛺orbit ≡ 1.0.
In this unified formulation:
𝐀 ⋅ 𝛷orbit(𝜏) represents the Deterministic Quasi-Periodic Skeleton (Theorem III), serving as an external gravitational pump driven by orbital eccentricity variations.
−𝐁 ⋅ 𝛹3(𝜏 − 𝛿fluid) denotes the Intrinsic Low-Dimensional Chaos Operators, parameterizing the delayed fluid-thermal convective feedback loop within the porous
core-ocean boundary, where 𝛿fluid is the characteristic relaxation time lag of internal hydrodynamic transport.
−𝜅𝛹 parameterizes the linear bulk dissipation representing continuous plume leakage and fault closure.
𝐂 ⋅ 𝜉(𝜏) injects a weak stochastic environmental noise tensor acting upon the fluid manifold.
This delayed negative feedback encapsulation aligns profoundly with recent high-resolution
planetary hydrodynamic models. Specifically, the pioneering multi-dimensional ocean
circulation simulations conducted by Kang and Flierl (2020) and Kang (2024) demonstrate that
Enceladus’ sub-surface ocean is characterized by intense hydrothermal plume convection and
robust salinity stratification dynamically driven by core heating. Within the geometry of the
kilometer-scale porous fault conduits beneath the South Polar Terrain (SPT), this complex fluid-thermal circulation naturally induces a significant macroscopic mass and energy transport
relaxation time lag.
By generalizing these complex, thousands-of-lines hydrodynamic grid constraints into a single, mathematically closed fixed-lag metric operator 𝛿fluid, the UMD framework successfully reduces the micro-scale fluid variables into a macroscopic state trajectory. Thus, the micro-convective instabilities identified by Kang et al. (2022) serve as the explicit physical foundation for the low-dimensional chaotic modulation metric 𝐁 utilized in our unified formulation.
Figure 5. Numerical simulation of the first-order fluid-dissipative UMD model
tailored for Enceladus. >
(Top Panel) The temporal evolution of the dimensionless state variable
𝛹(𝜏)representing the sub-surface chamber excess pressure over extended orbital cycles.
Driven by the external tidal skeleton, the system displays persistent localized envelope
fluctuations under non-linear fluid-thermal modulation, bounded safely beneath the
critical thermodynamic threshold (𝛹 → 1.0) without experiencing numerical divergence.
(Bottom Panel) The corresponding dimensionless power spectrum derived via Fast
Fourier Transform (FFT). The sharp, prominent order pillar at the primary tidal orbital
forcing frequency coexists with a continuous, broadband chaotic modulation field and a
structured background, successfully reproducing the physical mechanism underlying the
empirical nonstationary 1–2 year variability found in empirical plume index observations.
Numerical integration of the UMD state equation utilizing a (fixed-lag kernel) under standard
Euler-Maruyama constraints yields a profoundly rich dynamical topography (Fig. 5). In thetemporal domain (Fig. 5, top), the chamber pressure𝛹(𝜏) does not settle into a trivial static
equilibrium or disintegrate into white noise; instead, it tracks the orbital gravitational forcing
skeleton while experiencing strong non-linear modulations. Upon reaching the critical structural failure threshold (𝛹 → 1.0), a non-linear phase transition is triggered, successfully reproducing the intermittent, multi-scale burst behavior characteristic of the hydrothermal plumes.
Crucially, the resulting power spectrum (Fig. 5, bottom) exhibits a dual structure: a resilient,
sharp energy peak locked to the external orbital skeleton, flanked by a structured, continuous
broadband chaotic base. This specific topological configuration provides a direct physical
explanation for the empirical instability of the 1–2 year peak under bootstrap resampling: the
observed variability is not a static harmonic resonance, but rather a nonstationary emergent phenomenon born from the coupling of low-dimensional chaos and delayed fluid feedback. This successful mapping from second-order rigid fracture mechanics (the Moon) to
first-order dissipative hydrodynamics (Enceladus) firmly establishes the structural universality of the Five Fundamental Theorems of Predictability (Gu, 2026), proving that the Quasi-Periodic Deterministic Law governs the macrodynamics of complex planetary systems regardless of their specific chemical or material compositions.
It is worth noting that while classical planetary climate models focus primarily on steady-state or statically averaged thermal-orbital equilibrium profiles (e.g., Kang & Flierl, 2020; Kang et al.,
2022), they often bypass the long-term nonstationary time-series variations inherent in the
empirical observations. Our UMD integration bridge this exact theoretical vacuum.
While individual localized convective plumes resolve symmetric space-time profiles in steady-
state fluids, their multi-scale coupling under a rigid orbital skeleton yields an emergent
broadband chaotic spectrum. As illustrated in Figure 5, the system naturally generates
nonstationary low-frequency envelope fluctuations near the 1.3–1.8 year band, explicitly
explaining why these periodic signatures appear diffuse and unstable under empirical statistical re-sampling.
Mathematically, the transcendence of 𝐀，𝐁 and 𝐂 from linear algebraic coefficients into
dynamic metric operators directly establishes a profound epistemological link with the geometric formulation of Einstein's General Relativity (GR). In GR, the metric tensor 𝑔𝜇𝜈 dictates the curvature of spacetime and governs the geodesic trajectories of matter. Symmetrically, the UMD triad (𝐀, 𝐁, 𝐂) defines the non-Euclidean geometric metric of the system's dissipative phase space, framing the observed plume or seismic activity variations as deterministic trajectories along a complex thermodynamic manifold.
Crucially, when planetary bodies encounter external gravitational wave perturbations or high-
order tidal multipole tensors—which physically manifest as spacetime metric ripples ℎ𝜇𝜈—these perturbations are directly absorbed and converted through the Order Metric operator 𝐀. For a high-Q factor celestial manifold like the Moon or the rigid crustal interfaces of icy moons, this framework offers a deterministic mechanism wherein subtle relativistic spacetime metric fluctuations are resonantly amplified into macroscopic, phase-locked observable indexes.
Therefore, UMD provides a foundational bridge linking localized macroscopic continuum
mechanics with universal general relativistic field theories.
5.2 Notice
The present results do not support the existence of a statistically secure coherent periodicity in Enceladus plume activity.
However, the analysis consistently points toward a weaker and potentially more physically
realistic picture: the plume system may exhibit organized low-frequency variability superposed on persistent stochastic emission.
This interpretation differs substantially from a strict deterministic-clock model. Instead, the
observed behavior resembles weakly coherent variability commonly found in complex nonlinear natural systems, including climatic oscillations, geophysical release systems, and tidally modulated environments.
Several aspects of the present results are notable:
The low-frequency enhancement appears across multiple observables and composite formulations.
The signal persists across both FFT-based and Lomb–Scargle analyses, suggesting that it is not purely an interpolation artifact.
The peak position drifts with time and under resampling, favoring a quasi-periodic
interpretation.
The low-frequency structure does not clearly exceed red-noise expectations, emphasizing the limitations imposed by the short observational baseline.
Taken together, these results motivate a testable hypothesis rather than a definitive claim:
Enceladus plume activity may contain a weak, drifting, quasi-periodic low-frequency component on approximately 1–2 year timescales.
Future progress will likely require:
onger-duration monitoring;
improved phase-corrected brightness measurements;
additional plume-state observables;
wavelet or time-frequency analysis;
physically motivated tidal and fracture-flow models.
6. Conclusion
We constructed a first-pass composite plume activity index for Enceladus using normalized
proxies of plume mass flux, thermal output, and optical brightness.Across multiple independent spectral approaches, the analysis consistently reveals enhanced power within a broad 1–2 year variability band.
However, bootstrap, permutation, and AR(1) surrogate tests indicate that the evidence is
insufficient to establish a statistically significant coherent periodicity.
The present work therefore does not claim the discovery of a fixed Enceladus plume cycle.
Instead, it proposes a more cautious interpretation:
Enceladus plume activity may exhibit weak, nonstationary, quasi-periodic low-frequency
organization superposed on persistent stochastic variability.
Although preliminary, this framework provides a potentially useful direction for future
investigation of coupled plume dynamics in icy worlds.
References
Gu, J. (2026). Quasi-Periodic Deterministic Law of Seismic Systems: Five Fundamental
Theorems Supporting Earthquake Predictability. Zenodo.
Gu, J., Gu, Y., & Gu, J. Y. (2026). Detection of a 13.606-Day Periodicity in Lunar Seismicity
and Its Implications for Deterministic Tidal Forcing. Zenodo.
Gu, J. et al. (2026). Moon as a Natural Gravitational-Wave Antenna: Phase-Locked Lunar
Seismicity and the Possibility of Resonant Detection of Weak Spacetime Perturbations. Zenodo.
Khawaja et al. (2025). Detection of organic compounds in freshly ejected ice grains from
Enceladus’s ocean. Nature Astronomy
Souček et al. (2024). Variations in plume activity reveal the dynamics of water-filled faults on
Enceladus. Nature Communications.
Kang, W. (2024). Hydrothermal activity and ocean salinity stratification on Enceladus. Science
Advances, 10(14), eadj7785. https://doi.org/10.1126/sciadv.adj7785
Ershova et al. (2024). Dust-plume modeling based on Cassini in situ observations.
Villanueva et al. (2023). Webb observations of large-scale water vapor plumes from Enceladus.
Kang, W., Mittal, T., Bire, S., Campin, J. M., & Marshall, J. (2022). Across-latitude asymmetry
of Enceladus' ice shell explained by ocean circulation and sub-ice heating. Icarus, 374, 114808.
https://doi.org/10.1016/j.icarus.2021.114808Kang, W., & Flierl, G. (2020). Spontaneous formation of geysers at only one pole on Enceladus' ice shell. Proceedings of the National Academy of Sciences (PNAS), 117(26), 14764-14769.
Ingersoll, Ewald & Trumbo (2020). Time variability of the Enceladus plumes. Icarus. 344,
113735.
Preprint: This is a preprint (pre-submission version) and has not yet undergone peer review. The official version published in the journal will prevail if it is formally accepted and published.
Comments