Shock models governed by an inverse gamma mixed Poisson process

We study three classes of shock models governed by an inverse gamma mixed Poisson process (IGMP), namely a mixed Poisson process with an inverse gamma mixing distribution. In particular, we analyze (1) the extreme shock model, (2) the 𝛿 -shock model, and the (3) cumulative shock model. For the latter, we assume a constant and an exponentially distributed random threshold and consider different choices for the distribution of the amount of damage caused by a single shock. For all the treated cases, we obtain the survival function, together with the expected value and the variance of the failure time. Some properties of the inverse gamma mixed Poisson process are also disclosed.


Introduction
Nowadays, shock models are widely adopted for the description of lifetimes of systems working in random environments in various fields, including survival analysis, reliability engineering, insurance and finance (see Cha and Finkelstein [13]).In the literature, the existing shock models can be classified into five classes: (1) extreme shock models, where a system fails due to a single shock whose magnitude exceeds a fixed threshold (see, for instance, Cha and Finkelstein [9,11]); (2) cumulative shock models, in which a system fails when the cumulative damage caused by shocks exceeds a given threshold (see Gut [23] and Wang [38]); (3) run shock models, in which the failure of the system happens when the magnitudes of a given number of consecutive shocks are greater than a certain threshold (Mallor and Omey [28]; Ozkut [31] and Eryilmaz [17]); (4) -shock models, where we observe the failure of the device when the time between two consecutive shocks is less than a fixed threshold value (see Goyal, Hazra and Finkelstein [20]; Eryilmaz [16]; Jiang [25]; Li, Chan and Yuan [27]); (5) mixed shock models, in which different kinds of shock models are combined (see Wang and Zhang [39]).Other more general shock models have been studied over the last years: see, for instance, Belzunce et al. [4] for the case of multivariate shock models, Cha and Badía [8] for stochastically dependent dynamic shock models and Mallor and Santos [29] for a general shock model which includes a correlation structure.
A customary assumption in shock models is that the stochastic process governing the arrival of shocks is a homogeneous or a non-homogeneous Poisson process (see, for instance, Grandell [21]).The latter counting processes are both characterized by independent increments, this resulting in a high mathematical tractability of the model (see, for instance, Pellerey, Shaked, and Zinn [32] for some results concerning the epoch times and the inter-epoch intervals of a non-homogeneous Poisson process).However, this assumption seems unrealistic in many applications and, for this reason, some authors proposed different stochastic processes to model the occurrence of shocks.In particular, in the work by Eryilmaz [16], the -shock model is treated under the assumptions of shocks arriving according to a Pòlya process, i.e., a mixed Poisson process with mixing distribution given by a gamma probability density function.Furthermore, a generalized Pòlya process has been considered in Konno [26] and Cha [7], whereas in Di Crescenzo and Pellerey [15] the authors consider a geometric counting process governing the arrivals of shocks.
In the present manuscript we propose a mixed Poisson process with inverse gamma mixing distribution to describe the arrivals of shocks.As mentioned by Willmot [41], mixed Poisson processes are often used for modeling claim counts but, with the exception of few cases, they are difficult to be evaluated, this limiting their use.On the contrary, the counting process proposed in the present manuscript, named the inverse gamma mixed Poisson (IGMP) process, is characterized by a high mathematical tractability.In particular, the IGMP process shares with the Pòlya process the three characteristics stated in Beichelt [3]: both processes have (1) explicitly available finite dimensional distributions, (2) dependent increments, (3) free parameters allowing the model adaptation to a wide variety of data sets.Moreover, since the model properties arise from the choice of an inverse gamma mixing distribution, the IGMP process lends itself to being used in situations in which the underlying counting rate is better characterized by an inverse gamma distribution.
The inverse gamma distribution is one of the most useful distributions for describing different physical scenarios.For instance, in modeling shadow fading, the inverse gamma distribution is considered as an alternative to the gamma distribution since it admits a relatively simple mathematical formulation and exhibits heavy-tailed characteristics similar to the gamma density.The role of the inverse gamma distribution in describing the random variation of the mean signal power in composite fading channels is described in Yoo et al. [44], whereas an empirical validation of the inverse gamma distribution for modeling shadowing is provided in Yoo et al. [43].The use of the inverse gamma as a survival distribution is discussed in Glen [19], in which its upside-down bathtub (UBT) shaped hazard function is pinpointed.See also Ramirez-Espinosa and Lopez-Martinez [35] for a discussion concerning the benefits brought by the inverse gamma distribution to model shadowing, or Mead [30] for applications of such distribution in corrosion problems in new machines.The discretized form of the inverse gamma distribution has been recently proposed by Pundir et al. [34] as a suitable lifetime model.
Section 2 provides the main properties of the IGMP process.Similarly to other mixed Poisson processes, the inverse gamma mixed Poisson process is overdispersed.However, by analyzing the limit behavior of the coefficient of variation, the IGMP process, as time increases, tends to become qualitatively more regular than a Poisson process.
The distribution of the arrival time of the kth shock and that of the interarrival time between consecutive shock arrivals are also obtained in Section 2. In particular, we prove that the interarrival times follow the so-called K-distribution, which is characterized by heavy tails.This means that there is a high probability of observing long interarrival times and that the density does not decay rapidly.This property is shared by the Pòlya process, for which the interarrival times follow the Pareto distribution.However, the mean interarrival time of the Pòlya process can be infinite when the tail parameter is less than or equal to 1. On the contrary, the K-distribution typically has a finite mean and this property can be advantageous in applications where it is necessary to estimate or predict the average time between occurrences.
Throughout the manuscript, we consider three classes of shock models under the assumption of an IGMP process governing the arrival of shocks.In particular, in Section 3 we obtain the survival function for the extreme shock model, the -shock model and the cumulative shock model.In the latter case, we assume an exponentially distributed random threshold and consider different choices of the distribution of the random variable representing the amount of damage caused by the kth shock.We show that the IGMP process owns the weak positive dependence property, this meaning that the absence of shocks in a fixed time interval decreases the probability of a shock in the next short interval of time.Such property, which is shared by other processes such as the conditional Poisson process or the Geometric Counting Process (see Cha and Finkelstein [11]), eventually results in decreasing the failure rate.Indeed, both in the extreme and cumulative shock model (for an exponentially distributed random threshold) we prove that the system failure time is characterized by a decreasing failure rate.Hence, the probability of a failure in a fixed time interval in the future decreases over time, thus suggesting that the system improves as time goes on.This property provides a crucial distinction of the IGMP shock model from other well-known shock models for which the failure rate usually exhibits an increasing behavior.For instance, for an extreme shock model with constant system's failure probability, under the assumption of shocks governed by a generalized Pòlya process with constant baseline function, the failure rate function exhibits an increasing behavior (see Cha and Finkelstein [12]).
Finally, in order to highlight the practical use of the process in applications, in Appendix A.1 we propose some examples of real data fitting.In particular, we conduct a chi-square goodness of fit test to determine if the observed frequency distribution fits the IGMP distribution.The first dataset concerns the occurrences of earthquakes in Italy from 1,961 to 2,018; the second one describes the number of hourly vessel arrivals at the Hong Kong Port during May 2014.
It is worth pointing out that, as a byproduct of our investigations, in Sections 2 and 3 we obtain some results concerning the sum of series of modified Bessel functions of the second kind which have not been given before, to the best of our knowledge.

