The El Niño Southern Oscillation Recharge Oscillator with the Stochastic Forcing of Long-Term Memory

: The influence of the fast-varying variables that have a long-term memory on the El Niño Southern Oscillation (ENSO) is investigated by adding a fractional Ornstein–Uhlenbeck (FOU) process stochastic noise on the simple recharge oscillator (RO) model. The FOU process noise converges to zero very slowly with a negative power law. The corresponding non-zero ensemble mean during the integration period can exert a pronounced influence on the ensemble-mean dynamics of the RO model. The state-dependent noise, also called the multiplicative noise, can present its influence by reducing the relaxation coefficient and by introducing periodic external forcing. The decreasing relaxation coefficient can enhance the oscillation amplitude and shorten the oscillation period. The forced frequency is close to the natural frequency. The two mechanisms together can further amplify the amplitude and shorten the period, compared with the state-independent noise or additive noise, which only exhibits its influence by introducing non-periodic external forcing. These two mechanisms explicitly elucidate the influence of the stochastic forcing on the ensemble-mean dynamics of the RO model. It provides comprehensive knowledge to better understand the interaction between the fast-varying stochastic forcing and the slow-varying deterministic system and deserves further investigation.


Introduction
The El Ni ∼ no Southern Oscillation (ENSO) phenomenon is generated by a complex atmosphere-ocean coupling mechanism and is characterized by SST anomalies (SSTAs) that fluctuate significantly on a large scale in the tropical Pacific [1].ENSO has been considered to exert significant and profound influences on changes in the regional weather and even global climate patterns [2].Therefore, the fundamental dynamics of ENSO have been extensively investigated.Generally speaking, the ENSO cycle is controlled by positive feedback and negative feedback.This feedback determines the SSTA variations that can affect the trade wind, and hence the oceanic motion, which can further affect SSTA variation, in turn.The Bjerknes positive feedback emphasized that a warm SSTA in the eastern Pacific can weaken the trade wind, then reduce the oceanic upwelling, and then further amplify the warm SSTA there, and vice versa.There are two classic negative mechanisms.For example, a delayed oscillator model for ENSO is established to highlight the oceanic wave transit effect by introducing negative, delayed feedback [3,4].A recharge oscillator (RO) model is constructed to emphasize the recharge-discharge of the equatorial heat content as the phase-transition mechanism [5].Recently, a spatiotemporal oscillator model was proposed to characterize the complex temporal and spatial variations in ENSO [6].These studies considered ENSO as a slow deterministic oscillation of the atmosphere-ocean interaction.However, in fast atmospheric synoptic-intraseasonal variations, such as the Madden-Julian Oscillation (MJO) [7][8][9][10][11], the westerly wind bursts (WWBs) [12][13][14][15][16] may also have a significant influence on the variability and prediction of ENSO.This synoptic-intraseasonal variability, the source of the noise in the tropics, can be regarded as a stochastic forcing of ENSO [17].
Specific to the RO framework, the random wind stress forcing and random heating forcing can be added into SSTA and the thermocline depth anomaly (TDA) dynamics [18][19][20].For mathematical convenience, these random forcings are represented by the Langevin equation with a constant damping coefficient, that is, the Ornstein-Uhlenbeck (OU) process.The results suggest that the stochastic forcing can change the growth rate to enhance ENSO activity in the RO model.Moreover, the stochastic forcing can increase the ensemble spread and thus may reduce the effectiveness of the prediction of ENSO.However, the OU process is a Gaussian process with zero ensemble mean and exponential auto-correlation function, which denotes a short-term memory.It might not be a better prototype of the real world, where the auto-correlation functions of variables might be power functions, denoting a long-term memory [21,22].
Considering the fact that the ensemble-mean dynamics of the RO model under stochastic forcing with a short-term memory have been intensively investigated [18][19][20], we may explore the performance of stochastic forcing with a long-term memory.As an analog, the fractional OU (FOU) process that has a long-term memory [23][24][25] is a natural choice.Then, the influence mechanism of the FOU process of stochastic forcing on the ensemble-mean dynamics of the RO model, and the differences from the previous OU process stochastic forcing results [18], become the key scientific questions.Therefore, this paper attempts to answer these questions by applying a simple RO model [26] under the stochastic forcing of the FOU process.The investigation will be beneficial in allowing an in-depth understanding of the interaction between the fast-varying stochastic forcing and the slow-varying deterministic RO model.
Following this introduction section, the simple RO model forced by the FOU process forcing is presented in Section 2. The analytic solution of the FOU process is also discussed.The ensemble-mean dynamics of the RO model forced by stochastic forcing are derived and analyzed in Section 3. Finally, a conclusion and discussion end the paper in Section 4.

