Probability distribution functions for intermittent scrape-off layer plasma fluctuations

# Probability distribution functions for intermittent scrape-off layer plasma fluctuations

A. Theodorsen Department of Physics and Technology, UiT The Arctic University of Norway, N-9037 Tromsø, Norway    O. E. Garcia Department of Physics and Technology, UiT The Arctic University of Norway, N-9037 Tromsø, Norway
July 14, 2019
###### Abstract

A stochastic model for intermittent fluctuations in the scrape-off layer of magnetically confined plasmas has been constructed based on a super-position of uncorrelated pulses arriving according to a Poisson process. In the most common applications of the model, the pulse amplitudes are assumed exponentially distributed, supported by conditional averaging of large-amplitude fluctuations in experimental measurement data. This basic assumption has two potential limitations. First, statistical analysis of measurement data using conditional averaging only reveals the tail of the amplitude distribution to be exponentially distributed. Second, exponentially distributed amplitudes leads to a positive definite signal which cannot capture fluctuations in for example electric potential and radial velocity. Assuming pulse amplitudes which are not positive definite often make finding a closed form for the probability density function difficult, even if the characteristic function remains relatively simple. Thus estimating model parameters requires an approach based on the characteristic function, not the probability density function. In this contribution, the effect of changing the amplitude distribution on the moments, probability density function and characteristic function of the process is investigated and a parameter estimation method using the empirical distribution function is presented and tested on synthetically generated data. This proves valuable for describing intermittent fluctuations in the boundary region of magnetized plasmas.

## I Introduction

Radial propagation of filamentary structures is the main contributor to cross-field transport of particles and heat in the scrape-off layer (SOL) of magnetically confined plasmas LaBombard et al. (2001); Krasheninnikov (2001); Rudakov et al. (2002); D’Ippolito et al. (2004); Rudakov et al. (2005); Zweben et al. (2007); Krasheninnikov et al. (2008); GARCIA (2009); D’Ippolito et al. (2011); Kube et al. (2013); Zweben et al. (2015). This turbulence-driven transport results in broad plasma profiles and enhanced plasma-wall interactions Rudakov et al. (2005); Pitts et al. (2005); Lipschultz et al. (2007); Garcia et al. (2007a, b); Militello et al. (2013); Carralero et al. (2014, 2015); Militello and Omotani (2016a); Walkden et al. (2017a).

Statistical analysis of single-point measurements in the far-SOL of several tokamak experiments reveal skewed and flattened fluctuation probability density functions (PDFs), exponential auto-correlation functions and Lorenzian power spectra for positive definite variables such as ion saturation current, electron density and temperature, and gas puff imaging (GPI) intensity signals Garcia et al. (2007b, 2013); Militello et al. (2013); Kube et al. (2016); Garcia et al. (2017); Theodorsen et al. (2016, 2017a); Kube et al. (); Walkden et al. (2017b, a). Conditional averaging of large-amplitude fluctuations show that large structures exhibit fast exponential growth and slower exponential decay, with exponentially distributed peak amplitudes and exponentially distributed waiting times between them Garcia et al. (2007b, 2013); Kube et al. (2016); Garcia et al. (2017); Theodorsen et al. (2016); Walkden et al. (2017b); Kube et al. (). Measurements of the radial velocity is shown to have PDFs with exponential tails which are nearly symmetric around the mean value Theodorsen et al. (2016); Kube et al. (). Previously, PDFs with exponential tails have been investigated using the so-called instanton method Kim and Anderson (2008); Anderson and Xanthopoulos (2010).

In order to systematize and unify these observations, a well-known reference model for intermittent fluctuations Rice (1944, 1945); Fesce (1986); Kristensen et al. (1991); Jang (2004); Daly and Porporato (2006); Narasimha et al. (2007); Elter et al. (2015) has been introduced for SOL plasma fluctuations Garcia (2012); Garcia et al. (2016); Militello and Omotani (2016b). This model, called a shot noise process or filtered Poisson process (FPP), consists of a super-position of independent and identical pulse shapes with randomly distributed amplitudes, arriving according to a Poisson process. The predictions of this FPP have been shown to be in excellent agreement with experimental measurements of PDFs, auto-correlation functions and frequency power spectra, conditional averaging, and higher order statistics such as threshold level crossings and excess time statistics Garcia et al. (2013); Kube et al. (2016); Theodorsen et al. (2016); Garcia et al. (2017); Theodorsen et al. (2017a); Walkden et al. (2017b); Kube et al. ().

This model provides a framework for comparing measurements of SOL fluctuations. For example, it has been demonstrated from GPI data in Alcator C-Mod that far-SOL fluctuations have highly skewed and flattened PDFs, while PDFs close to the separatrix more closely resemble normal distributions. Garcia et al. (2013); Theodorsen et al. (2017a) At the same time, Lorenzian power spectra with the same time scale are observed at all radial positions Theodorsen et al. (2017a). Interpreting the PDFs by model, blobs are numerous and close together in the near SOL, while they are further apart in the far SOL. The blobs retain their basic shape while traveling through the SOL, however, as indicated by the universality of the spectra Theodorsen et al. (2017a). By comparing PDFs and power spectra, and by estimating model parameters, measurements in the SOL of different fusion experiments, in various confinement modes and for a range of plasma parameters can be compared. We note, however, that due to the time invariance of the Poisson process, the model inherently only describes statistically stationary turbulence in the SOL.