The inverse gamma mixed Poisson process
Let N(t) be a mixed Poisson process (MPP) with marginal distribution: where N ( ) (t) is a Poisson process with intensity > 0 and U (•) is an inverse gamma distribution defined on R + 0 with density: We shall call N(t) the inverse gamma mixed Poisson (IGMP) process and, for t > 0, the expression of its probability distribution is given by (see Willmot [41]) where K v (•) is the modified Bessel function of the second kind, defined as: Some plots of the probability distribution k (t) are shown in Figure 1 for different values of t.
Remark 2.1.We recall that, according to well-known properties of mixed Poisson processes (see for instance, Rolski et al. [36]), the probabilities k (t) for k = 1, 2, . . ., can be expressed in terms of the zero-state probability 0 (t).Indeed, the function 0 (t) is infinitely often differentiable such that: where https://doi.org/10.1017/S0269964823000232Published online by Cambridge University Press Recalling Eq. ( 1), it immediately follows that, for each t ≥ 0, As a byproduct of this result, the following corollary shows a simple relation involving the gamma function and a series of modified Bessel functions, which does not appear to have been given before.
Corollary 2.1.For > 0 and z > 0 we have The following Proposition provides some results concerning the IGMP process, which immediately follow from the properties of the mixed Poisson process (see, for instance, Rolski et al. [36] or Grandell [21]).To the best of our knowledge, the IGMP process has never been studied in detail in the literature.Proposition 2.1.Let N(t) be an inverse gamma mixed Poisson process.Then, for t > 0, the following properties hold: (i) N(t) admits stationary increments, i.e., for , > 0, n = 1, 2, . . .and each sequence {k r ; r = 1, . . ., n} of nonnegative integers, it is: where k = n r=1 k r and 0 .
In particular, the expected value and the variance of N(t) are given by: (iv) For , > 0, the probability generating function of N(t) is given by: Remark 2.2.From Eq. ( 7) we note that, for all t > 0, Var[N (t)] > E[N (t)], thus implying that the IGMP process is overdispersed.However, for > 2, the coefficient of variation of N(t) at time t > 0 is given by: Hence, the coefficient of variation is decreasing in t, and it is less than one for large t.This suggests that the IGMP process, as time increases, tends to become qualitatively more regular than a Poisson process.
Remark 2.3.From Eq. ( 8) we can deduce that adjacent increments of N(t) are positively correlated for > 2. Hence, for the IGMP process, a large number of events in a given time interval tends to lead to a large number of events in the next interval.
In the following Theorem we prove that the IGMP process owns the weak positive dependence property.We recall that a counting process {M (t), t ≥ 0} satisfies the weak positive dependence property if where s, t > 0 and I(E) represents the indicator function for the event E (see, for instance, Cha and Finkelstein [11]).
Theorem 2.1.The IGMP process has the weak positive dependence property.
https://doi.org/10.1017/S0269964823000232Published online by Cambridge University Press Proof.Due to Eq. ( 10), the IGMP satisfies the weak positive dependence property if the inequality: holds for any s, t ∈ R + .Due to Proposition 6 and recalling Eq. ( 1), Eq. ( 11) reads: In order to prove this inequality, we define the functions: where recalling that (see Gaunt [18]): it results 0 < f (x) < 1 and g(x) > 1, ∀x > 0. Due to the above assumptions, the proof of Eq. ( 12) can be obtained by showing that, for x, h > 0, one has with g(0) = 1, being (see, for instance, Abramowitz and Stegun [1]): For small h we have that , which means that the function ln g(x) is concave, or, equivalently, ln f (x) is convex.This latter condition can be proved noting that due to a Turan-type inequality (see Baricz [2]), for any positive it results Hence, being the proof easily follows.
The following Proposition provides a further characterization of the IGMP process in terms of a suitable stochastic ordering.We recall that, given two discrete random variables Y and W, with discrete densities h and g, respectively, Y is said to be smaller than W in the likelihood ratio order, Y ≤ lr W, if increases in t over the union of the supports of Y and W, (15) (see, for instance, Belzunce, Martinez Riquelme and Mulero [5]).We recall that stochastic orders allow to compare probability distributions or measures and represent a powerful tool used in many areas such https://doi.org/10.1017/S0269964823000232Published online by Cambridge University Press as reliability, survival analysis, risks, finance, and economics.In particular, the likelihood ratio order is stronger than other orderings, such as the stochastic order, this justifying the goal of exploring its validity.
Proposition 2.2.Let N(t) be a IGMP process.Then, for , > 0, N(t) increases in the likelihood ratio order when t increases.
Proof.We aim to show that for 0 < s < t, it is N (s) ≤ lr N (t), i.e., due to Eq. ( 15), the ratio P(N (t)=n) is increasing with respect to n.
Being (see Eq. ( 3)) we have to show that, for n ∈ N, it is i.e. the function is increasing with respect to z.Noting that, due to Eq. (9.6) of Abramowitz and Stegun [1], the proof immediately follows from Eq. (2.2) of Yang and Zheng [42].