The Model Framework
Previous studies [18][19][20] have added multiplicative noise forcing onto the simple conceptual RO model to describe the influence of the fast atmospheric fluctuations that associate with the MJO/WWB activity.This theoretical framework has the following form: where T and h are the sea surface temperature anomaly (SSTA) in the eastern Pacific and the thermocline depth anomaly (TDA), respectively.λ is the relaxation of SSTA toward climatology.ω denotes the thermocline upwelling process.σ is the noise amplitude.ξ(t) is a stochastic process.W(t) is a standard Brownian motion (or the Wiener process), and its increment dW is the Gaussian white noise.r > 0 is a constant damping coefficient.G = 1 + BT represents the state-dependent noise forcing due to the interannual modulation of WWB/MJO activity via SSTAs.Particularly, if G = 1 (or B = 0), then the system is forced by the additive noise, which is state-independent.The first two equations in Equation (1) represent the simple RO model with a stochastic forcing.The third equation in Equation ( 1) is the OU Langevin equation.As Jin et al. [18] suggested, the stochastic forcing term σξ(t)G represents the effects of the fast atmospheric wind anomalies (short de-correlation timescales).However, as stated previously, the fast atmospheric anomalies might own a long-term memory that cannot be expressed by the stochastic OU process.To take the long-term memory into consideration, a simple approach is to generalize the OU process to a fractional one; that is, where W α (t) is a fractional Brownian motion, and its increment dW α can be called a fractional Gaussian noise.Here, the superscript α (0 < α < 1) denotes the derivative order.Note that the FOU process Equation (2) has the same form as the OU process but with different type of noise (they will be the same if α = 1).Therefore, it can be called a fractional Langevin equation.According to [27], the fractional derivative ξ (α) of order α satisfies the usual OU process; that is, dξ and its solution is where ξ Here, the Caputo-type derivative is adopted; that is, where ξ ′ (t) is the first derivative.Applying the fractional integral of order α, we can derive the solution of Equation (2): where E t (α, −r) = (−r) −α e −rt γ(α,−r) is the Mittag-Leffler function [28].Although the analytic solution can be derived, it introduces an extra parameter that is not easy to specify.Therefore, we may just propose an alternative equation for Equation (2) or Equation ( 3) by replacing the first derivative of ξ(t) with a fractional derivative of an order α.It can be written as and its solution [29] is where ξ 0 = ξ(0) is the initial value.Compared with Equation ( 6), it only needs the initial value of the stochastic process.Moreover, if α = 1, then Equation ( 7) will reduce to the standard OU process.We still call Equation ( 7) an FOU process for convenience.
Replacing the third equation in Equation (1) with Equation ( 7), the RO model with the fractional stochastic forcing becomes Equation ( 9) is the theoretical framework that this investigation applies.Since we mainly emphasize the influence of stochastic forcing with long-term memory, we do not introduce the seasonal cycle in this investigation.It can be easily conducted to set the seasonal varying relaxation coefficient λ = λ 0 + λ ac cos(ω ac t − ϕ), where ω ac , ϕ, and λ ac are the frequency, phase, and amplitude of the annual cycle, respectively.

The Stochastic Forcing of Long-Term Memory
The only difference between this investigation and previous studies is the introduction of the fractional stochastic forcing.Therefore, it is necessary to firstly compare the two types of stochastic forcings.As shown in Figure 1, there is a significant discrepancy between them.For the OU process (Figure 1a), the stochastic noise basically declines exponentially at the beginning time, during which the influence of the Gaussian white noise is relatively weak and the system is mainly determined by its deterministic damping rate.The stochastic noise begins to oscillate around the steady state (zero value, since the exponential law declines toward zero rapidly) with a continuously increasing time, so that the Gaussian white noise begins to exhibit its influence.