Previous theoretical work has revealed the convergence of the lowest order moments for the process Kube and Garcia (2015), extended the model to include additive noise Theodorsen et al. (2017b), revealed the auto-correlation and power spectrum for different pulse shapes Garcia and Theodorsen (2017a) and for randomly distributed pulse durations Garcia and Theodorsen (2017b) and discussed the rate of threshold crossings and average time above a given threshold Theodorsen and Garcia (2016). It has also been demonstrated that radial propagation of filament structures with these statistical properties results in exponential profiles in the SOL, consistent with experimental observations Garcia et al. (2016); Militello and Omotani (2016b). In this contribution, we focus on the PDF, the characteristic function and the lowest order moments of the process for various relevant amplitude distributions.

While conditional averaging has demonstrated exponentially distributed amplitudes for large-amplitude fluctuations, the statistical properties of small amplitudes has not been revealed. Assuming a positive definite time series (as is the case for ion saturation current, electron density, electron temperature or GPI intensity), exponentially distributed amplitudes is an obvious assumption. Another candidate is the Gamma distribution with shape parameter . This distribution is unimodal and decays exponentially for large amplitudes, but has vanishing probability for amplitudes approaching zero. In this contribution we will compare the distributions for the FPP given exponentially and Gamma distributed amplitudes with shape parameter , since this captures the essential differences between the two distributions while allowing for analytical treatment of the PDF and the characteristic function of the FPP.

While the ion saturation current, electron density and temperature, and GPI intensity all are positive definite variables, electric potential and radial velocity are not. Thus, in order to correctly model fluctuations in these quantities, amplitude distributions with non-zero probability for negative amplitudes are required. The asymmetric Laplace distribution fulfills this requirement while still having exponential tails. The PDF of the FPP with symmetric Laplace distributed amplitudes has been favorably compared to measurements of radial velocity of filaments in the SOL Theodorsen et al. (2016); Kube et al. (). In the asymmetric case, however, the resulting FPP does not have a closed form solution for the PDF. Thus methods for estimating model parameters based on or requiring the PDF of the process are not applicable. However, the characteristic function for the model can still be found in closed form. This allows for a method based on the empirical characteristic function, which is general enough to allow for any asymmetry in the Laplace distributed amplitudes (of which the exponential distribution is a special case) and noise level. We present this method and its application to the FPP.

## Ii Moments of the filtered Poisson process

In this section, the FPP is introduced as a model for SOL fluctuations. The general form of the moments are presented for two-sided, exponential pulse shapes and three different pulse amplitude distributions: the exponential distribution, the Gamma distribution with shape parameter and the asymmetric Laplace distribution. We also investigate how normally distributed additive noise affects the moments of the FPP.

The FPP is given by a super-position of identical pulse shapes with randomly distributed amplitudes arriving at times restricted to the range . All random variables are assumed independent. The pulses arrive according to a Poisson process with intensity , where is the average waiting time between pulses. Thus the arrival times are uniformly distributed on the interval and the waiting time between pulses is exponentially distributed with mean value .

We express the FPP as Rice (1944, 1945); Parzen (1999); Pécseli (2000)

 ΦK(t)=K(T)∑k=1Akφ(t−skτd), (1)

where is the fixed pulse duration time. In general, the duration times may be randomly distributed, but only the mean value of this variable (that is, ) plays a role for the moments and distribution of , see App. A. Thus, for simplicity of notation and without loss of generality, we will in the following consider a constant duration time. The pulse shape , where is a unitless variable, is normalized according to

 ∫∞−∞dθ|φ(θ)|=1. (2)