Distribution of arrival and interarrival times
For the IGMP process, let us denote by T k , k ∈ N + , the arrival time of the kth event and by , the interarrival time between consecutive events.The explicit expression of the distribution of T k is provided in the following theorem.
Some plots of the probability density ( 16) are shown in Figure 2 for different values of the parameters and and for different choices of k.
In the next proposition we obtain the joint distribution of the arrival epochs T k : Proposition 2.3.The joint density of the vector T = (T 1 , . . ., T m ), for (t 1 , . . ., t m ) ∈ R m such that 0 < t 1 < . . .< t m , is given by: Proof.Recalling Theorem 5.2 of Chapter VIII in Rolski et al. [36], we have that Hence, due to Eq. (3), we easily obtain Eq. (17).
Some properties of the interarrival times are disclosed in the following propositions.Let us denote by X m = (X 1 , X 2 , . . ., X m ) the vector of the interarrival times of the IGMP process.Proposition 2.4.For x 1 , x 2 , . . ., x m > 0 we have Proof.Due to Theorem 5.3 of Chapter VIII in Rolski et al. [36], it is: where (m) 0 is defined in Eq. ( 5).Hence the proof follows by making use of straightforward calculations.
As an immediate consequence of Proposition 2.4, we obtain the marginal density of the interarrival times X k , k ∈ N + .Proposition 2.5.For all x > 0, k ∈ N + , the marginal density of the interarrival times X k is given by: The density (18) is the so-called K-distribution (see Ward, Watts and Tough [40]) where is the shape factor, is the scale factor, and K (•) is the modified Bessel function of second kind defined in Eq. ( 4).
Starting from Eq. ( 18), we provide the expression of the expected value and the variance of the interarrival times Aiming to detect the dependence structure of the interarrival times X k , in the following Proposition we evaluate the autocorrelation function, defined as: Var(X j )Var(X j+i ) , i = 1, 2, . . . .Proposition 2.6.For > 0, the autocorrelation function of the interarrival times is given by: Proof.Recalling Eq. (4.6) of Cox and Lewis [14], it is possible to evaluate the autocorrelation function as: with where G N (t) (•) is the probability generating function of the IGMP process, provided in Eq. ( 9).Due to Eq. ( 9), and making use of straightforward calculations, we have that https://doi.org/10.1017/S0269964823000232Published online by Cambridge University Press which identifies with the probability generating function of a mixed Poisson process with mixing distribution given by an inverse gamma density with parameters ( + 1, ).Due to Eq. ( 21) and recalling Eq. ( 3), we thus have Hence, the proof immediately follows from Eq. ( 20) and recalling Corollary 2.2.
Note that the autocorrelation function (19) does not depend on the lag; this is in agreement with the stationarity of the interarrival times series.