The Influence of FOU Process on the RO Model
Following Jin et al. [18], the parameter values are approximately set to  1/6 month −1 ,  2/48 month −1 ,  1/1.5 month −1 ,  1/√24 month −1 , and the values of For the FOU process (Figure 1b-d), the stochastic noise represents a more pronounced oscillation, characterized by a wider amplitude range and general larger amplitudes.Furthermore, the smaller the fractional order, the more significant the oscillation.For example, the significant decay range at the beginning time may even become implicit (Figure 1b) when α = 1/3.Moreover, there are marked differences in the probability distribution functions (PDFs) between those two types of noise sequences (Figure 2).For the OU process, the maximum probability exhibits a peak at the value of zero and decreases when the values are far away from zero.Actually, this PDF is a standard normal distribution, since the time is long enough that the exponential decay signal at the beginning can be neglected and the OU process can be regarded as Brownian motion.For the FOU process, in contrast, the maximum probability is not located on the value of zero but a positive increment.The smaller the fractional order becomes, the larger the increment will be.This can be interpreted via the solution in Equation (8) or Equation ( 6), the first part of which denotes a deterministic damping between the stretched exponential and the negative power law.The former has a very fast decay, while the latter has a very slow decay for a large time [30].Therefore, even though it tends to be zero when the time tends to be infinity, it does have a relatively large value when the time is large.In the integration time range that is long enough for the RO model, the stochastic process shows a heavy tail, denoting the long-term memory.Thus, it can be predicted that compared with the original OU process, the FOU stochastic noise will bring markable change to the RO model.

The Influence of FOU Process on the RO Model
Following Jin et al. [18], the parameter values are approximately set to  1/6 month −1 ,  2/48 month −1 ,  1/1.5 month −1 ,  1/√24 month −1 , and the values of the initial system state are set as  0, ℎ 2. According to a recent study [31], the values

The Influence of FOU Process on the RO Model
Following Jin et al. [18], the parameter values are approximately set to λ = 1/6 month −1 , ω = 2π/48 month −1 , r = 1/1.5 month −1 , σ = 1/ √ 24 month −1 , and the values of the initial system state are set as T = 0, h = 2.According to a recent study [31], the values of RO parameters can be improved by objectively optimizing the RO equations fit to observations.
For the additive noise case (B = 0), the ensemble-mean solution of the linear RO model shows the dependence on the stochastic forcing (Figure 3).The numerical ensemble-mean SSTA obtained from a 5000-member ensemble (Figure 3a) suggests that the FOU process noise forcing can amplify the amplitude.This is totally different from the OU process stochastic forcing, which is indistinguishable from the deterministic solution with the same initial conditions but without stochastic forcing, according to Jin et al. [18].Compared with the deterministic SSTA amplitude (or the OU process noise), the SSTA amplitude forced by the FOU process noise can increase by around 15% compared with the initial value.The smaller the fractional derivative order, the larger the strength that the amplitude can approach.The situation for TDA (Figure 3b) is quite similar.The ensemble-mean curves do not represent a complete ENSO cycle due to the positive relaxation coefficient λ.
For the multiplicative noise case (B = 1), the ensemble-mean solution shows similar variations to the previous additive case but with a more significant oscillation amplitude (Figure 4).The SSTA (Figure 4a) can increase to exceed 1.5 • C within eight months, with an increase of around 50% (compared with the initial value) in the case of α = 1/3, far larger than the increase of around 10% in the case of α = 1, the corresponding OU-process multiplicative noise.Then, the SSTA declines to around 0.35 • C within around 34 months in the case of α = 1/3 and to around −0.16 • C within around 37 months in the case of α = 1.It can be seen that this phase transition time has been advanced for around 3 months.Similar to SSTA variations, TDA (Figure 4b) also exhibits a significantly amplified amplitude and earlier phase transition time in the case of FOU process noise.
larger than the increase of around 10% in the case of  1, the corresponding OU-process multiplicative noise.Then, the SSTA declines to around 0.35 °C within around 34 months in the case of  1/3 and to around −0.16 °C within around 37 months in the case of  1.It can be seen that this phase transition time has been advanced for around 3 months.Similar to SSTA variations, TDA (Figure 4b) also exhibits a significantly amplified amplitude and earlier phase transition time in the case of FOU process noise.

