Level crossings and excess times due to a super-position of uncorrelated exponential pulses

TL;DR

This paper analyzes a stochastic model of intermittent fluctuations using superimposed exponential pulses, deriving excess time statistics and exploring different intermittency regimes with analytical and numerical methods.

Contribution

It provides new analytical expressions for excess time statistics in a super-position of exponential pulses, extending previous results to various intermittency levels.

Findings

Derived expressions for level crossing rates and times above thresholds.

Analytical distribution of times above threshold for strongly intermittent processes.

Numerical verification shows good agreement with known models like Ornstein-Uhlenbeck.

Abstract

A well-known stochastic model for intermittent fluctuations in physical systems is investigated. The model is given by a super-position of uncorrelated exponential pulses, and the degree of pulse overlap is interpreted as an intermittency parameter. Expressions for excess time statistics, that is, the rate of level crossings above a given threshold and the average time spent above the threshold, are derived from the joint distribution of the process and its derivative. Limits of both high and low intermittency are investigated and compared to previously known results. In the case of a strongly intermittent process, the distribution of times spent above threshold is obtained analytically. This expression is verified numerically, and the distribution of times above threshold is explored for other intermittency regimes. The numerical results compare favorably to known results for the distribution of times above the mean threshold for an Ornstein-Uhlenbeck process. This contribution generalizes the excess time statistics for the stochastic model which find applications in a wide diversity of natural and technological systems.

Equations

$$X(\Phi )=T\underset{0}{\overset{\infty }{\int }}\text{d}\dot{\Phi }\dot{\Phi }{P}_{\Phi \dot{\Phi }}(\Phi ,\dot{\Phi }).$$

$$X(\Phi )=T\frac{{\dot{\Phi }}_{\text{rms}}}{2\pi {\Phi }_{\text{rms}}}\exp \left(-\frac{{\left(\Phi -⟨\Phi ⟩\right)}^2}{2{\Phi }_{\text{rms}}^2}\right),$$

$$\tilde{\Phi }=\frac{\Phi -⟨\Phi ⟩}{{\Phi }_{\text{rms}}},$$

$$X(\tilde{\Phi })=T\frac{{\dot{\Phi }}_{\text{rms}}}{2\pi {\Phi }_{\text{rms}}}\exp \left(-\frac{{\tilde{\Phi }}^2}{2}\right).$$

$$⟨\Delta T⟩\left(\Phi \right)=\frac{T\left[1-C_{\Phi }(\Phi )\right]}{X\left(\Phi \right)}.$$

$$⟨\Delta T⟩(\tilde{\Phi })=\pi \frac{{\Phi }_{\text{rms}}}{{\dot{\Phi }}_{\text{rms}}}\mathrm{erfc}\left(\frac{\tilde{\Phi }}{\sqrt{2}}\right)\exp \left(\frac{{\tilde{\Phi }}^2}{2}\right),$$

$${\Phi }_K(t)=\sum _{k=1}^{K(T)}{A}_k\phi \left(\frac{t-t_k}{{\tau }_{\text{d}}}\right),$$

$$P_K(K)=\frac{1}{K!}{\left(\frac{T}{{\tau }_{\text{w}}}\right)}^K\exp \left(-\frac{T}{{\tau }_{\text{w}}}\right),$$