Shock models
Let us consider a system subjected to external shocks that arrive according to the IGMP process, defined in Section 2. Due to the weak positive dependence (see Theorem 2.1), in this model the absence of shocks in a fixed time interval decreases the probability of a shock in the next short interval of time, this being common for some kind of shock processes.To the best of our knowledge, shock models governed by IGMP processes are new in the scientific landscape and they represent an alternative to the well-known Pòlya shock models.
Denoting by S the failure time of the system, we are interested in the study of the survival function and the rate function of S in three different cases: (i) extreme shock models (systems fail due to one single large shock), (ii) cumulative shock models (systems fail because of some cumulative effect) and (iii) -shock models (systems fail when the time between two consecutive shocks is less than a fixed value ).

Extreme shock model
Let us consider the case in which the shocks occur independently and each one causes the failure of the system with the same probability p ∈ (0, 1).Denoting by M the geometric random variable counting the number of shocks which cause the failure of the system, and by assuming M independent of N(t), for t ≥ 0 we have: where k (t) is defined in Eq. ( 1).Note that the assumption of independence between M and the process N(t) could be seen as restrictive in the description of real-world phenomena.Indeed, the occurrence of previous shocks could weaken the system, making it susceptible to further shocks (see, for instance, Cha and Finkelstein [10], for a shock model based on a history-dependent approach).Such consideration could be regarded as the starting point for a future investigation.
The explicit forms of the survival probability and the expected value and variance of the failure time S are provided in the following theorem.In addition, the corresponding failure rate is also obtained.Theorem 3.1.Let N(t) be the IGMP process.In the extreme shock model, the survival function of the system failure time S and the corresponding failure rate are given by: Moreover, the expected value and the variance of S are Proof.Due to Eq. ( 22) and recalling Eq. ( 1), for t ≥ 0 and 0 < p < 1, we have where G N (t) (z) is the probability generating function of the IGMP process (Eq.( 9)).
The corresponding failure rate function is given by: Being K +1 (z) = 2 z K (z) + K −1 (z) (see, for instance, Abramowitz and Stegun [1]), Eq. ( 24) follows after straightforward calculations.Finally, the expected value and the variance of S can be obtained from Eq. ( 23).
Remark 3.1.The failure rate ( 24) is decreasing for all t > 0 and tends to zero as t goes to infinity.This result immediately follows from the Hankel's expansion for the modified Bessel function of the second kind (see Abramowitz and Stegun [1]), which, for sufficiently large z, reads: As a byproduct of Theorem 3.1.,recalling Eqs. ( 22) and (3), in Corollary 3.1 we obtain the explicit value of a series of modified Bessel functions of the second kind.Corollary 3.1.For > 0 and x > 0 we have The following corollary provides the explicit expressions of the survival function and the rate function for a particular choice of .