The Ensemble-Mean Dynamics
It becomes obvious that the FOU process forcing in both additive and multiplicative cases has a significant influence on the amplitude of SSTA and TDA in the ensemble-mean sense.It means that the stochastic forcing, no matter whether it is state-dependent or not, can affect the ensemble-mean state of ENSO.To interpret this phenomenon, let us discuss the ensemble-mean dynamics by separating the variables into the ensemble mean and the departure, namely,  〈〉  and ℎ 〈ℎ〉 ℎ , where T and h are the en- semble means for T and h and T  and h denote the departures, respectively.The ensemble-mean equations for the RO model are The third term on the right-hand side in the first equation of Equation ( 10) is associated with the ensemble mean for the stochastic process.The ensemble mean for the OU process converges rapidly (exponential law) to zero so that it can be neglected.However,

The Ensemble-Mean Dynamics
It becomes obvious that the FOU process forcing in both additive and multiplicative cases has a significant influence on the amplitude of SSTA and TDA in the ensemble-mean sense.It means that the stochastic forcing, no matter whether it is state-dependent or not, can affect the ensemble-mean state of ENSO.To interpret this phenomenon, let us discuss the ensemble-mean dynamics by separating the variables into the ensemble mean and the departure, namely, T = ⟨T⟩ + T ′ and h = ⟨h⟩ + h ′ , where ⟨T⟩ and ⟨h⟩ are the ensemble means for T and h and T ′ and h ′ denote the departures, respectively.The ensemble-mean equations for the RO model are The third term on the right-hand side in the first equation of Equation ( 10) is associated with the ensemble mean for the stochastic process.The ensemble mean for the OU process converges rapidly (exponential law) to zero so that it can be neglected.However, according to Equation (6) or Equation (8) the ensemble mean for the FOU process converges slowly (negative power law) to zero so that it can manifest its influence during the integration time period.In reality, the ensemble mean determined by the first term on the right-hand side of Equation (6) or Equation ( 8) has an asymptotic representation f (α)t −α when t → ∞ , where f (α) is a bounded function [30].Therefore, we may set the ensemble mean for the FOU process ⟨ξ⟩ ∼ t −α .This means that even if we only consider an additive noise (B = 0), as an external forcing, it will still exert its influence on the ensemble-mean states.Furthermore, its influence will be more significant at the beginning time when ⟨ξ⟩ has a relatively large value.This is totally different from the previous investigations [18][19][20], which suggested that the additive OU process noise had no effect on the ensemble-mean dynamics.Strictly, the additive OU process noise also plays a role, which decreases rapidly to zero.Now we know that whether a stochastic process can exert influence on the ensemble-mean dynamics or not depends on the evolution of its ensemble mean.For a stochastic process, such as the OU process, it has little effect on the ensemble-mean dynamics, since its ensemble mean converges to zero very rapidly.For a stochastic process, such as the FOU process that this investigation applies, it can affect the ensemble-mean dynamics due to its slowly convergent ensemble mean.Of course, its influence will be insignificant when the time is very large.The situation will be more complex if we consider the multiplicative noise, which means that the ensemble mean for ⟨ξT ′ ⟩ should be considered.Since it will be inconvenient when we apply the fractional derivative form of the FOU process (Equation ( 7)), we turn to its first derivative form (Equation ( 2)) to derive the ensemble-mean equations: Ignoring the higher-order term ξ 2 T ′ , the system can be enclosed and it has the same form of the RO model.We can easily derive its solution if we assume that ⟨ξ⟩ is independent of time.It is where λ + = (λ + ⟨ξ⟩σB), ω 2 = 1 2 4ω 2 − λ 2 + , A 1 , and A 2 are constants.Supposing ⟨ξ⟩ > 0, we can find out that the oscillation frequency ω 2 < ω 0 = 1 2 √ 4ω 2 − λ 2 , the frequency of the RO model.Meanwhile, the damping coefficient λ + > λ, the relaxation coefficient of the RO model.This means that the multiplicative noise will introduce a slower and stronger decayed oscillation to force the RO model.With increasing time, ⟨ξ⟩ will be smaller and smaller and tend toward zero.It will lead λ + ≈ λ and hence ω 2 ≈ ω 0 , that is, the first order term ⟨ξT ′ ⟩ has the same frequency as the RO model.Moreover, the larger the σB, the more significant the influence that the stochastic process can exert.
Substituting Equation (12) into Equation ( 10), the ensemble-mean equations can finally be written as d⟨T⟩ dt = −λ − ⟨T⟩ + ω⟨h⟩ + g(t) d⟨h⟩ dt = −ω⟨T⟩ (13) where λ − = λ − ⟨ξ⟩σB, g(t) = σ⟨ξ⟩ + σBe −λ + t (A 1 cos ω 2 t + A 2 sin ω 2 t).And its solution is where and C 2 are constants; ⟨T⟩ * is the special solution determined by g(t).Similarly, we can find out that ω 1 > ω 0 and λ − < λ, denoting a natural oscillation with higher frequency and weaker damping.It provides a theoretical explanation for previous calculation results, as Figure 4 portrays.Note that the special solution has a frequency ω 2 , which is close to the frequency of the RO model (ω 0 ) and also close to the frequency of the ensemble-mean oscillation (ω 1 ), namely, ω 2 ≤ ω 0 ≤ ω 1 .Furthermore, we may assume that ω 1 ≈ ω 0 + δ and ω 2 ≈ ω 0 − δ, where δ is a small positive frequency increment.Now, let us consider a simple case, where ⟨T⟩ ∼ sin ω 1 t + sin ω 2 t, (15) in which the natural oscillation and the forced oscillation have the same amplitude.Equation ( 15) can be rearranged as ⟨T⟩ ∼ sin ω 0 t cos δt.(16) It is clear that the ensemble-mean SSTA oscillation will show a wave train shape.The sin ω 0 t part is the carrier wave, denoting the phase propagation.The cos δt part is the wave packet, the envelop line of the carrier wave, denoting the amplitude propagation.Due to the fact that the relaxation coefficient is positive and the ensemble mean for stochastic forcing converges to zero, the ensemble-mean SSTA is a decayed oscillation.Therefore, we do not observe the wave train feature in Figure 4.