$$\phi (x)=\{\begin{array}{ll}\exp \left(x/\lambda \right),& x<0,\\ \exp \left(-x/(1-\lambda )\right),& x\ge 0,\end{array}$$

$$\gamma =\frac{{\tau }_{\text{d}}}{{\tau }_{\text{w}}},$$

$$P_A(A)=\frac{1}{⟨A⟩}\exp \left(-\frac{A}{⟨A⟩}\right),A>0,$$

$$P_{\Phi }(\Phi )=\frac{1}{⟨A⟩\Gamma (\gamma )}{\left(\frac{\Phi }{⟨A⟩}\right)}^{\gamma -1}\exp \left(-\frac{\Phi }{⟨A⟩}\right),\Phi >0,$$

$$1-C_{\Phi }(\Phi )=Q\left(\gamma ,\gamma \Phi /⟨\Phi ⟩\right),$$

$$⟨\exp (i\Phi u)⟩=(1-i⟨A⟩u{)}^{-\gamma }.$$

$${\kappa }_n=\gamma ⟨A^n⟩I_n,$$

$$I_n=\underset{-\infty }{\overset{\infty }{\int }}\text{d}x{\left[\phi (x)\right]}^n.$$

$${\kappa }_n=(n-1)!\gamma {⟨A⟩}^n.$$

$$⟨\Phi ⟩$$

$${\Phi }_{\text{rms}}^2$$

$$S_{\Phi }$$

$$F_{\Phi }$$

$${\Theta }_K(t)=\frac{{\tau }_{\text{d}}}{2}\frac{\text{d}{\Phi }_K}{\text{d}t}=\sum _{k=1}^{K(T)}{A}_k\vartheta \left(\frac{t-t_k}{{\tau }_{\text{d}}}\right),$$

$$\vartheta (s)=\frac{1}{2}\frac{\text{d}\phi }{\text{d}s}=\frac{1}{2}\{\begin{array}{ll}(\lambda {)}^{-1}\exp \left(s/\lambda \right),& s<0,\\ -(1-\lambda {)}^{-1}\exp \left(-s/(1-\lambda )\right),& s\ge 0.\end{array}$$

$$⟨\Phi \Theta ⟩=\frac{{\tau }_{\text{d}}}{4}\frac{\text{d}}{\text{d}t}⟨{\Phi }^2⟩=0.$$

$$⟨\Theta ⟩$$

$${\Theta }_{\text{rms}}^2$$

$$S_{\Theta }$$

$$F_{\Theta }$$

$$⟨\exp (i\Theta v)⟩={\left(1-i⟨A⟩\frac{v}{2\lambda }\right)}^{-\lambda \gamma }{\left(1+i⟨A⟩\frac{v}{2(1-\lambda )}\right)}^{-(1-\lambda )\gamma }.$$

$$P_{\tilde{\Theta }}(\tilde{\Theta })=\sqrt{\frac{2\gamma }{\pi }}\frac{2^{-\gamma /2}}{\Gamma (\gamma /2)}{\left|\sqrt{\gamma }\tilde{\Theta }\right|}^{(\gamma -1)/2}{\mathcal{K}}_{(\gamma -1)/2}\left(\left|\sqrt{\gamma }\tilde{\Theta }\right|\right),$$

$$P_{\Phi \Theta }\left(\Phi ,\Theta \right)=\frac{1}{{\left(2\pi \right)}^2}\underset{-\infty }{\overset{\infty }{\int }}\text{d}u\underset{-\infty }{\overset{\infty }{\int }}\text{d}v\exp \left(-i\Phi u-i\Theta v\right)⟨\exp \left(iu\Phi +iv\Theta \right)⟩.$$

$$⟨\exp \left(iu\Phi +iv\Theta \right)⟩=\exp \left(\gamma \underset{-\infty }{\overset{\infty }{\int }}\text{d}A{P}_A(A)\underset{-\infty }{\overset{\infty }{\int }}\text{d}s\left[\exp \left(iuA\phi (s)+ivA\vartheta (s)\right)-1\right]\right).$$

$$⟨\exp (iu\Phi +iv\Theta )⟩={\left[1-i⟨A⟩\left(u+\frac{v}{2\lambda }\right)\right]}^{-\gamma \lambda }{\left[1-i⟨A⟩\left(u-\frac{v}{2(1-\lambda )}\right)\right]}^{-\gamma \left(1-\lambda \right)}.$$

Full text

I Introduction

A stochastic process given by a super-position of uncorrelated pulses can be considered as a reference model for intermittent fluctuations in physical systems. It has found applications in a broad range of fields, including economics, electronics, fission chambers, magnetically confined fusion plasmas, meteorology, oceanography and optics.Campbell (1909); Rice (1944); Segal et al. (1985); Fesce (1986); Jang (2004); Claps, Giordano, and Laio (2005); Lefebvre (2008); Garcia (2012); Elter et al. (2015); Dalmao and Mordecki (2015) In many of these applications, the failure or survival of the system depends sensitively on the frequency of large-amplitude fluctuations and the duration of times spent above a critical threshold level. Accordingly, much work has been done in order to calculate the rate of level crossings and average excess times above a threshold level.Rice (1945); Bar-David and Nemirovsky (1972); Barakat (1988); Leadbetter and Spaniolo (2004); Alili, Patie, and Pedersen (2005); Daly and Porporato (2010); Yi (2010); Biermé and Desolneux (2012); Theodorsen and Garcia (2016); Yura and Hanson (2010)

This contribution is primarily motivated by turbulent flows in the boundary region of magnetically confined plasmas. Evidence points towards these fluctuations being caused by filamentary structures transporting particles and heat towards main chamber walls.D’Ippolito, Myra, and Zweben (2011); Zweben et al. (2007) Experimental results provide strong evidence that large-amplitude plasma fluctuations in the boundary region can be described as a super-position of uncorrelated pulses with fixed, exponential pulse shape of constant duration and exponentially distributed pulse amplitudes, with exponentially distributed waiting times between the pulse arrivals.Garcia et al. (2013); Garcia, Horacek, and Pitts (2015); Theodorsen et al. (2016); Kube et al. (2016); Theodorsen et al. (2016); Garcia et al. (2016a); Theodorsen et al. (sion); Kube et al. (shed) A stochastic model with these properties has Gamma distributed amplitudes, a parabolic relation between the skewness and flatness moments, an exponential autocorrelation function and a Lorenzian power spectrum.Garcia (2012); Garcia et al. (2016b); Theodorsen, Garcia, and Rypdal (2017)

This stochastic model can be extended in several ways, including adding a noise term,Theodorsen, Garcia, and Rypdal (2017) using different pulse shapesLowen and Teich (2005); Pécseli (2000); Garcia and Theodorsen (2017a, b) or distributions of amplitudes,Daly and Porporato (2010); Theodorsen et al. (2016) or allowing for a distribution of pulse durations.Garcia and Theodorsen (2017a, b) In this contribution, which is an extended version of LABEL:theodorsen-pop-2016, the rate of threshold crossings and average time above a given threshold are derived and discussed in the case of exponential pulses with fixed duration and shape.

Given the joint probability density function (PDF) $P_{\Phi\dot{\Phi}}(\Phi,\dot{\Phi})$ for a stationary random variable $\Phi(t)$ and its derivative $\dot{\Phi}=\text{d}\Phi/\text{d}t$, the number of up-crossings of the level $\Phi$ in a time interval of duration $T$ is given by integrating over all positive values of the derivative,Rice (1945); Kristensen et al. (1991); Sato, Pécseli, and Trulsen (2012); Fattorini et al. (2012)

[TABLE]

For independent, normally distributed $\Phi$ and $\dot{\Phi}$, this gives the celebrated result known as the Rice formula,Rice (1945); Kristensen et al. (1991); Sato, Pécseli, and Trulsen (2012); Fattorini et al. (2012); Dalmao and Mordecki (2015)

[TABLE]

where $\langle{\Phi}\rangle$ is the mean value of $\Phi$ and $\Phi_{\text{rms}}$ and $\dot{\Phi}_{\text{rms}}$ are the standard deviation or root mean square (rms) values of $\Phi$ and $\dot{\Phi}$, respectively. Here and in the following, $\langle{\bullet}\rangle$ denotes an average over all random variables. The number of level crossings is clearly largest for threshold values close to the mean value of $\Phi$. In this contribution, we will frequently use the normalization

[TABLE]

giving

[TABLE]

The average time $\langle\Delta T\rangle$ spent above a threshold value $\Phi$ by the stationary process is given by the ratio of the total time spent above the level $\Phi$ to the number of up-crossings $X$ in an interval of duration $T$. The former is by definition given by $T[1-C_{\Phi}(\Phi)]$, where $C_{\Phi}(\Phi)$ is the cumulative distribution function (CDF) of $\Phi$. This gives the average excess time as

[TABLE]

For jointly normally distributed $\Phi$ and $\dot{\Phi}$ with zero correlation (that is, the processes are independent), the average excess time is given by Kristensen et al. (1991); Sato, Pécseli, and Trulsen (2012); Fattorini et al. (2012)

[TABLE]

where $\operatorname{erfc}(x)$ denotes the complementary error function (and $\operatorname{erf}(x)$ is the error function).Olver et al. Here and in the following, $x$ denotes a real, unitless variable, used in definitions of special functions. It should be noted that the standard deviation of $\dot{\Phi}$, which appears in both Eqs. (4) and (6), is challenging to estimate from measurement data, thus limiting the usefulness of the expressions above.

The goal of this contribution is to derive expressions for level crossing rates and excess times for a filtered Poisson process (FPP). In Sec. II, the FPP with a two-sided exponential pulse shape with fixed pulse duration time and exponentially distributed amplitudes is introduced, and some of its statistical properties are reviewed. The derivative of the process is discussed and the joint PDF between the process and its derivative is derived. In Sec. III, expressions for the rate of level crossings and the average excess time for the FPP are given. Limits of a one-sided exponential pulse shape, the normal limit and the limit of strong intermittency are discussed in detail. In Sec. IV, we discuss the distribution of excess times in the strong intermittency limit and in the normal limit. Sec. V gives numerical results for the distribution of excess times in the general case, and compares this to the analytic expressions from Sec. IV. The convergence of the rate of level crossings to its analytic expression is also considered in Sec. IV. Concluding remarks are given in Sec. VI.

II The filtered Poisson process

In this section, the FPP is introduced and its general features are discussed. First, we present the distribution and moments of the FPP. Secondly, the derivative of the FPP is derived, and its distribution and moments are presented. Lastly, we derive and discuss the joint PDF between the FPP and its derivative.

II.1 Super-position of pulses

The FPP can be described as a super-position of uncorrelated pulses, Rice (1944); Garcia (2012); Garcia et al. (2016b); Kube and Garcia (2015); Kube et al. (2016); Theodorsen and Garcia (2016); Theodorsen et al. (2016); Garcia et al. (2016a); Theodorsen et al. (sion); Kube et al. (shed); Campbell (1909); Parzen (1999); Lowen and Teich (2005); Pécseli (2000)

[TABLE]

where for event $k$, $t_{k}$ is the pulse arrival time and $A_{k}$ is the pulse amplitude. The pulse duration time $\tau_{\text{d}}$ and the pulse shape $\varphi(x)$ are assumed to be the same for all events. We will assume the waiting time between pulses to be uncorrelated and exponentially distributed with mean waiting time $\tau_{\text{w}}$. From this it follows that $K(T)$ is Poisson distributed with constant rate $1/\tau_{\text{w}}$,

[TABLE]

and therefore that the pulse arrival times $t_{k}$ are uniformly distributed on $[0,T]$.

In the following, the pulse shape is described by a two-sided exponential function

[TABLE]

where $\lambda$ is a pulse asymmetry parameter restricted to the range $0<\lambda<1$. The ratio between the pulse duration and average waiting time,

[TABLE]

is called the intermittency parameter, and determines the degree of pulse overlap. It is the most fundamental parameter of the stochastic model. Realizations of this process for various values of $\gamma$ are shown in Fig. 1, using the normalization given in Eq. (3). For small $\gamma$, the pulses are well separated and the process is strongly intermittent. For large $\gamma$, there is significant pulse overlap and realizations of the process resembles random noise, with relatively small and symmetric fluctuations around the mean value. For intermediate $\gamma$, large-amplitude bursts can be constructed from one separate large-amplitude pulse, or several smaller amplitude pulses. Because of this, the parameter $\gamma$ can be interpreted as an intermittency parameter for the process, with low values of $\gamma$ giving a highly intermittent process and high values of $\gamma$ giving a weakly intermittent process.

Assuming the PDF of the pulse amplitudes $A$ is an exponential distribution,

[TABLE]

the stationary distribution of the random variable $\Phi(t)$ can be shown to be a Gamma distribution with shape parameter $\gamma$ and scale parameter $\left<A\right>$;Bondesson (1982); Daly and Porporato (2010); Garcia (2012); Garcia et al. (2016b); Garcia and Theodorsen (2017b)

[TABLE]

where $\Gamma(x)$ is the gamma function.Olver et al. The complementary CDF of $\Phi$ is then given by

[TABLE]

where $Q(a,x)$ is the regularized upper incomplete gamma function with parameter $a$.Olver et al. In this contribution, we will also use the upper incomplete gamma function $\Gamma(a,x)$.Olver et al. $Q(a,x)$ is defined as $Q(a,x)=\Gamma(a,x)/\Gamma(a)$.

The complementary CDF of $\Phi$ as a function of $\widetilde{\Phi}$ for various values of $\gamma$ is presented in Fig. 2. This function can be interpreted as the fraction of time a signal spends above the threshold $\widetilde{\Phi}$. As $\gamma$ increases, the PDF approaches a normal distribution. In the normal regime $\gamma\gg 1$, the fraction of time above threshold falls rapidly with increasing threshold level since the fluctuations in the signal are concentrated around the mean value. In the strong intermittency regime, $\gamma\ll 1$, the signal spends long periods of time close to zero value as few pulses overlap significantly. Thus, the total time above threshold increases rapidly as the threshold approaches zero. Also note that for large values of $\Phi$, the total time above threshold is orders of magnitude higher for a process with high intermittency than for a process with low intermittency. For $\gamma=1$, the PDF of $\Phi$ is an exponential distribution.

Correspondingly, the characteristic function of $\Phi$ is

[TABLE]

It can likewise be shown that the cumulants of the process for arbitrary pulse shape and amplitude distribution are given byGarcia (2012); Garcia et al. (2016b)

[TABLE]

where

[TABLE]

For the pulse shape given in Eq. (9), $I_{n}=1/n$ and using that $\langle{A^{n}}\rangle=n!\,\left<A\right>^{n}$ for exponentially distributed amplitudes, the cumulants for $\Phi$ areGarcia (2012); Garcia et al. (2016b)

[TABLE]

Note that the cumulants, and therefore also the PDF, are independent of the pulse asymmetry parameter $\lambda$. Given the cumulants, we can find the lowest order moments of the process:Rice (1944); Bondesson (1982); Garcia (2012); Garcia et al. (2016b)

[TABLE]

Here, $S_{\Phi}$ is the skewness of the random variable $\Phi$, and $F_{\Phi}$ is its flatness. The relative fluctuation level is $\Phi_{\text{rms}}/\left<\Phi\right>=1/\gamma^{1/2}$. There is a parabolic relation between skewness and flatness: $F_{\Phi}\left(S_{\Phi}\right)=3+3S_{\Phi}^{2}/2$. It can be shown that the distribution of the normalized process $\widetilde{\Phi}=(\Phi-\left<\Phi\right>)/\Phi_{\text{rms}}$ resembles a standard normal distribution (that is, a normal distribution with zero mean and unit standard deviation) in the limit $\gamma\to\infty$, independent of pulse shape and amplitude distribution.Rice (1944) In this case, both the skewness $S_{\Phi}$ and the excess kurtosis $F_{\Phi}-3$ vanish.Rice (1944); Garcia (2012) Conversely, for $\gamma\to 0$, the skewness and kurtosis moments both tend to infinity. Note that from the definition of skewness and flatness, it follows that $S_{\Phi}=S_{\widetilde{\Phi}}$ and $F_{\Phi}=F_{\widetilde{\Phi}}$.

Note that $\Phi$ is non-negative, giving $\widetilde{\Phi}\geq-\gamma^{1/2}$. By contrast, a normally distributed random variable has infinite support. The difference between the PDF of $\widetilde{\Phi}$ and a standard normal distribution due to this discrepancy is negligible in practice, since values of $-\gamma^{1/2}$ or less are highly unlikely for a standard normal distribution in the case of $\gamma\gg 1$.

II.2 The derivative of the filtered Poisson process

In order to calculate the joint distribution of the process and its derivative, the normalized time derivative is defined by

[TABLE]

where the pulse shape is given by

[TABLE]

Here, we have divided by a factor 2 in order for the pulse shape to fulfill $\int_{-\infty}^{\infty}\text{d}s\,{\left|\vartheta(s)\right|}=1$.111This is a correction to $\vartheta$ given in LABEL:theodorsen-pop-2016. This leads to corrections in the rms-value of $\Theta$ and the joint PDF between $\Phi$ and $\Theta$, but not the excess statistics. This is another stochastic process of the same type as that given in Eq. (7), but with a different pulse shape. Since the process $\Phi(t)$ is stationary, it follows that $\left<{\Theta}\right>=0$, which is easily verified from Eq. (17). The processes $\Phi(t)$ and $\Theta(t)$ are evidently dependent yet also uncorrelated,

[TABLE]

In Appendix A, the joint PDF between $\Phi$ and $\Theta$ is used to demonstrate that $\Phi$ and $\Theta$ become independent in the limit $\gamma\to\infty$.

The lowest order moments of $\Theta$ are readily calculated as

[TABLE]

In the limit of $\lambda\to 0$ or $\lambda\to 1$, the moments $\Theta_{\text{rms}}$, $S_{\Theta}$ and $F_{\Theta}$ diverge, meaning the PDF of $\Theta$ does not exist in this case. In these limits, the pulse shape in $\Phi$ is discontinuous and the derivative of the pulse shape contains delta functions. Thus we require the two-sided exponential pulse shape to calculate the rate of level crossings. It will later be shown that these limits exist for the rate of level crossings and are consistent with other methods starting from the one-sided exponential pulse shape. Thus, while this method cannot be used to calculate the rate of level crossings for a discontinuous signal, the rate still exists.Bar-David and Nemirovsky (1972); Daly and Porporato (2010); Biermé and Desolneux (2012)

Using the same approach as in Refs. Garcia, 2012 and Garcia et al., 2016b, the characteristic function of $\Theta$ is given by

[TABLE]

This characteristic function can be interpreted as originating from the sum of two independent gamma distributed variables, one over positive values with shape parameter $\gamma\lambda$ and scale parameter $\left<A\right>/(2\lambda)$, and the other over negative values with shape parameter $\gamma(1-\lambda)$ and scale parameter $\left<A\right>/[2(1-\lambda)]$. The PDF of this compound process is a convolution of the two gamma distributions, which to the best of the authors knowledge does not have a closed form. Still, the argument in Refs. Rice, 1944; Garcia, 2012; Garcia et al., 2016b applies here as well, and the PDF of $\Theta$ resembles a normal distribution in the limit $\gamma\to\infty$. In Fig. 3, realizations for $\widetilde{\Theta}$ are presented for $\lambda=1/4$ and various values of $\gamma$. Arrival times and pulse amplitudes are the same as in Fig. 1. Again, the process is strongly intermittent for low values of $\gamma$, and resembles random noise for high values of $\gamma$.

By choosing $\lambda=1/2$, the pulse shape $\varphi(x)$ is symmetric. In this case, the characteristic function in Eq. (23) has an inverse transformation in closed form, and the corresponding PDF is given by

[TABLE]

where $\mathcal{K}_{a}(x)$ is the modified Bessel function of the second kind.Olver et al. This PDF is presented in Fig. 4 for various values of $\gamma$. For small values of $\gamma$, this PDF has exponential tails and is sharply peaked at the mean value, while it resembles a normal distribution for large values of $\gamma$. The same PDF for $\gamma=2$ and various values of $\lambda$ is presented in Fig. 5. As the asymmetry parameter approaches [math], the skewness and flatness of $\Theta$ increases. It can be seen from Eq. (23) that in the case $\lambda=1/2$ and $\gamma=2$, $\Theta$ is symmetrically Laplace distributed with zero mean and standard deviation $\Theta_{\text{rms}}=2\left<A\right>^{2}$.Kozubowski and Podgórski (2000)

II.3 The joint PDF of the filtered Poisson process

The joint PDF of $\Phi$ and $\Theta$ is generally given by

[TABLE]

Using that individual events are uncorrelated and that the number of pulses is Poisson distributed, the characteristic function of $\Phi$ and $\Theta$ can be calculated as

[TABLE]

This expression is given in Refs. Pécseli, 2000; Sato, Pécseli, and Trulsen, 2012 for the case of fixed (degenerately distributed) pulse amplitudes, although the generalization is straightforward. Exchanging the order of integration, we find that

[TABLE]

We note that we recover the expression for the characteristic function of $\Phi$ in Eq. (14) by setting $v=0$ in this equation, and we recover the characteristic function of $\Theta$ in Eq. (23) by setting $u=0$. Substituted into Eq. (25), the stationary joint PDF can be obtained in closed form. We change variables to $x=\left<A\right>[u+v/(2\lambda)]$ and $y=\left<A\right>\{u-v/[2(1-\lambda)]\}$, and use the notation

[TABLE]

The joint PDF can now be written as

[TABLE]

The integrals can be performed separately, and we get the closed form expression

[TABLE]

This is non-zero only for the limited range $-\Phi/[2(1-\lambda)]<\Theta<\Phi/(2\lambda)$, which follows from the fact that the signal $\Phi(t)$ cannot decrease faster than the rate of decay of individual pulse structures, nor increase slower than the rate of growth of individual pulses, since the pulse amplitudes are positive definite. The dependence between $\Phi$ and $\Theta$ is evident from Eq. (29), since the joint PDF is not separable into a product of the marginal PDFs. As expected, $P_{\Phi}(\Phi)$ can be recovered by integrating over $\Theta$. Also note that the expression for the joint PDF diverges in the limits $\lambda\to 0$ and $\lambda\to 1$, as was the case for the moments and PDF of $\Theta$. As the PDFs of both $\Phi$ and $\Theta$ resemble normal distributions in the limit $\gamma\to\infty$ and they are uncorrelated, the joint PDF for $\Phi$ and $\Theta$ resembles the product of two normal distributions, that is, a joint normal distribution with vanishing correlation coefficient. This is demonstrated explicitly in Appendix A. Thus, in the normal limit $\gamma\to\infty$, the classical Rice formula given by Eq. (2) is recovered. As in the case of $P_{\Phi}$, there is a discrepancy between $P_{\Phi\Theta}$ and a joint normal distribution due to the limited region of non-zero values of $P_{\Phi\Theta}$. The domain of non-zero values can be written as $-(\widetilde{\Phi}+\gamma^{1/2})/(1-\lambda)<\widetilde{\Theta}/\sqrt{\lambda(1-\lambda)}<(\widetilde{\Phi}+\gamma^{1/2})/\lambda$, where $\widetilde{\Theta}=\Theta/\Theta_{\text{rms}}$. For standard normally distributed variables, values outside of this domain are highly unlikely in the case of $\gamma\gg 1$, and this discrepancy is in practice negligible.

The joint distribution $P_{\Phi\Theta}(\Phi,\Theta)$ is presented in Fig. 6 for $\gamma\in\{10^{-1},1,10\}$ and $\lambda\in\{1/4,1/2\}$. It should be noted that logarithmic scaling is used for $\gamma=10^{-1}$ and $1$, while linear scaling is used for $\gamma=10$. The white area in all figures are the regions where $P_{\Phi\Theta}$ vanishes, as given by Eq. (29). The joint distribution for $\gamma\leq 1$ diverges at $\Phi=0$ and $\Theta=0$, corresponding to $\widetilde{\Theta}=0$, $\widetilde{\Phi}=-\gamma^{1/2}$, since the pulses arrive rarely enough for the signal to fall close to zero value for long time durations. In this case, the signals are very likely to decay or grow undisturbed at the rate of individual pulses, explaining the increased value of the joint distribution near the lines $\Theta=-\Phi/[2(1-\lambda)]$, $\Theta=\Phi/(2\lambda)$.

III Excess time statistics

In this section we present the rate of threshold crossings and average time above threshold for the FPP. Limits of one-sided exponential pulse shape, and weak and strong intermittency are explored and compared to previous works.

III.1 Formulation of excess time statistics

The rate of up-crossings above a threshold level $\Phi$ is now readily calculated from Eqs. (1) and (29) as

[TABLE]

which, together with the complementary CDF in Eq. (13), gives the average time above the threshold for each threshold crossing,

[TABLE]

Note that both Eqs. (30) and (31) can be written as a pre-factor depending on $\gamma$ and $\lambda$, multiplied by a function of $\gamma$ and the variable

[TABLE]

This indicates that the functional shape of both equations with threshold level depend only on the intermittency parameter $\gamma$, while the function value depends on both $\gamma$ and $\lambda$. By contrast, the complementary CDF given by Eq. (13) does not depend on $\lambda$.

From the joint PDF in Eq. (29), it is clear that the dependency between $\Phi$ and $\Theta$ is important for the rate of threshold crossings. In order to investigate the effect of this dependency, we calculate the rate of threshold crossings divided by the PDF of $\Phi$:

[TABLE]

On the other hand, starting from Eq. (30) and assuming $\Phi$ and $\Theta$ are independent gives

[TABLE]

which is independent of $\Phi$. Thus an assumption of independence will always give the wrong algebraic factor, although this is not very relevant for large $\Phi$ where the exponential term dominates. Also note that Eq. (34) gives the correct result in the limit $\gamma\to\infty$, where the process and its derivative are indeed independent. However, inserting the PDF of $\Theta$ from Sec. II.2 into Eq. (34) gives a surprisingly complicated result, presented in Appendix B. There is significant discrepancy between this expression and the prefactor in Eq. (33). Thus accurately accounting for the dependency between $\Phi$ and $\Theta$ is necessary for correctly predicting the rate of threshold crossings.

The rate of up-crossings as function of the threshold level for various values of $\gamma$ is presented in Fig. 7. Full lines show the case of $\lambda=1/2$, while dashed lines show the rate of level crossings in the limit $\lambda\to 0$. The analytical expression in this limit will be discussed further in Sec. III.2. The total number of crossings is evidently proportional to the length of the time series $T$ and inversely proportional to the pulse duration $\tau_{\text{d}}$. The rate of threshold crossings is highest for thresholds close to the mean value of the process in all cases. In the normal regime $\gamma\gg 1$, there are comparatively few crossings for threshold levels much smaller or much larger than the mean value due to the low probability of large-amplitude fluctuations. The rate of level crossings is therefore a narrow Gaussian function in this limit. In the strong intermittency regime, $\gamma\ll 1$, the signal spends most of the time close to zero value, and virtually any pulse arrival will give rise to a level crossing for finite threshold values. As seen in Fig. 7, the rate of level crossings approaches a step function in this limit. For $\lambda=1/2$, the rate of level crossings at the mean value, $\widetilde{\Phi}=0$, approaches a definite value. In Sec. III.3 this value is shown to be $1/\pi$. In contrast, there is no limiting value for $\lambda=0$. In this case $X(\widetilde{\Phi}=0)\to\infty$ as $\gamma\to\infty$, as will be demonstrated in Sec. III.3.

The average time above threshold is presented in Fig. 8 for various values of $\gamma$. Full lines show the case of $\lambda=1/2$, while dashed lines show the average time above threshold in the limit $\lambda\to 0$. While both the rate of threshold crossings and the fraction of time above threshold vary qualitatively as $\gamma$ changes, the shape of the average time above threshold is fairly similar. In all cases the average excess time decreases monotonically with the threshold level, with a fast drop for small threshold values. This is followed by a slow tapering off for large threshold values. For the range of intermittency parameters considered here, the average excess time is of the order of the pulse duration or shorter for large threshold values. For $\lambda\to 0$, the average time above threshold decreases by about half a decade for each tenfold increase in $\gamma$, but the functional shape varies little. For $\lambda=1/2$, the average time above threshold converges to the Rice result, as will be shown in Sec. III.3. It can be shown that for given $\gamma$ and $\lambda$, $\left<\Delta T\right>/\tau_{\text{d}}$ scales as $1/\widetilde{\Phi}$ in the limit $\widetilde{\Phi}\to\infty$. As the threshold value increases above the mean signal value, up-crossings of the threshold become fewer while the signal spends less time in total above the threshold. Evidently these two effects nearly cancel, and the average excess time decreases slowly with increasing threshold level.

III.2 Limit of the one-sided pulse shape

As stated in Sec. II.3, the limit of the one-sided exponential pulse shape does not exist for $P_{\Theta}$ or $P_{\Phi\Theta}$. This is due to the fact that the pulse shape $\varphi(x)$ is discontinuous in this case, and therefore second and higher order moments of its derivative do not exist. However, the rate of level crossings for the discontinuous process still exists, and has been discussed in for example Refs. Bar-David and Nemirovsky, 1972; Biermé and Desolneux, 2012; Daly and Porporato, 2010. Taking either of the limits $\lambda\to 0$ and $\lambda\to 1$ give the same result, and yield

[TABLE]

This result was also obtained in LABEL:bierme-2012 by considering the Fourier transform of the number of level crossings. Since the complementary CDF of $\Phi$ does not depend on $\lambda$, the total time the signal spends above threshold remains unchanged, and the average time above threshold is simply

[TABLE]

The functional shape of Eqs. (35) and (36) are the same as in the more general expressions given by Eqs. (30) and (31), since $\lambda$ only appears in the prefactor of these equations. The approach discussed in LABEL:daly-2010 also leads to the results presented in this section, although they are not explicitly given in the reference.

III.3 The normal limit

In the limit of large $\gamma$, the expression for $X(\Phi)$ can be simplified and shown to be equal to the case for a normally distributed process. Using Stirling’s approximation for the Gamma functions in Eq. (30), we have in the normal limit:

[TABLE]

Inserting this result into Eq. (30), and using the normalized threshold in Eq. (32), the rate of crossings in the weak intermittency case $\gamma\gg 1$ can be written as

[TABLE]

In Appendix C, we show that

[TABLE]

and the rate of level crossings in the limit $\gamma\to\infty$ can be written as

[TABLE]

This expression is equal to Eq. (4), when using $\Phi_{\text{rms}}$ from Eq. (18b) and $\dot{\Phi}_{\text{rms}}=2\Theta_{\text{rms}}/\tau_{\text{d}}$ from Eq. (22b). As mentioned in the discussion of Fig. 7, in the case of $\lambda=1/2$, we have that $\lim_{\gamma\to\infty}\tau_{\text{d}}X(\widetilde{\Phi}=0)/T=1/\pi$.

In Appendix D, it is shown that

[TABLE]

and the expression for the average time above threshold in Eq. (31) can be shown to be equivalent to the expression given by Eq. (6) in the case $\gamma\to\infty$. Note that for $\lambda=1/2$, we have the limit $\lim_{\gamma\to\infty}\left<\Delta T\right>(\widetilde{\Phi}=0)/\tau_{\text{d}}=\pi/2$.

Starting from Eq. (35) and going through the same procedure as above, we have in the cases $\lambda=0$ and $\lambda=1$

[TABLE]

There is a clear discrepancy between Eqs. (40) and (42), suggesting a qualitative difference in the level crossing rate for a continuous and discontinuous pulse shape. This result is in agreement with the careful analysis in LABEL:bierme-2012. The rate of level crossings is much higher for a process with jumps in the pulse shape (and continues to increase with the square root of $\gamma$ as $\gamma$ increases). No matter how strong the pulse overlap is, the discontinuous pulses are much more likely to trigger threshold crossings than the continuous pulses.

We further note that the average time above threshold for $\lambda\in{0,1}$ can be written as

[TABLE]

Just as the rate of level crossings increases without bound for increasing pulse overlap in the cases $\lambda=0$ and $\lambda=1$, the average time above threshold decreases with increasing $\gamma$. Thus, in the normal limit, the process is characterized by frequent threshold crossings but short excess times. In the case of a discontinuous pulse shape, the derivative of the process does not exist, and the method we have used to find the rate of threshold crossings is not valid (but still gives results in agreement with other methods). In this case, Rice’s formula, Eq. (6) does not exist for the process (as $\Theta_{\text{rms}}$ does not exist). Thus, the rate of pulse arrivals will always play a role in the expressions for the rate of threshold crossings and average excess times.

III.4 The strong intermittency limit

We will now investigate the limit of $\gamma\to 0$, where we can neglect overlap of individual pulses, such that each pulse appears as one isolated burst in realizations of the process. In this section, we will use $\Phi/\left<A\right>$ instead of the expressions in Eq. (32), to avoid $\gamma$ where possible. In the previous section, $\widetilde{\Phi}$ approached a standard, normally distributed variable. Here, $\widetilde{\Phi}$ approaches a random variable with infinite skewness and flatness, and the advantage of normalizing the signal to remove the dependence on $\left<A\right>$ is diminished. In the limit $\gamma\to 0$, we can find the number of threshold crossings, the average time above threshold and even the distribution of time above threshold for each up-crossing without going through the joint PDF of $\Phi$ and $\Theta$.

For non-overlapping pulses, the total number of upward crossings of the threshold must be the same as the total number of pulses with amplitude higher than the threshold value. Therefore, the total number of up-crossings can be written as

[TABLE]

where $\langle{K}\rangle=T/\tau_{\text{w}}=\gamma T/\tau_{\text{d}}$ and $\varphi_{\text{max}}$ is the largest positive value of $\varphi$. For the exponential pulse shape in Eq. (9), $\varphi_{\text{max}}=\varphi(0)=1$. This expression can also be reached by taking the limit $\gamma\to 0$ in either Eqs. (30) or (35), suggesting that the number is the same for a continuous and a discontinuous pulse. This can be explained by the fact that each sufficiently large-amplitude pulse triggers one crossing above the threshold, and this is independent of the pulse shape.

Using the complementary CDF from Eq. (13), we have the total time above the threshold level $\Phi$ in the strong intermittency limit,

[TABLE]

Estimating $\left<\Delta T\right>$ by $T\,(1-C_{\Phi})/X$, given by Eqs. (44) and Eq. (45), we find that the average time above threshold for each level crossing is given by

[TABLE]

The rate of level crossings, given by Eq. (44), and the fraction of time above threshold, given by Eq. (45), both decay as $\gamma$ in the limit $\gamma\to 0$. Since the dependency of these two expressions on $\gamma$ is the same, the average time the signal spends above the threshold is independent of $\gamma$ in the strong intermittency limit.

IV The distribution of excess times

In this section, we will investigate the PDF of the times spent above threshold. In the strong intermittency limit, there is a closed analytical expression for this distribution. In the normal limit, with $\lambda\to 0$, an analytical expression can also be found for crossings above the mean threshold value, but it depends explicitly on the intermittency parameter $\gamma$. In the following, we will use $L$ to denote the threshold value.

IV.1 The strong intermittency limit

In this section, we will derive the PDF of the time above threshold in the case when overlap of pulses can be neglected, that is, the strong intermittency limit $\gamma\to 0$. For brevity of notation, we will not include the limit in the following. We will also assume $\varphi_{\text{max}}=\varphi(0)=1$. Generalization to arbitrary $\varphi_{\text{max}}$ is done by replacing the threshold $L$ by $L/\varphi_{\text{max}}$.

For a given pulse with amplitude $A>L$, the signal spends a time ${\Delta T}$ above the threshold. With the two-sided exponential pulse shape, ${\Delta T}$ can be divided into a time before the peak, ${\Delta T}_{-}$, and a time following the peak, ${\Delta T}_{+}$. Assuming the pulse has peak amplitude at time $t=0$, the pulse crosses the threshold $L$ upwards at time ${\Delta T}_{-}$, given by $L=A\exp({\Delta T}_{-}/\lambda\tau_{\text{d}})$, which gives

[TABLE]

Similarly, the pulse crosses the threshold downwards at time ${\Delta T}_{+}$, given by $L=A\exp(-{\Delta T}_{+}/[(1-\lambda)\tau_{\text{d}}])$, which gives

[TABLE]

Thus, the total time that the pulse spends above the threshold is

[TABLE]

and the pulse asymmetry plays no further role. Note that ${\Delta T}$ is always positive, since $A>L$ by assumption. Using that $A$ is exponentially distributed with mean value $\left<A\right>$, the conditional PDF of $A$ given that $A>L$, is given by the truncated exponential distributionStark, H.: Woods (2012)

[TABLE]

Changing the random variable from $A$ to ${\Delta T}$ and ensuring proper normalization for the PDF of excess times gives

[TABLE]

This is the so-called Gompertz distribution with parameters $L/(\tau_{\text{d}}\left<A\right>)$ and $1/\tau_{\text{d}}$, see e.g. Ref. Olive, 2014, Ch. 10.20. It is presented in Fig. 9 for various values of $L/\left<A\right>$. For $L\geq\left<A\right>$, the PDF decays monotonically from ${\Delta T}=0$, while for $L<\left<A\right>$, the PDF has a maxima at ${\Delta T}/\tau_{\text{d}}=\ln(\left<A\right>/L)$.

The mean value of the Gompertz distribution can be calculated as

[TABLE]

which is equivalent to the expression in Eq. (46). The PDF of ${\Delta T}$ will be compared to synthetic data in Sec. V.1.

IV.2 The normal limit

It is well known that the distribution of a random variable given by a superposition of uncorrelated pulses approaches a normally distributed process in the normal limit $\gamma\to\infty$.Rice (1944); Garcia (2012); Garcia et al. (2016b) In the case of a one-sided exponential pulse shape, $\lambda\to 0$, the rescaled process $\widetilde{\Phi}$ is in the normal limit characterized by a Gaussian PDF and an exponential auto-correlation function. The statistical properties of a normally distributed random process are completely described by its PDF and auto-correlation function, and the process is thus statistically identical to any process with a standard normal distribution and exponential auto-correlation function generated by different means.

Much work has been done to elucidate the level crossing statistics and time above or below threshold for an Ornstein–Uhlenbeck (OU) process.Yi (2010); Alili, Patie, and Pedersen (2005) We give the OU-process in our notation as

[TABLE]

where $\text{d}W$ is a standard Wiener process and the initial value is given by $X(0)=x_{0}>0$. This process is normally distributed with mean $\langle{X}\rangle(t)=x_{0}\exp(-t/\tau_{\text{d}})$ and variance $X_{\text{rms}}^{2}(t)=1-\exp(2t/\tau_{\text{d}})$, and has an exponential auto-correlation function with e-folding time $\tau_{\text{d}}$. In Appendix E, it is shown that the moments of $\widetilde{\Phi}(t)$ can be written in the same way. $X(t)$ is thus identical to $\widetilde{\Phi}(t)$ with the conditions described above.

For the case of zero threshold, LABEL:yi-2010 gives the PDF of time above threshold ${\Delta T}$ (and a discussion of relevant references) as

[TABLE]

The initial value $x_{0}>0$ can be identified as the normalized value of the signal below the threshold plus the value of the pulse which brought the signal above the threshold. If the un-normalized initial value is $\Phi_{0}$, the relationship between $x_{0}$ and $\Phi_{0}$ is

[TABLE]

We show in Appendix F that for a threshold value $L$, $\Phi_{0}$ has a truncated exponential distribution

[TABLE]

With $x_{0}$ given above and the threshold being the zero crossing of $\widetilde{\Phi}$, which corresponds to crossing the mean value of $\Phi$ (as was also commented in LABEL:yi-2010, crossing any stationary mean value is statistically equivalent to crossing the stationary mean value [math]), we have

[TABLE]

Thus the full PDF of ${\Delta T}$ can be shown to be

[TABLE]

Changing variables to $\tau=\gamma\left[\exp(2{\Delta T}/\tau_{\text{d}})-1\right]/2$, this PDF can be written more compactly as

[TABLE]

which is independent of $\gamma$. The mean value of the excess time ${\Delta T}$ can also be found,

[TABLE]

where $\operatorname{erfi}(x)=-i\operatorname{erf}(ix)$ and $\operatorname{Ei}(x)$ is the exponential integral.Olver et al. We note that

[TABLE]

in agreement with the result in Eq. (43) for the threshold $\widetilde{\Phi}=0$. We can also find that

[TABLE]

suggesting that $P_{\Delta T}({\Delta T})$ has an exponential tail for large ${\Delta T}$.

V Monte–Carlo studies

In this section, we will investigate some properties of excess time statistics for which we do not have analytical results. Firstly, we will employ a Monte–Carlo approach for investigating the PDF of ${\Delta T}$ for general $\gamma$. Secondly, the question of how quickly the rate of threshold crossings converges to the analytical value will be investigated.

V.1 PDF of excess times

The PDF of excess times in the case where pulse overlap can be neglected was investigated in Sec. IV, and the special case of crossings over the mean value in the case $\gamma\gg 1$ was discussed in Sec. IV.2. The search for an expression for the distribution of time until a process crosses a given threshold is not new, and is frequently referred to as the distribution of first passage time. The Laplace transform for the time until a FPP crosses a given threshold from below is given in Refs. Tsurui and Osaki, 1976; Perry, Stadje, and Zacks, 2001; Novikov et al., 2005. The related problem of the first passage time for an Ornstein–Uhlenbeck process has been investigated by for example Refs. Madec and Japhet, 2004; Alili, Patie, and Pedersen, 2005; Yi, 2010.

To the best of the authors knowledge, there is no closed form expression for the distribution of times above threshold, and discussion of numerically computed PDFs are rare. In this section we therefore present a simulation study of the complementary CDF of ${\Delta T}$ in the case of a one-sided exponential pulse shape, $\lambda=0$. Determining the PDF of times above threshold by simulating the process, with some examples presented in Fig. 1, and estimating $P_{\Delta T}({\Delta T})$ from the realization is computationally prohibitive, in particular for large $\gamma$ and threshold values. We will therefore use a more direct algorithm, according to the following procedure:

At time $t=0$, a pulse arrives, taking the signal from below to above the threshold $L$. The signal takes on the value $\Phi(0)>L$ immediately after the pulse arrival. How $\Phi(0)$ is computed is discussed below. 2. 2.

This arrival ensures that the signal at least spends a time $t_{0}=\tau_{\text{d}}\ln(\Phi(0)/L)$ above the threshold, which is the excess time in the case of no other pulse arrivals in this time interval. 3. 3.

Draw a waiting time $\tau_{1}$ from the exponential waiting time distribution. If $\tau_{1}>t_{0}$, the signal decays below the threshold before the next pulse arrives, and the excess time is $t_{0}$. If $\tau_{1}<t_{0}$, the signal now spends a time

[TABLE]

above the threshold, where $A_{1}$ is the exponentially distributed amplitude associated with the pulse arriving at $\tau_{1}$. 4. 4.

Draw a new waiting time $\tau_{2}$, and compare $\tau_{1}+\tau_{2}$ to $t_{1}$. If $\tau_{1}+\tau_{2}<t_{1}$, make $t_{2}$ in the same way as above. 5. 5.

Continue until the sum of the waiting times would place the arrival of the $n$’th pulse after the signal has decayed below the threshold. The time above threshold is then $t_{n-1}$ for this iteration. 6. 6.

Repeat as often as necessary, and estimate $P_{\Delta T}({\Delta T})$ from all times above threshold found in steps 1-5 above.

Step 1 requires calculating $\Phi(0)$, which consists of two parts. Assume a stationary FPP takes the value $\Phi_{-}<L$ just before time $t=0$. A pulse with amplitude $A_{0}$ arrives and takes the signal above the threshold, $\Phi(0)=\Phi_{-}+A_{0}>L$. It is shown in Appendix F that the PDF of $\Phi(0)$ is

[TABLE]

independent of the intermittency parameter $\gamma$. Samples from this distribution are readily drawn using inverse random sampling. The algorithm presented above is reasonably fast, and allows for accurate computation of the empirical CDF.

In Fig. 10, we present plots of $1-C_{{\Delta T}/\langle{{\Delta T}}\rangle}({\Delta T}/\langle{{\Delta T}}\rangle)$ as a function of ${\Delta T}/\langle{{\Delta T}}\rangle$ for $\gamma\in\{10^{-3},10^{-2},10^{-1},1,10,10^{2}\}$ and various values of the rescaled threshold value $\widetilde{\Phi}$. The full lines give the empirical complementary CDF for $10^{7}$ excess time simulations. In Figs. 10(a), 10(b) and 10(c), the dashed lines give the complementary CDF for ${\Delta T}$ in the limit $\gamma\to 0$ given by Eq. (51). This expression matches the simulated results for short times above threshold, but underestimate the result for longer excess times. This is due to the fact that for small but finite $\gamma$, pulse overlap is significant enough to make longer times above threshold more likely. There is a clear bump in the complementary CDF for $\gamma=10^{-3}$, which is also visible for $\gamma=10^{-2}$. This bump signifies the departure of the simulated distribution from the analytic result in the limit $\gamma\to 0$, and is due to the breakdown of the assumption of negligible pulse overlap, caused by the arrival of a second pulse after the original one.

In Figs. 10(e) and 10(f), the dashed line represents the complementary CDF in the case of $\gamma\gg 1$, from Eq. (IV.2). This is calculated from $1-C_{\tau}(\tau({\Delta T}/\left<\Delta T\right>))$ given by Eq. (59), where $\tau({\Delta T}/\left<\Delta T\right>)=\gamma[\exp(2\left<\Delta T\right>{\Delta T}/\tau_{\text{d}})-1]/2$ and $\left<\Delta T\right>$ is taken from Eq. (IV.2). The $\gamma$-values of the respective figures have been used in this calculation. It is evident that the simulated PDF approaches the analytical one in the limit $\gamma\to\infty$.

For $\gamma<1$, the distribution is concave (on the logarithmic scale) and transitions to a convex distribution for $\gamma>1$. As seen in Fig. 10(d), the distribution for $\gamma=1$ is an exponential distribution for all values of the threshold level. Exponential tails for large ${\Delta T}$ are seen for $\gamma=10^{-1}$ and larger. In the limit $\gamma\to\infty$, this was already suggested by Eq. (62). The exponential tails are not a universal trait of this PDF; the Gompertz distribution for ${\Delta T}$ in the case $\gamma\to 0$ decays as $\exp(-\exp({\Delta T}))$.

V.2 Convergence of excess statistics

In this section, we will quantify how fast the rate of level crossings converges to the analytical value. The process is as follows:

Choose the duration $T/\tau_{\text{d}}$ of a realization of the process, the imtermittency parameter $\gamma$ and the pulse asymmetry parameter $\lambda$. Generate a realization of the process. 2. 2.

Choose $N=200$ threshold values $L_{n},n=1,2,...,199,200$ evenly spaced between $\widetilde{\Phi}=2$ and $\widetilde{\Phi}=10$, and estimate the rate of level crossings $\widehat{X}_{n}$ for each $L_{n}$. 3. 3.

Find the mean squared logarithmic error $\mathcal{E}=\frac{1}{N}\sum_{n=1}^{N}\left[\ln(\widehat{X}_{n})-\ln(X(L_{n}))\right]^{2}$. We use the logarithmic error instead of the linear error since the rate of threshold crossings falls exponentially with increasing threshold for large threshold values, and we wish to emphasize large threshold values. 4. 4.

Repeat as often as necessary to estimate the mean of $\mathcal{E}$ for different $T/\tau_{\text{d}}$, $\gamma$ and $\lambda$.

In Fig. 11, we present the estimated mean squared error of synthetic data for $\gamma\in\{10^{-1},1,10,10^{2}\}$ and $\lambda\in\{10^{-1},1/2\}$. The algorithm described above was repeated 100 times for each set of parameters. In all cases, the mean squared error is inversely proportional to $T/\tau_{\text{d}}$. In Figs. 11(a)-11(d), we see that the error for $\lambda=10^{-1}$ is larger than the error for $\lambda=1/2$ in all cases. This is most likely a side effect of the algorithm used, where the pulses are forced to arrive at integer multiples of $\triangle_{t}$. This introduces a slight bias in the synthetic data, which becomes larger the more asymmetric the pulse shape is. It is also evident from Figs. 11(e) and 11(f) that the error decreases with increasing $\gamma$. Higher $\gamma$ for equal $T/\tau_{\text{d}}$ signifies more pulses, which may lead to quicker convergence, as the samples $\{A_{k}\}_{k=1}^{K}$ and $\{t_{k}\}_{k=1}^{K}$ more closely reflect their underlying distributions for larger K.

VI Discussion and conclusion

In this contribution, a reference model for intermittent fluctuations in physical systems has been investigated. The model consists of a super-position of uncorrelated pulses with a fixed, exponential pulse shape and exponentially distributed pulse amplitudes arriving according to a Poisson process. The PDF and moments of the process were reviewed, and the moments and distribution of its derivative were discussed. The joint PDF between the process and its derivative was derived and used to obtain predictions for level crossing rates and average excess times for fluctuations above a given threshold level. These predictions depend on two model parameters, the intermittency parameter $\gamma$ and the pulse shape asymmetry parameter $\lambda$. It was shown that the functional shape of the rate of level crossings with the threshold level is strongly dependent on the intermittency parameter $\gamma$ of the process, while the functional shape of the average excess time varies little with the parameter $\gamma$. In both cases, the functional shape is independent of $\lambda$, as this parameter only appears in the pre-factor. The limit of $\lambda\to 0$ was considered, and was shown to be in agreement with previous works using different methods.Biermé and Desolneux (2012); Daly and Porporato (2010) The limits of highly intermittent signals as well as the normal limit were investigated. The normal limit was shown to be in agreement with the well-known Rice’s formulaRice (1945) for $0<\lambda<1$, and was shown to have qualitatively different behaviour for $\lambda\in\{0,1\}$.

The PDF of the time the stochastic process spends above the threshold was found analytically in both the limit of strong intermittency for general threshold level and in the normal limit for threshold equal to the mean value, adapted from studies of Ornstein-Uhlenbeck-processes.Yi (2010); Alili, Patie, and Pedersen (2005) In the strong intermittency limit, the time above threshold was shown to be Gompertz distributed. Both limits were in agreement with a Monte-Carlo study of synthetically generated time series, and the shape of the complementary CDF of time above threshold from synthetic data was presented for various values of the intermittency parameter $\gamma$. In order to investigate the convergence of the rate of level crossings to the analytical expression, another Monte-Carlo study was performed. The convergence was shown to be proportional to $\tau_{\text{d}}/T$.

The model presented here has previously been shown to be a good description of intermittent fluctuations in the boundary region of magnetically confined plasmas,Garcia et al. (2013); Garcia, Horacek, and Pitts (2015); Kube et al. (2016); Theodorsen et al. (2016); Garcia et al. (2016a); Kube et al. (shed); Theodorsen et al. (sion) and the rate of level crossings has compared favorably to large-amplitude fluctuations in the SOL.Garcia et al. (2016a); Theodorsen et al. (sion) In comparing the model to experimental data, the results presented here provide two major improvements over the classical Rice’s formula in the case of intermittent fluctuations. Firstly, any discrepancy between the normal limit for excess time statistics and measurement data has previously been interpreted as a signature of intermittency in the process. The formulas derived here quantifies the level of intermittency by the model parameters $\lambda$ and $\gamma$. Secondly, Rice’s formula requires the rms-value of the derivative of the signal, which is difficult if not impossible to reliably estimate for discretely sampled data containing measurement noise. In contrast, estimates for $\lambda$ and $\gamma$ can be found from the signal using the lowest order moments of $\Phi$ and its correlation function or frequency power spectrum.Theodorsen et al. (2016); Kube et al. (2016) The variation of rate of level crossings and excess time statistics shows that the average time above threshold varies little with the intermittency parameter, suggesting that the rate of level crossings might be a more useful tool in comparing the model to experimental data in order to assess intermittency effects, although the time above threshold may be the more relevant statistic for applications in terms of failure or stability.

Even though the total time above a given threshold level may be the same for realizations of two different intermittent processes, this can be realized through either many short bursts or few but long lasting bursts events. This may have profound implications for systems where long lasting, large amplitude events can lead to severe damaging while the system can recover from the impacts of shorter burst events, depending on their frequency of occurrence. An example would be intermittent plasma-wall interactions in magnetically confined plasmas.Sato, Pécseli, and Trulsen (2012); Fattorini et al. (2012) Thus accurately predicting the rate of level crossings and average excess times for an intermittent process is of considerable interest to statistical modelling of fluctuations in the boundary region of magnetically confined plasmas. In future work, the novel predictions presented here will be further compared to experimental measurement data from the scrape-off layer of magnetically confined plasmas.

Acknowledgements

This work was supported with financial subvention from the Research Council of Norway under grant 240510/F20. Discussions with M. Rypdal are gratefully acknowledged. The authors acknowledge the generous hospitality of the MIT Plasma Science and Fusion Center where this work was conducted.

Appendix A The joint PDF of $\Phi$ and $\Theta$ in the normal limit

We will here demonstrate that the joint PDF of $\Phi$ and $\Theta$ given by Eq. (29) is a joint normal distribution with zero correlation coefficient in the limit $\gamma\to\infty$. We begin by changing variables to the normalized

[TABLE]

where the moments of $\Phi$ and $\Theta$ are given in Eqs. (18) and (22), respectively. Then we have

[TABLE]

By Stirling’s formula, both fractions in the pre-factor are equal to $1/\sqrt{2\pi}$. Using the notation

[TABLE]

we have that $\widetilde{\Phi}=\lambda\alpha+(1-\lambda)\beta$ and

[TABLE]

Thus, in the limit of $\gamma\to\infty$, the joint PDF of $\widetilde{\Phi}$ and $\widetilde{\Theta}$ approaches a joint normal distribution of two independent variables.

Appendix B An integral connected to the rate of threshold crossings

In Eq. (34), the integral

[TABLE]

was presented. In Sec. II.2, this PDF was shown to be a convolution between a gamma distribution $P_{+}(\Theta)$ over positive values of $\Theta$ with shape parameter $\gamma\lambda$ and scale parameter $\left<A\right>/2(1-\lambda)$ and a gamma distribution $P_{-}(\Theta)$ over negative values of $\Theta$ with shape parameter $\gamma(1-\lambda)$ and scale parameter $\left<A\right>/2\lambda$. The PDF of $\Theta$ is therefore

[TABLE]

where the integration limits are due to the domain of non-zero values for the gamma functions. Inserting this into Eq. (71) and exchanging the order of integration lets us compute the integral. The result is

[TABLE]

where ${}_{2}F_{1}(a,b;c;z)$ is the hypergeometric function.Olver et al. This expression contains the prefactor of Eq. (33), but the terms inside the curly brackets give this expression a very different behavior.

Appendix C The rate of level crossings in the normal limit

In this Appendix, we derive a necessary result in order to go from Eq. (38) to Eq. (40). The derivation is analogous to Eq. (A):

[TABLE]

Appendix D Upper incomplete Gamma function to error function

In this Appendix, an asymptotic limit of the upper incomplete Gamma function is derived. We have

[TABLE]

By substituting $u=(t-\gamma)/\sqrt{\gamma}$ and using that $\gamma\Gamma(\gamma)=\Gamma(\gamma+1)$, this expression becomes

[TABLE]

The fraction is $1/\sqrt{2\pi}$ by Stirling’s formula, and using Eq. (C), we have that

[TABLE]

This result is used to show the equivalence between the average excess time in Eq. (31) and Rice’s result in Eq. (6) in the limit $\gamma\to\infty$.

Appendix E Time-dependent moments of the FPP

In this Appendix, we derive the first two time-dependent moments of a normalized FPP. In LABEL:tsurui-1976, the time-dependent characteristic function of the FPP is given as

[TABLE]

We explicitly demand a pulse arriving at $t=0$ with value $\Phi_{0}$. Thus we can write a modified version of the FPP as

[TABLE]

We assume $\Phi_{0}$ is given, such that $\Psi(t)$ has the characteristic function

[TABLE]

The two first moments of $\Psi$ are $\langle{\Psi}\rangle(t)=\gamma\left<A\right>(1-\exp(-t/\tau_{\text{d}}))+\Phi_{0}\exp(-t/\tau_{\text{d}})$ and $\Psi_{\text{rms}}^{2}(t)=\gamma\left<A\right>^{2}(1-\exp(-2t/\tau_{\text{d}}))$, and the stationary moments are $\langle{\Psi}\rangle=\left<{\Psi}\right>(t\to\infty)=\gamma\left<A\right>$ and $\Psi_{\text{rms}}^{2}=\Psi_{\text{rms}}^{2}(t\to\infty)=\gamma\left<A\right>^{2}$. Normalizing $\Psi$ by

[TABLE]

it is straightforward to show that

[TABLE]

Writing $\Phi_{0}$ as in Eq. (55), this equation can be written as

[TABLE]

whose first two moments are $\langle{\widetilde{\Psi}}\rangle(t)=x_{0}\exp(-t/\tau_{\text{d}})$ and $\widetilde{\Psi}_{\text{rms}}^{2}(t)=1-\exp(-2t/\tau_{\text{d}})$. These moments are independent of $\gamma$ and are equal to the moments of the OU process in Eq. (53).

Appendix F The truncated exponential distribution

In this Appendix, we derive the result presented in Eq. (64). Consider a stationary stochastic process $\Phi$ consisting of a super-position of uncorrelated random pulses. Assume the pulses have a positive jump at the arrival time, and only are non-zero after the arrival time. Just before $t=0$, the value of $\Phi$ is below the threshold $L$;

[TABLE]

A pulse with amplitude $A$ arrives at $t=0$, taking the signal above the threshold:

[TABLE]

It is assumed that $A$ is exponentially distributed with mean value $\left<A\right>$. $\Phi_{-}$ can in principle have an arbitrary distribution. The distribution of $\Phi_{0}$ is then found from integrating the joint distribution of $A$ and $\Phi_{-}$ over the region $A+\Phi_{-}<\Phi_{0}$, under the conditions in Eqs. (84) and (85). Since $A$ and $\Phi_{-}$ are independent, we have

[TABLE]

where $C_{\Phi_{0}}$ is the CDF of $\Phi_{0}$ and

[TABLE]

and

[TABLE]

The truncated distributions are calculated by using the method given in LABEL:stark-psrpe. This gives

[TABLE]

The derivative with respect to $\Phi_{0}$ can be brought inside the first integral, and we have that

[TABLE]

where we have used that $A$ is exponentially distributed. This does not depend on $\Phi_{-}$, so we have

[TABLE]

The fraction is unity by definition, and the distribution is

[TABLE]

Thus, if we assume exponentially distributed pulse amplitudes and pulses with jumps as described above, the distribution of the stationary process $\Phi$ plays no role for the distribution of the jumps above the threshold. In particular, for a FPP with $\lambda\to 0$ and $\gamma\to\infty$, the process is normally distributed, but the value of the signal just after the threshold is crossed is exponentially distributed.