Corollary 3.2. Under the assumptions of Theorem 3.1, in the special case
Proof.It immediately follows from Theorem 3.1 recalling that

Cumulative shock model
Let us denote by W k , k ∈ N + , the random variable representing the amount of damage caused by the kth shock, whose arrival is governed by the IGMP process N(t).We assume that {W k , k ∈ N + } is a sequence of independent and identically distributed random variables with the distribution of the r.v.W, and that they are independent of N(t).
The total accumulated damage of the system from all shocks in [0, t) is given by: where we make the customary assumption 0 k=1 W k = 0.The failure of the system occurs when the cumulative effect of shocks reaches a certain deterministic or random level R. Hence, the survival function of the system is given by: In the sequel, we shall analyze two choices for the random level R, by considering the case of (i) constant and (ii) exponentially distributed threshold.
In the case of constant threshold, the following Proposition provides the expression of the survival function of the system failure time S, together with the expected value and variance.Proposition 3.1.Let us assume that N(t) is an IGMP process and suppose that the threshold is a deterministic boundary, R = d.Then, the survival function of the system failure time S is given by: where F (k) W (•) is the k-fold convolution of the distribution function F W (•).Moreover, the expected value and the variance of S are given by: https://doi.org/10.1017/S0269964823000232Published online by Cambridge University Press In the particular case of exponentially distributed random variables W k , k ∈ N with parameter 1/ , > 0, Proposition (3.5) can be specified in the following result.
Proposition 3.2.Let us assume that N(t) is an IGMP process and suppose that the random variables W k , k ∈ N, are exponentially distributed with parameter 1/ , > 0. For a deterministic threshold R = d, the survival function of the system failure time S is given by: whereas for the expected value and the variance of S we have , , , > 0.
In the sequel, we shall analyze the case of an exponentially distributed random threshold R with mean > 0. We denote by L W (s) := E e −sW , the Laplace transform of the probability density function of W.
Theorem 3.2.Let N(t) be the IGMP process.In the cumulative shock model, for an exponentially distributed random threshold R with mean > 0 and by assuming L W 1 < 1, the survival function of the failure time of the system can be expressed as: Furthermore, the failure rate function is given by: Proof.For t ≥ 0, due to Eq. ( 26), we have: where G N (t) (z) is the probability generating function of the IGMP process N(t) and it is provided in Eq. ( 9).Hence, Eqs. ( 29) and (30) follow by means of straightforward calculations.
https://doi.org/10.1017/S0269964823000232Published online by Cambridge University Press Figure 3 provides some plots of the expected value of the failure time when the random variables representing the amount of damage are Erlang (LHS) and uniformly distributed (RHS).In both cases the mean failure time is increasing as function of , which is the expected value of the threshold R.Moreover, in the Erlang case, E(S) is decreasing as the parameter n increases.Similarly, in the uniform case, the expected value of the failure time decreases as the width of the interval [0, b] increases.This confirms the intuition that large values of the damage result in earlier system failure.