Conclusions and Discussion
To investigate the influence of fast atmospheric activities that may have a long-term memory on ENSO, we introduce stochastic forcing with long-term memory that is represented by an FOU process in a simple conceptual RO model.The FOU process shows quite different features compared with those of the OU process.The process decreases to zero very slowly with a negative power law when the time is large.During a long integration period that is long enough for the RO model, the FOU process forcing does not converge to zero rapidly like the Gaussian OU process.Therefore, the FOU process noise seems to have a non-zero mean value in the integration period.
This new feature can significantly influence the ensemble-mean dynamics of the simple RO model.The additive noise can exert is influence through introducing a non-periodic external forcing that obeys the negative power law.The external forcing can enhance the oscillation amplitude but has little effect on the oscillation period, which can be seen from Figure 3.The multiplicative noise plays its role by reducing the relaxation coefficient and by introducing a periodic external forcing, which is solved using the second-order enclosed ensemble-mean RO equations.The frequency of the periodic forcing is higher than, and close to, the frequency of the first-order enclosed ensemble-mean RO model.The two close frequencies can cause the ensemble-mean SSTA oscillation to present a wave train shape.On the other hand, a decreasing relaxation coefficient can strengthen the oscillation amplitude.Moreover, a decreasing relaxation coefficient can increase the oscillation frequency, denoting a faster oscillation.With these two mechanisms, the multiplicative noise can exert a more significant influence on the ensemble-mean dynamics than the additive noise, which can be clearly seen from Figures 3 and 4.
The finding that both additive and multiplicative noise can affect the ensemble-mean dynamics of ENSO further develops the previous knowledge that emphasized the importance of the multiplicative noise and deemed additive noise ineffective.Therefore, we may conclude that there is no essential discrepancy between the FOU process additive and multiplicative noises, although there is indeed a difference in the degree of impact.
Stochastic noises with long-term memory present extensively in nature.The simple RO framework forced by FOU-process stochastic noise can be further extended to consider more physical processes.The ensemble spread can also be discussed by applying the enclosed ensemble-mean equations.Therefore, continuous explorations may be beneficial for not only deeply understanding the influence of the fast atmospheric activities on slow ENSO variations, but also reducing the ensemble spread and improving the predictability.ENSO is the most dominant interannual signal, but not the only signal.To improve our physical understanding of internal climate variability across the global oceans, a conceptual model hierarchy that consists of the classic stochastic climate model formulated by Klaus Hasselmann, the classic RO model, and the effect of the seasonal cycle on both of these models can be also proposed to capture the essence of observed SST variability from sub-seasonal to decadal timescales [32,33].Under this model framework, the classic RO model is just a building block, only contributing to the deterministic interannual variability (labelled as T ENSO ).The total interannual SSTA variability (labelled as T) is also affected by a damping term and a stochastic forcing term, the classic Hasselmann climate model.Moreover, the seasonable cycle can be introduced to both the damping term and the RO model.Although this stochastic-deterministic model is different from the theoretical framework that this investigation applies, there is indeed linkage between them.The stochastic-deterministic model can be seen as an extension of the classic stochastic Hasselmann climate model, which is driven by deterministic interannual oscillation, seasonable cycles, and stochastic noise.The RO framework with stochastic noise also means that ENSO SSTA variability is driven by deterministic interannual oscillation and stochastic noise.Actually, if we introduce the seasonal cycle to the RO model, as previous investigations [34][35][36] have successfully demonstrated, the SSTA variability can also be concluded to these three terms, that is, deterministic interannual oscillation, seasonal cycles, and the stochastic noise.In this sense, the stochastic-deterministic model and the stochastically forced RO model are actually performing the same function if we mainly focus our discussion in the Pacific.Note that the stochastic-deterministic model still uses Gaussian white noise.This inspires us to extend the classic Hasselmann climate model to a fractional one to further develop the model.It will undoubtedly be of importance and deserve further investigations.