Here and in the following, only the two-sided exponential pulse shape will be considered. This pulse shape has an asymmetry parameter , and is given by

 φ(θ)={exp(θ/λ),θ<0,exp(−θ/(1−λ)),θ≥0. (3)

We define the integral of the ’th power of the pulse function, which for the two-sided exponential pulse is independent of the asymmetry parameter ,

 In=∫∞−∞dθ[φ(θ)]n=1n. (4)

While in principle could be randomly distributed, as discussed in App. A, we will in this contribution assume all pulses to be identical, with the same, fixed and .

Under the assumptions given above, the characteristic function of the FPP has been derived in App. A, and is given by Eq. (69). Inserting the integral of the ’th power of the pulse function given in Eq. (4) into Eq. (69), the logarithm of the characteristic function of the FPP is given by the sum

 lnCΦ(u)=γ∞∑n=1(iu)nn!n⟨An⟩, (5)

where is the so-called intermittency parameter of the process. This parameter determines the degree of pulse overlap, and thereby the intermittency of the process. For low , each pulse duration is short compared to the average time between pulses, and the process is strongly intermittent. For high , many pulses arrive in the duration of one pulse event and pulse overlap becomes significant.

The cumulants of the process are given by the coefficients in the expansion of the logarithm of the characteristic function,

 lnCΦ(u)=∞∑n=1κn(iu)nn!, (6)

which according to Eq. (28) are given by

 κn=γn⟨An⟩. (7)

The mean value of the process is , the variance is , where rms denotes the root mean square value, and the skewness and flatness moments are related to the cumulants by Rice (1944); Garcia (2012)

 SΦ =κ3κ3/22, (8a) FΦ =3+κ4κ22. (8b)

According to Eq. (7), each cumulant is proportional to for any amplitude distribution. Thus, the mean value is proportional to , the rms-value is proportional to , the skewness is proportional to and the flatness to . For increasing , any FPP will have mean and rms tending to infinity (for finite and ) and vanishing skewness and flatness. For the FPP there is a universal parabolic relationship between and , independent of the intermittency parameter Garcia (2012); Garcia et al. (2016)

 FΦ=3+κ2κ4κ23S2Φ=3+I2I4I23⟨A2⟩⟨A4⟩⟨A3⟩2S2Φ=3+98⟨A2⟩⟨A4⟩⟨A3⟩2S2Φ. (9)

The physical basis for a parabolic relationship between skewness and kurtosis has been explored previously Sattin et al. (2009); Bergsaker et al. (2015). We note that while many relationships between skewness and kurtosis based on Eq. (8) are possible, only the one presented in Eq. (9) is independent of the intermittency parameter .

We will consider 3 different amplitude distributions for the FPP; the exponential distribution, the Gamma distribution with shape parameter and the asymmetric Laplace distribution. These all give closed form expressions for the characteristic function of . The PDF of the exponential distribution has a finite value for and is monotonically decreasing. The Gamma distribution with shape parameter is unimodal and tends to for . Since both have exponential tails, comparing the distribution of the FPP with exponentially and Gamma distributed amplitudes will highlight the importance of small-amplitude pulses while keeping the effect of large-amplitude pulses equal. The Laplace distribution allows for both positive and negative values of , whereas the exponentially and Gamma distributed amplitudes are strictly positive. Thus the Laplace distribution is the only one of these capable of describing measurement data which is not positive definite.

### ii.1 Exponentially distributed amplitudes

The exponential distribution is a one parameter distribution with scale parameter . The distribution and its moments are given by

 PA(A;α) =1αexp(−A/α),A>0, (10a) ⟨An⟩ =αnn!, (10b)

for integer values of .

For exponentially distributed amplitudes, the first four moments of are given by Garcia (2012); Garcia et al. (2016)

 ⟨Φ⟩ =γα, (11a) Φ2rms =γα2, (11b) SΦ =2γ1/2, (11c) FΦ =6γ+3, (11d)

and we have the parabolic relationship between the skewness and flatness moments,

 FΦ=3+32S2Φ, (12)

where the pre-factor is simply , as the scale parameter is cancelled out.

### ii.2 Gamma distributed amplitudes

The Gamma distribution has a shape parameter and a scale parameter . It is given by

 PA(A;α,β) =Aβ−1αβΓ(β)exp(−A/α),A>0, (13a) ⟨An⟩ =αnΓ(β+n)Γ(β). (13b)

For , this is equivalent to the exponential distribution with scale parameter . For , we have

 PA(A;α,β) =Aα2exp(−A/α),A>0, (14a) ⟨An⟩ =αnΓ(2+n). (14b)

For large amplitudes , the Gamma distribution has an exponential tail with the same decay rate as the exponential distribution with equal . The moments are not equal, however, as , giving a factor more for all moments in the case of Gamma distributed pulse amplitudes.

For Gamma distributed amplitudes with general shape parameter , the first four moments are given by

 ⟨Φ⟩ =γαβ, (15a) Φ2rms =γα2β(β+1)2, (15b) SΦ =23/23γ1/2β+2[β(β+1)]1/2, (15c) FΦ =1γ(β+2)(β+3)β(β+1)+3. (15d)

The parabolic relationship between the skewness and flatness moments depends on the shape parameter for the amplitude distribution,

 FΦ=3+98β+3β+2S2Φ. (16)

Setting gives the same moments and parabolic relation as for the exponentially distributed amplitudes. For the special case , the moments are

 ⟨Φ⟩ =2γα, (17a) Φ2rms =3γα2, (17b) SΦ =833/2γ1/2, (17c) FΦ =103γ+3. (17d)

and the parabolic relationship simplifies to

 FΦ=3+4532S2Φ. (18)

This relationship is very close to the case of exponentially distributed amplitudes, Eq. (12), (it would be equal for the prefactor ).

### ii.3 Asymmetrically Laplace distributed amplitudes

The asymmetric Laplace distribution can be formulated in a few different ways, see for example Refs. Kozubowski and Podgórski, 2000; Reed, 2006. We will use a different formulation which easily admits the exponential distribution as a limiting case. With as a scale parameter and as a shape parameter, we have

 PA(A;α,β) =12α⎧⎪⎨⎪⎩exp(−A2α(1−β)),A>0,exp(A2αβ),A<0, (19a) ⟨An⟩ =(2α)nn![(−1)nβn+1+(1−β)n+1]. (19b)

This distribution is symmetric for , and is equivalent to the exponential distribution in the limit . In the limit , the distribution is an exponential distribution mirrored around , with zero probability for positive -values and finite probability for negative -values.

For Laplace distributed pulse amplitudes, the first four moments are given by

 ⟨Φ⟩ =2γα(1−2β), (20a) Φ2rms =4γα2[β3+(1−β)3], (20b) SΦ =(1−β)4−β4[β3+(1−β)3]3/22γ1/2, (20c) FΦ =[β5+(1−β)5][β3+(1−β)3]26γ+3. (20d)

and the parabolic relationship between skewness and flatness is,

 FΦ=3+32[β5+(1−β)5][β3+(1−β)3][(1−β)4−β4]2S2Φ. (21)

For , we have the same expressions as for the exponentially distributed amplitudes. For , the Laplace distribution is symmetric, so all odd moments of vanish, see Eq. (19b), giving

 ⟨Φ⟩ =0, (22a) Φ2rms =γα2, (22b) SΦ =0, (22c) FΦ =6γ+3. (22d)

In this case, there is no parabolic relationship between skewness and flatness, since . As approaches (from either side), the prefactor in Eq. (21) tends to infinity.

### ii.4 Comparisons

In Fig. 1, realizations of the process are presented for and the amplitude distributions given above. The bottom (blue) lines give realizations with amplitudes distributed according to an exponential distribution. The middle (orange dashed) lines are for Gamma distributed amplitudes with shape parameter , and the top (green dotted) realizations are computed using the Laplace distribution with shape parameter . In all cases, the realizations have been normalized to have zero mean and unit standard deviation, in order to remove the dependency on the scale parameter in the amplitude distribution. In Fig. (a)a, and all processes are strongly intermittent, alternating between periods of activity and inactivity. In Fig. (b)b, , and the large degree of pulse overlap leads to weaker intermittency and makes individual pulses harder to discern. While the signals with exponentially and Gamma distributed amplitudes are easy to separate visually from the Laplace case for , this is not true for the case with , where all realizations look like a random walk around a mean value. Indeed, it can be shown that the PDF of the normalized process approaches a standard normal distribution as , independent of the amplitude distribution and the pulse shape Rice (1944); Pécseli (2000). Visually, it is very difficult to separate the filtered Poisson process with exponentially distributed amplitudes from the one with Gamma distributed amplitudes with shape parameter .

In Fig. 2, the inverse of the prefactor in the parabolic relationship between the skewness and flatness moments as a function of is shown for exponentially distributed amplitudes (blue), Gamma distributed amplitudes (orange dashed) and Laplace distributed amplitudes (green dotted). The inverse is used, since the prefactor itself tends to infinity in the case of Laplace distributed amplitudes for , as discussed above. This prefactor is constant for exponentially distributed amplitudes, since the exponential distribution has no shape parameter. For Gamma distributed amplitudes, the inverse of the prefactor is smaller than for exponentially distributed amplitudes for and larger for , thus the prefactor itself is larger for and smaller for . From Eq. (16) we see that the prefactor approaches as . The FPP with Laplace distributed amplitudes and or has the same prefactor in the parabolic relation as the FPP with exponentially distributed amplitudes. The skewness is equal in magnitude but with different sign for these two cases, giving the same prefactor in the parabolic relationship.

In many applications of the model, there may be some normally distributed additive noise to the process, either as white noise connected to measurements or as noise with the same power spectrum as the FPP, connected to the dynamics. This situation has been explored in detail for exponentially distributed amplitudes and one-sided exponential pulses Theodorsen et al. (2017b). We assume the noise process, denoted by , to be normally distributed and independent of the FPP, and denote the signal with additive noise as

 Ψ=Φ+X. (23)

The basic properties of the distribution of a sum of independent random variables are reviewed in App. 70. For a normal distribution, only the first two cumulants are non-zero. Using that the variance of the process is given by the second cumulant, we define the noise ratio as

 ϵ=X2rmsΦ2rms=κX2κΦ2. (24)

Using this noise ratio and that the cumulants of a sum of independent random variables is the sum of their cumulants, see Eq. (73), we have

 ⟨Ψ⟩ =κΦ1+κX1, (25a) Ψ2rms =(1+ϵ)κΦ2, (25b) SΨ =κΦ3[(1+ϵ)κΦ2]3/2, (25c) FΨ =3+κΦ4[(1+ϵ)κΦ2]2. (25d)

The parabolic relationship between the skewness and flatness moments is

 FΨ=3+κΦ4[(1+ϵ)κΦ2](κΦ3)2S2Ψ. (26)

The effect of noise on the moments is to increase the variance and decrease the skewness and flatness, leading to more closely resembling a normally distributed process than . As increases, the prefactor in the parabolic relationship increases as well.

## Iii Probability distributions

Under the assumptions that we have an FPP with fixed and and independent amplitudes and arrival times, the characteristic function of the FPP has been derived in App. A, and is given in various forms by Eqs. (64), (65) and (69). For the two-sided exponential pulse function, we can split the integral in Eq. (65) into two parts, one over negative values of and one over positive values of , and substitute the integration for the pulse shape in Eq. (3). This gives

 lnCΦ(u)=γu∫0dvCA(v)−1v, (27)

where is the characteristic function of the amplitudes. Inserting the integral of the ’th power of the pulse function given in Eq. (4) into Eq. (69), we can alternatively give the characteristic function of the FPP as the sum

 lnCΦ(u)=γ∞∑n=1(iu)nn!n⟨An⟩, (28)

as in the previous section. The PDF of is given by

 PΦ(Φ)=12π∫∞−∞duexp(−iΦu)CΦ(u). (29)

We will in the following frequently use the normalization

 ˜Φ=Φ−⟨Φ⟩Φrms, (30)

which removes the dependence on the scale parameter in the amplitude distribution. The PDF of is in general

 P˜Φ(˜Φ)=ΦrmsPΦ(Φ%rms˜Φ+⟨Φ⟩), (31)

while its characteristic function is

 C˜Φ(u)=exp(−i⟨Φ⟩Φrmsu)CΦ(uΦrms). (32)

### iii.1 Exponential amplitude distribution

In the case of exponentially distributed amplitudes, it is well known Bondesson (1982); Pécseli (2000); Daly and Porporato (2006, 2010); Garcia (2012) that the distribution of is a gamma distribution with shape parameter and scale parameter ,

 PΦ(Φ) =Φγ−1αγΓ(γ)exp(−Φα),Φ>0, (33a) CΦ(u) =(1−iαu)−γ. (33b)

Using the normalization in Eq. (30), the dependence explicit on disappears, and the distribution becomes

 P˜Φ(˜Φ) =γγ/2Γ(γ)(˜Φ+γ1/2)γ−1exp(−γ1/2˜Φ−γ), (34a) C˜Φ(u) =exp(−i√γu)(1−iu√γ)−γ. (34b)

In Fig. 3, the effect of the parameter is illustrated by presenting the complementary cumulative distribution function (cCDF) of for various parameter values. The cCDF at a given function value can be interpreted as the probability that the random variable takes the value or a larger value. It can also be interpreted as the fraction of time a signal spends above a threshold value . is positive definite, so the lowest possible value for is . For small , the cCDF falls off slowly with increasing , indicating high probability of large amplitude fluctuations. As increases, the probability of large values of decreases, as the signal transitions from an intermittent signal to one resembling random motion around a mean value. As stated earlier, for , the process approaches a standard normal distribution, presented by the black diamonds in Fig. 3. For the -values presented here, the Gamma-distributed signals all have higher probability of fluctuations with amplitude compared to a normally distributed signal, highlighting the importance of intermittency for threshold phenomena.

### iii.2 Gamma amplitude distribution

Using gamma distributed amplitudes, we have

 lnCΦ(u)=γβiαu3F2[111+β22;iαu], (35)

where is the generalized hypergeometric function Olver et al. (). For , this simplifies to the characteristic function for an exponential amplitude distribution, as discussed above. For , the characteristic function is

 lnCΦ(u)=γ[11−iαu−1−ln(1−iαu)], (36)

giving

 CΦ(u)=exp(−γ)exp(γ1−iαu)(1−iαu)−γ. (37)

In App. C, the corresponding PDF is shown to be Daly and Porporato (2006, 2010)

 PΦ(Φ)=1α(Φγα)(γ−1)/2exp(−γ−Φα)Iγ−1(2√γΦα). (38)

Here, is the modified Bessel function of the first kind Olver et al. (). Using the mean and rms-values given by Eq. (17), we have that the distribution of the normalized variable is again independent of , and we have

 P˜Φ(˜Φ) =√3γ(√3γ˜Φ+2)(γ−1)/2exp(−√3γ˜Φ−3γ)Iγ−1⎛⎜⎝2γ ⎷√3γ˜Φ+2⎞⎟⎠, (39a) C˜Φ(u) =exp(−γiu√3γ1+2iu/√3γ1−iu/√3γ)(1−iu√3γ)−γ. (39b)

A comparison between the PDF in Eq. (39a) and the PDF for the case of exponentially distributed amplitudes, given by Eq. (34a), is presented in Fig. 4 for various values of the intermittency parameter. It is evident that the PDF of the FPP with Gamma distributed amplitudes can be very well approximated by the PDF of a FPP with exponentially distributed amplitudes and slightly larger . Thus the PDF seems a poor choice for differentiating Gamma distributed amplitudes with from exponentially distributed amplitudes in a given realization of the process.

### iii.3 Asymmetric Laplace amplitude distribution

With asymmetrically Laplace distributed amplitudes, the characteristic function of the filtered Poisson process is

 CΦ(u)=(1+i2αβu)−γβ(1−i2α(1−β)u)−γ(1−β). (40)

Note that this can be seen as the characteristic function of a sum of two independent gamma distributed variables, one over positive values with shape parameter and scale parameter , and one over negative values with shape parameter and scale parameter . While the PDF is in general not possible to find in closed form, a numerical estimate can be found by noting that the PDF of a sum of independent random variables is the convolution of their respective distributions. Estimating the two gamma distributions and convolving them numerically gives PDFs as illustrated in Fig. 5. In Fig. (a)a, the intermittency parameter is and the shape parameter varies. For , this is a symmetric Laplace distribution. As , the PDF approaches a Gamma distribution. Fig. (b)b shows the distribution for and various values of . As discussed below, one can find the PDF in closed form in this case. The PDF is symmetric around for all . For small , it is sharply peaked and convex around , while for large , it is concave and approaches a normal distribution as increases. For all combinations of and (finite) , the PDF has exponential tails.

In the limit , the Laplace distribution approaches the exponential distribution and the standard Gamma distribution for is recovered. In the case , the Laplace distribution is symmetric and we can find the PDF in closed form. The distribution is Daly and Porporato (2006, 2010); Theodorsen et al. (2016); Garcia and Theodorsen (2017b)

 PΦ(Φ;γ,α,β=1/2) =1π1/2αΓ(γ/2)(|Φ|2α)(γ−1)/2K(γ−1)/2(|Φ|α), (41a) CΦ(u) =(1+α2u2)−γ/2, (41b)

where is the modified Bessel function of the second kind Olver et al. (). The normalized variable has the distribution

 P˜Φ(˜Φ;γ,β=1/2) =γ1/2π1/2Γ(γ/2)(γ1/2|˜Φ|2)(γ−1)/2K(γ−1)/2(γ1/2|˜Φ|), (42a) C˜Φ(u) =(1+u2γ)−γ/2. (42b)

This PDF is presented in Fig. (b)b for various values of the intermittency parameter .

Adding noise to the FPP is straightforward as long as only the characteristic function is considered. Using the FPP with asymmetrically Laplace distributed amplitudes and additive noise as an example, we have

 Ψ(t)=Φ(t)+X(t), (43)

where is a FPP with Laplace distributed amplitudes and is normally distributed noise with vanishing mean and standard deviation . The characteristic function of is the product of the characteristic functions of and , see App. 70. We have

 CΨ(u)=(1+i2βαu)−γβ(1−i2(1−β)αu)−γ(1−β)exp(−12X2rmsu2). (44)

Using the noise parameter from Eq. (24), we have . The moments for the FPP with additive noise from Eq. (25) give and . Normalizing the process by and using Eq. (32), we have

 C˜Ψ(u)=(1+iβu√γ(1+ϵ)B(β))−γβ(1−i(1−β)u√γ(1+ϵ)B(β))−γ(1−β)exp(−ϵu22(1+ϵ)−i√γ(1−2β)u√1+ϵB(β)), (45)

where . Performing the inverse Fourier transform on this expression to get the PDF of does in general not lead to a closed analytical expression. However, it can be done for the FPP with exponentially distributed amplitudes (that is, in the limit ) Theodorsen et al. (2017b).

## Iv Parameter estimation

In Sec. III.2, it was shown that the distribution of an FPP with Gamma distributed amplitudes with shape parameter does not differ significantly from an FPP with exponentially distributed amplitudes and a slightly larger . The exponential distribution is a special case of the asymmetric Laplace distribution, and so we consider the FPP with asymmetrically Laplace distributed amplitudes to describe the most general form of the distribution of the FPP presented in this contribution. Adding noise as well gives the process described in Sec. III.4. This process does not in general have a closed form PDF, and so standard methods for estimating the process parameters which rely on the PDF does not work. We do, however, have a closed form for the characteristic function. Estimating parameters using the characteristic function has been discussed in Refs. Feuerverger and Mureika, 1977; Feuerverger and McDunnough, 1981; Yu, 2004. The main problem is finding a reasonable set of variables for the characteristic function. This has been considered by Refs. Tran, 1998; Zhang and He, 2016. All consider dividing the characteristic function into real and imaginary components explicitly. The approach taken here is more compact, and should be equivalent.

We have a set of independent and identically distributed observations . We assume we know the distribution family, but not the parameters; the PDF is given by with the corresponding characteristic function , where is a vector of parameters. The goal is to estimate these parameters from the observations. Define the empirical characteristic function

 CN(u)=1NN∑n=1exp(iuYn), (46)

and the error function

 εN(u,θ)=N1/2(CN(u)−C(u,θ)). (47)

Given a discrete set of sampling points for the characteristic function , the error vector will be given by (here and in the following, is suppressed for simplicity of notation)

 εN=[εN(u1),εN(u2),…,εN(uM)]T. (48)

It is known that is asymptotically normal with zero mean and an covariance matrix Feuerverger and Mureika (1977); Feuerverger and McDunnough (1981),

 Ωkl=⟨εN(uk)εN(ul)⟩=C(uk+ul)−C(uk)C(ul). (49)

An estimator of can be found by minimizing Tran (1998)

 εNˆΩ−1εTN, (50)

where is an estimate of . We estimate the covariance matrix by defining an -matrix with

 Ekl=exp(iukYl)−CN(uk). (51)

Then we have

 ˆΩ=11−NEET. (52)

Thus, in order to estimate the parameters , we choose the set of (as described below) and calculate . With an initial guess , can be constructed and Eq. (50) can be iteratively minimized to find the estimator .

Choosing the bin size for a histogram can have large effects on the resulting PDF. In the same way, choosing the sampling points has a large effect on the estimation procedure. Various approaches are discussed in Refs. Tran, 1998; Zhang and He, 2016. Not many points are needed; both agree on around 10 points as sufficient. Choosing is harder and requires the derivative of with respect to (and preferably, the PDF of , which we do not have). As a low-complexity, high-cost brute force method, we note that is insensitive to the initial values for a good choice of . Thus, we do the estimation for a large range of and many different initial values. The results are chosen where many initial values lead to the same result.

### iv.1 Examples of parameter estimation

In this section, we attempt to estimate the parameters from realizations of the process introduced in Sec. III.4. In the absence of noise, the method works well as long as the degree of pulse overlap is not excessive (). We therefore present particularly challenging examples from the full 3-parameter model. Three sets of parameters have been chosen. To describe an FPP with experimentally relevant intermittency level, exponentially distributed amplitudes and low noise level, the parameters were used. Both pulse overlap and high noise level lead to a distribution more closely resembling a normal distribution. Therefore, in order to reveal the sensitivity to intermittency in a process with high noise level, the parameters were used. Lastly, as strong intermittency leads to a more symmetric distribution, revealing the presence of moderate asymmetry in a process with high intermittency was tested with the parameters .

In Fig. 6, parameter estimation has been performed for parameters (top), (middle) and (bottom), using the L-BFGS-B - algorithm wrapped by the scipy.optimize.minimize - package Byrd et al. (1995); Zhu et al. (1997); Jones et al. (). The synthetic time series have data points, and we use and . For each , parameters have been estimated for all combinations of the initial values , and , giving 18 estimated values for each and each of , and . These estimated values are shown by the blue dots and orange stars in Fig. 6. Large scatter for a specific signifies high sensitivity to initial values in the estimation procedure. For , there was high dependence on initial values in all cases, so we only show results for . In all cases, the black line gives the true value and the green dashed line gives the estimated value. The estimated value is found by taking the mean of all orange stars, while ignoring blue dots. The criteria for deciding which values are used in the estimate are as follows. Top, we assume is known or suspected. Many points cluster around the line, and these all correspond to the same value of . Thus, these points are chosen for the parameter estimation, using for the orange stars. Middle, the estimates diverge for and in the case , so these values are not used, and marked with blue dots. Bottom, we have two fixed points for different initial values, one with high and , and low , and one with low and , and high . In this case, we choose to explain as little of the variation in the signal as possible with the additional noise level described by . Low corresponds to high and , so was used as a condition to mark estimated values with orange stars.

In Table 1, the estimated values are presented. This is the true parameter value subtracted from the mean value of the accepted estimated parameters with uncertainty equal to one standard deviation. For the two lowest example intermittency parameters, all estimated parameters are very close to the true value. For the example parameters , the procedure underestimates and overestimates . This is most likely due to the fact that as increases, the FPP approaches a normal distribution, making it difficult to separate the FPP from the normally distributed noise in the PDF. Still, all estimated parameters are within two standard deviations of the true parameter values.

## V Discussion and conclusions

In this contribution, we have investigated the filtered Poisson process given three different distributions of the pulse amplitudes; exponentially distributed amplitudes, Gamma distributed amplitudes and asymmetrically Laplace distributed amplitudes. For all of these, the mean, variance, and skewness and flatness moments of the resulting process were presented, as well as the parabolic relationship between the skewness and flatness moments. In addition, it was discussed how normally distributed noise affects the moments and the parabolic relation. In all cases, the characteristic function of the filtered Poisson process has a closed form solution, while the probability density function only has closed form solutions for exponentially distributed amplitudes, Gamma distributed amplitudes with shape parameter and symmetrically Laplace distributed amplitudes.

It was furthermore shown that exponentially distributed amplitudes and Gamma distributed amplitudes with shape parameter lead to PDFs for the FPP which are practically indistinguishable. In previous work Garcia and Theodorsen (2017b), it was shown that the amplitude distribution does not influence the power spectrum or the auto-correlation function of the FPP. Thus, a realization of an FPP with Gamma distributed amplitudes with shape parameter cannot be easily distinguished from a realization of a FFP with exponentially distributed amplitudes and a slightly larger intermittency parameter . The assumption of Gamma distributed amplitudes therefore leads to more complicated derivations and expressions, but practically equivalent predictions, for the PDF, auto-correlation and power spectrum, and can safely be simplified to the standard assumption of exponentially distributed amplitudes.

Laplace distributed amplitudes lead to expressions which are qualitatively different from exponentially or Gamma distributed amplitudes, since the two latter do not admit negative function values. Only for the case of symmetrically Laplace distributed amplitudes does the process have a closed form PDF. This makes parameter estimation methods requiring the PDF not applicable in general. However, the characteristic function has a closed form expression for any value of the asymmetry parameter. It has been demonstrated that a method based on the empirical characteristic function can be used to reliably estimate the model parameters in realizations of the process. This method can also handle additional noise to the process. The only problem in applying this method is deciding on the sampling points for computing the empirical characteristic function. In this contribution we have disregarded complex, iterative procedures for a simple brute-force method relying on the fact that the estimation procedure should be insensitive to initial parameter guesses for a good choice of sampling points. This method was capable of finding the correct noise ratio and asymmetry parameters in situations complicated by low intermittency (leading to more symmetric and Gaussian-like distributions) in the signal.

The FPP is a reference model for intermittent fluctuations in physical systems, where large amplitude bursts of similar shape dominate the fluctuations. Due to the Poisson process driving the model, it only considers turbulence which is statistically stationary in time. In magnetized plasmas, the model has been successfully applied to measurements of SOL fluctuations, where it is used to systematize and unify measurements. Finding the correct assumptions for the amplitude distributions and having good methods for estimating the model parameters is vital in being able to compare and contrast data for varying plasma parameters and machine configurations. In the future, comparisons between different machines will also be carried out using the framework presented here.

## Acknowledgements

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

## Appendix A Derivation of the characteristic function

In general, it has been assumed that the pulse amplitudes are independently and identically distributed, that the pulse arrival times are independently uniformly distributed, and that the pulses have a fixed shape. The characteristic function of the FPP (and its cumulants) has in this case been discussed in Refs. Rice, 1944, 1945; Parzen, 1999; Pécseli, 2000; Lowen and Teich, 2005.

In this section, we will derive the characteristic function of the FPP in a form as general as possible, and investigate exactly which assumptions are necessary in order to obtain a closed form expression. In its most general form, the FPP is given by

 ΦK(t)=K(T)∑k=1Akφ(t−skτk,λk), (53)

where the amplitudes , arrival times , duration times and asymmetries all are random variables. We begin by assuming that all random variables are independent and identically distributed across pulses, that is for all :

 pAk,λk,τk,sk(Ak,λk,τk,sk)=pAl,λl,τl,sl(Al,λl,τl,sl), (54)

and

 pAk,λk,τk,sk,Al,λl,τl,sl(Ak,λk,τk,sk,Al,λl,τl,sl)=pAk,λk,τk,sk(Ak,λk,τk,sk)pAl,λl,τl,sl(Al,λl,τl,sl). (55)

One could easily imagine an alternative FPP where for example the amplitude of the pulse depends on the waiting time from the last pulse. In this case, would depend on and , and a very different treatment from the present one would be required. In the following we drop the index on the random variables, since each pulse is statistically identical. The characteristic function of is the product of all characteristic functions of . Fixing for the moment, we have

 CΦ(u;K,t)=K∏k=1Cϕ(u;t)=Cϕ(u;t)K, (56)

where the variables after the semicolon are parameters in the characteristic function. By definition,

 Cϕ(u;t)=⟨exp(iuϕ)⟩=∞∫−∞dA∞∫−∞dλ∞∫−∞dτ∞∫−∞dspA,λ,τ,s(A,λ,τ,s)exp(iAuφ(t−sτ,λ)). (57)

The PDF of for fixed is

 PΦ(Φ|K)=12π∞∫−∞duexp(iuΦ)Cϕ(u;t)K. (58)

Using that is Poisson distributed,

 PK(K|T)=1K!(Tτw)Kexp(−Tτw), (59)

we have

 PΦ(Φ|T)=∞∑K=0PΦ(Φ|K)PK(K|T)=12π∞∫−∞duexp(iuΦ)exp{Tτw[Cϕ(u;t)−1]}. (60)

The expression inside the last exponential function can be identified as the logarithm of . Since the joint PDF of , , and integrates to by definition, we have that this expression is

 lnCΦ(u;T,t)=Tτw∞∫−∞dA∞∫−∞dλ∞∫−∞dτ∞∫−∞dspA,λ,τ,s(A,λ,τ,s)[exp(iAuφ(t−sτ,λ))−1], (61)

Since is Poisson distributed, is uniformly distributed on . Assuming that is independent of all the other random variables, we have

 lnCΦ(u;T,t)=1τw∞∫−∞dA∞∫−∞dλ∞∫−∞dτpA,λ,τ(A,λ,τ)T∫0ds[exp(iAuφ(t−sτ,λ))−1]. (62)

Exchanging the variable to , assuming stationarity and ignoring end effects for (that is, setting the integration limits back to ), we have

 lnCΦ(u)=1τw∞∫−∞d% A∞∫−∞dλ∞∫−∞dττpA,λ,τ(A,λ,τ)∞∫−∞dθ[exp(iAuφ(θ,λ))−1]. (63)

This is the most general form of the characteristic function of , where we only assume that pulses are independent of each other and that the arrivals follow a Poisson process and are independent of the other properties of the pulses.

Assuming that , and are independent, and that takes on a specific value, we arrive at

 lnCΦ(u)=γ∞∫−∞dAPA(A)∞∫−∞dθ[exp(iuAφ(θ))−1]. (64)

Changing the order of integration in Eq. (64) and using the definition of the characteristic function, we have Garcia and Theodorsen (2017b)

 lnCΦ(u)=γ∞∫−∞dθ⎡⎢⎣−1+∞∫−∞dAPA(A)exp(iuAφ(θ))⎤⎥⎦=γ∞∫−∞dθ[CA(uφ(θ))−1], (65)

where is the characteristic function for the amplitude distribution .

If we instead expand the exponential function in Eq. (63) into a sum, we have

 lnCΦ(u;T,t)=1τw∞∑n=1(iu)nn!∞∫−∞dA∞∫−∞dλ∞∫−∞dττpA,λ,τ(A,λ,τ)An∞∫−∞dθφ(θ,λ)n. (66)

The last integral gives , and

 lnC