-Shock model
In the extreme (Section 3.1) and in the cumulative shock model (Section 3.2), the failure of the system depends on the magnitude of shocks.However, in some cases, the system fails if the time between one shock and the next is shorter than a fixed threshold .This model is known in the literature as the -shock model (see, for instance, Eryilmaz [16]).
In the present Section we aim to analyze the -shock model under the assumption that the occurrence of shocks is governed by the IGMP process, this implying that the interarrival times between consecutive shocks X 1 , X 2 , . . .are identically distributed (see Section 2.1) but not independent.Recalling that S is the random variable representing the failure time of the system, the following Proposition provides the explicit expression of the survival function of S in the case of a -shock model.Proposition 3.3.In the -shock model, under the assumption of shocks governed by an IGMP process, for t > 0 the survival function of the system is given by: where x denotes the integer part of x.
Proof.By conditioning on the number of shocks occurred in (0, t), we have https://doi.org/10.1017/S0269964823000232Published online by Cambridge University Press Hence, being (see, for instance, Theorem 5.4 of Chapter VIII in Rolski et al. [36]) where x + = max(x, 0), the proof immediately follows from Eq. ( 3).
Let us consider the random variable M counting the number of shocks which determine the failure of the system.The probability law of M is disclosed in the following Proposition.Proposition 3.4.For m = 1, 2, . .., the probability law of the random variable M is given by: Moreover, for the expected value we have Proof. Being conditioning on and recalling that the random variables X 1 , X 2 , . . .are conditionally independent, we have Hence, due to Eq. ( 2), Eq. ( 37) immediately follows.Finally E(M) can be obtained starting from Eq. ( 37) making use of straightforward calculations.
Remark 3.4.The distribution of M can be expressed in terms of the probability law of the IGMP process.Indeed, recalling Eq. ( 3), for m = 1, 2, . .., we have In the following Proposition we provide the expression of the expected value of the failure time S. The proof proceeds along the lines of Theorem 2 of Eryilmaz [16] and so it is omitted.Proposition 3.5.In the -shock model, if N(t) is a IGMP process, for , > 0 the expected value of the failure time S is given by: where the random variable M * is distributed as M (Eq.37) with parameters ( + 1, , ).  Figure 4 shows some plots of the expected failure time in the -shock model scenario.We note that E(S) is always decreasing as increases and, for fixed ( ), it takes larger (smaller) values for smaller values of ( ).