Figure 2 .
Figure 2. The PDF of the ensemble-mean value of the OU/FOU process noises from 50-member simulations.

Figure 2 .
Figure 2. The PDF of the ensemble-mean value of the OU/FOU process noises from 50-member simulations.

Figure 2 .
Figure 2. The PDF of the ensemble-mean value of the OU/FOU process noises from 50-member simulations.

Figure 3 .
Figure 3.The ensemble-mean evolution of SSTA (a) and TDA (b) that are derived from 5000-member simulations with additive OU/FOU process noise.

Figure 3 .
Figure 3.The ensemble-mean evolution of SSTA (a) and TDA (b) that are derived from 5000-member simulations with additive OU/FOU process noise.Fractal Fract.2024, 8, x FOR PEER REVIEW 7 of 11

Figure 4 .
Figure 4.The ensemble-mean evolution of SSTA (a) and TDA (b) that are derived from 5000-member simulations with multiplicative OU/FOU process noise.

Figure 4 .
Figure 4.The ensemble-mean evolution of SSTA (a) and TDA (b) that are derived from 5000-member simulations with multiplicative OU/FOU process noise.

Author Contributions:
Conceptualization, Y.L.; methodology, Y.L. and X.L.; software, X.L.; validation, Y.L. and X.L.; formal analysis, Y.L. and X.L.; investigation, Y.L. and X.L.; resources, X.L.; writing-original draft preparation, X.L.; writing-review and editing, Y.L.; visualization, X.L.; supervision, Y.L.; funding acquisition, Y.L.All authors have read and agreed to the published version of the manuscript.Funding: This research was jointly funded by the National Natural Science Foundation of China (Grants 42275051) and the Fundamental Research Funds for the Central Universities.Data Availability Statement: Data are contained within the article.