Conclusions
In the present manuscript, we analyzed three different classes of shock models governed by an inverse gamma mixed Poisson process.Generalizations of the Poisson process play a key role in modeling since, in many applied contexts, the assumption of exponentially distributed and independent interarrival times may be unrealistic.On the contrary, under the assumption of IGMP process, the interarrival times between shocks are dependent and represent a stationary series.The interarrival times follow the K-distribution which is characterized by heavy tails, this meaning that long interarrival times are observed with high probability.This property is shared by the Pòlya process, for which the interarrival times follow the Pareto distribution.However, the mean interarrival time of the Pòlya process can be infinite whereas the K-distribution has typically a finite mean.This property can be advantageous in applications where it is necessary to estimate or predict the average time between occurrences.
Shock models governed by IGMP processes are characterized by explicitly available finite dimensional distributions and free parameters allowing the model adaptation to a wide variety of data sets.For this reason they represent an alternative to shock models based on Pólya processes.Among the different classes of shock models we have analyzed: (i) the extreme shock model, (ii) the cumulative shock model and (iii) the -shock model and we have obtained the explicit expression of the survival function and the mean failure time.In particular, both in the extreme and cumulative shock models we prove that the system failure time is characterized by a decreasing failure rate.This property provides a crucial distinction of the IGMP shock model from other shock models characterized by an increasing failure rate.For instance, an extreme shock model with constant systems failure probability, under the assumption of shocks arriving according to a generalized Pólya process with constant baseline function, exhibits an increasing failure rate function.
In the proposed shock models certain independence assumptions could be weakened in order to describe more realistic scenarios.Indeed, in a future research, in order to show potential applications in engineering and actuarial sciences, we aim to study a cumulative shock model in which the damages W k depend on the inter-arrivals X k (see, for instance, Sumita and Shanthikumar [37]).Moreover, other possible future researches involve (i) an extreme shock model when there are m possible sources of shocks acting on the system (see, for instance, Bozbulut and Eryilmaz [6]) and (ii) a generalization of the cumulative shock model with a threshold R following a three parameter generalized exponential distribution (see Gupta and Kundu [22]).

• A case study on seismology data
We provide the observed frequencies and compare them with the expected frequencies evaluated for (i) the IGMP process and (ii) the Poisson process.In all the cases, the parameters involved in the distribution are estimated by means of the Method of Moments and by using one year as the time unit.The critical value 2 0.05,n is calculated at a significance level equal to 0.05 in the case of n degrees of freedom.The observations are classified into the mutually exclusive classes shown in Table A1 for the case (i) and in Table A2 for the case (ii).
-Case (i): the estimated values of parameters of the IGMP process are = 5.03 and = 3.61.Since the 2 -statistic takes value 3.07 < 2 0.05,1 = 3.84, we cannot reject the null hypothesis of the distribution of the IGMP process.
-Case (ii): the estimated value of the parameter of the Poisson distribution (A.1) is = 0.90.Hence, the 2 -statistic takes value 5.02 > 2 0.05,1 = 3.84 and we have to reject the null hypothesis of the distribution of the Poisson process.

• A case study on vessel arrival data
As in the previous case, we carry out a chi-square test concerning (i) the IGMP process and (ii) the Poisson process.Since the dataset (cf.Table A1 of Huang [24]) counts how many times a fixed number of ships arrived within an hour, we set the time unit equal to one hour.The parameters in the distributions are estimated by the Method of Moments, and the critical value 2 0.05,n is calculated at a significance level equal to 0.05 in the case of n degrees of freedom.The classes in which the data are divided, together with the observed and expected frequencies, are shown in Table A3.
-Case (i): the estimated values of parameters of the IGMP process are = 14.8 and = 41.38.Since the 2 -statistic takes value 5.5 < 2 0.05,6 = 12.59, we cannot reject the null hypothesis of the distribution of the IGMP process.
-Case (ii): the estimated value of the parameter of the Poisson distribution (A.1) is = 2.99.The 2statistic takes value 16.61 > 2 0.05,7 = 14.07, so that we have to reject the null hypothesis of the distribution of the Poisson process.

Table A1 .
Observed and expected frequencies of earthquakes in Italy for the IGMP process.

Table A2 .
Observed and expected frequencies of earthquakes in Italy for the Poisson process.