Exponentially Weighted Moving Averages of Counting Processes When the Time between Events Is Weibull Distributed

There are control charts for Poisson counts, zero-inflated Poisson counts, and over dispersed Poisson counts (negative binomial counts) but nothing on counting processes when the time between events (TBEs) is Weibull distributed. In our experience the in-control distribution for time between events is often Weibull distributed in applications. Counting processes are not Poisson distributed or negative binomial distributed when the time between events is Weibull distributed. This is a gap in the literature meaning that there is no help for practitioners when this is the case. This book chapter is designed to close this gap and provide an approach that could be helpful to those applying control charts in such cases.


Introduction
Statistical process control and monitoring (SPCM) methods originally arose in the context of industrial/manufacturing applications, developed during and after World War II. Since then, it has become a popular way of monitoring all processes. Today, large volumes of data are often available from a variety of sources, in a variety of environments that need to be monitored. This means one needs to make sense of these data and then be able to make efficient monitoring decisions based on them. While constructing and applying monitoring tools, a fundamental assumption, necessary to justify the end results, involves the assumption about the distribution that the data have been generated from. This is the heart of outbreak detection of events particularly for estimating the in-control false discovery rates.
Selecting an appropriate probability distribution for the data is one of the most important and challenging aspects of data analysis. The estimates of outbreak false discovery rates often hinge on this crucial selection. We focus on the distribution for the time between events (TBEs) because there are often many of these events as compared to counts. Therefore, it is easier to fit an appropriate distribution using these many TBE values. The most commonly assumed distribution in the application of TBE is the Weibull distribution which is asymmetric and sometimes severely skewed. However, depending on the context, other distributions may also be used, such as the exponential distribution (which is a special case of a Weibull distribution) or the gamma distribution. If the TBE distribution is exponentially distributed, then the related counts are Poisson distributed. This book chapter focusses on counting processes when the distribution of TBE values is known to be Weibull distributed.
The challenge of making and meeting the distributional assumption is faced by all practitioners and data analysts. In many monitoring settings, event data are collected in a nearly continuous stream, and it is often more meaningful to monitor the individual TBE data [1][2][3][4] when outbreaks are of large magnitude. This individual event data are aggregated over fixed time intervals (e.g., daily) to form counts. In this chapter, these counts are monitored to detect outbreaks resulting from small changes in the incident of events. The focus is on the steady-state situation because this is the most common situations in event monitoring. Note that we cannot stop the process and investigate the out-of-control situation because often in nonmanufacturing settings it is not under our control. Events may include warranty claims of a product, health presentation at emergency departments, sales of an online products, etc. Here the term "quality of the process" is used in a general sense, which is context dependent. In the case of sales, an outbreak would represent an increased sales opportunity, provided the inventory stock can support this outbreak and the products are not sold out before the next order arrives. However, for warrantee claims, this would represent an undesirable outbreak of increased claims which may require a failure mode effect analysis [5]. Monitoring of in-control nonhomogeneous counting processes has traditionally been carried out using either Poisson or the negative binomial distributions for the counts. Many statistical tools used in SPCM for counts data are documented in Sparks et al. [1][2][3], Sparks et al. [6], Weiß [7,8], Weiß and Testik [9,10], Yontay et al. [11], and Albarracin et al. [12]. These control charts are perhaps the most well-known count monitoring methodologies. The control chart graphic is a time series plot of a signal-to-noise ratio designed for the user to make decisions about the outbreak of events.
In any case, before defining an "in-control" process, we need the information about the probability distribution of the events being monitored. When this information is available, it is possible to calculate the probabilities (or the chance) from the event in-control distribution, which could be defined by the TBEs or the counting of events within a fixed time interval. Deciding on whether event distribution constitutes an event outbreak is based on whether the event distribution is extreme compared to what is usual, i.e., counts are higher than expected. This is usually gauged by some upper threshold for the signal-to-noise ratios.
In a vast majority of SPCM applications, it is common to assume that the underlying probability distribution is of a (given) known form. In this chapter we assume that the TBEs are distributed as a Weibull distribution. However, we explore approaches to monitoring the counts of these events over fixed periods of varying length to find the period width that leads to earlier detection of outbreaks in terms of the average time to signal (ATS). The distribution of these counts when the TBEs are Weibull distributed is neither Poisson nor negative binomial distributed. Therefore, this chapter offers a different approach compared to others in the literature. In addition, the appropriate period of aggregation for the counts is explored in terms of how it influences early detection of outbreak events.
In practice event data are often collected in a continuous stream defined by the TBEs. Besides, wherever outbreaks are of large magnitude, it is more meaningful to monitor these individual TBE values [4]. In such situations, one also needs to deal with the issue of autocorrelation which may be thought of as the effect of time, as data values in close proximity with respect to time and space are likely dependent. This violates one of the basic assumptions in process monitoring, and a common way to deal with this issue is by fitting a time series model and monitoring the standardized residuals. The assumption of a distribution is also an important part of this analysis. In this chapter, our focus is on the Weibull distributional assumption for the TBE values. However, rather than monitoring the TBEs, we monitor the counts over a fixed time interval, because this improves the early detection of smaller outbreaks (Sparks et al., [4]).

Monitoring homogeneous counts
Considering a fixed distribution for TBEs during the monitoring period, we assume that the time of the day and date stamp of all events are available. For a series comprised of n events, let the day numbers for events be denoted as with the event times within days which are defined as Note that these times are measured in fractions of a day, e.g., 0 ≤ τ i < 1. Daily counts are given by counting the number d i values that are the same and then adding the days with zero counts for those days that are missing from the d i values. Then TBEs are given by where w i represents the time between i þ 1 ð Þth and ith events. We flag an outbreak wherever w i for any i ¼ 1, … , n À 1 are consistently lower than an expected value. These TBEs are assumed to be Weibull distributed with fixed scale and shape parameters. Using R code, we define the counting process, for say daily counts, as Denote these counts as c i for the ith day. We define the exponentially weighted moving average (EWMA) statistic for these daily homogeneous counts as follows: where e 0 ¼ E c i ð Þ is assumed constant when in-control for days i ¼ 1, 2, … :. The smoothing parameter,α 0 < α < 1 ð Þ , determines the level of memory of past observations included in e i . Smaller values of α retain the more memory of past counts than do the larger values. Therefore, smaller values of α are more efficient at detecting smaller changes in mean counts, while larger values of α are more efficient at detecting larger changes. EWMA statistic, e i , has a minimum value of zero, and we do not allow it to fall below zero because we are only interested in outbreaks, i.e., counts that are higher than expected. This EWMA statistic controls the worst-case scenarios, whereas the traditional EWMA statistic can wonder well below zero before going out of control causing the EWMA statistic to take a long time to signal the outbreak. An outbreak is flagged when e i exceeds a predetermined threshold, h c ARL 0 , scale, shape, α ð Þ . For a given pair of shape and scale parameters, h c is determined so that a desired in-control average run length (ARL 0 ) is achieved. Appendix A provides models for establishing the thresholds for various values for 0:2 ≤ scale ≤ 0:46, 0:6 ≤ shape ≤ 1:4, α ¼ 0:1 and ARL 0 ¼ 100, 200, 300 or 400: Given a h c and the process is out of control, we want the ARL to be as low as possible. Monitoring daily counts, the ARL is defined as the number of days before a signaled outbreak. Later we will use the ATS to assess the relative performance of control chart plans, because the ARL may vary according to the aggregation period we select for the counts. Note that a false outbreak is flagged by this EWMA statistic being significantly larger than expected given it is in-control.

Monitoring in-control nonhomogeneous counts
For most of the real-world cases, TBE distribution may vary due to different circumstances while not experiencing an outbreak. For instance, if we aggregate the TBEs in daily intervals, occurring intervals for the events may vary based on the day of the week, and even working and non-working days may affect the distribution of TBEs. As a result, we often face nonhomogeneous count processes. Hence, we define an adaptive exponentially weighted move average (AEWMA) statistic for nonhomogeneous daily counts as where ae i is the AEWMA statistic at time t, and ae 0 ¼ The other notations are as defined earlier. To control the false discovery rate, we set h c ARL 0 ¼ c, scale, shape, α ð Þ so that a desiredARL 0 is achieved. An outbreak is flagged approximately whenever ae i > 1.
The α 0 < α < 1 ð Þin Eq. (1) determines the level of memory of past observations in this average ae i . Smaller α values retain more memory of past counts; therefore, small α values are efficient at retaining enough memory of past counts to have the power to flag smaller outbreaks. However, larger values of α are needed when there is a larger size outbreak, because a shorter range of memory is adequate to build enough power to detect the outbreak. In addition, the shorter memory for the EWMA average with larger values of α is less inclined to be influenced by too many past in-control counts.

Simulation results
To assess the validity and applicability of our proposed adaptive method, we employ simulations studies. We restrict our attention to plans that have 10 < E c i ð Þ < 50 which limits the volume of simulations that are required to make a realistic judgment. This is also the range where the EWMA counts become competitive for small changes in the scale parameter (see [4] and Sparks et al., 2020). To do so, we consider four counting processes with daily aggregations. The shape parameters for the TBE distributions when the processes are in-control are equal to 0.85, 0.95, 1.15, and 1.25. The in-control scale parameter for each process is equal to 0.02, 0.025, 0.03, and 0.035, respectively. We intentionally imposed outbreaks in the simulated data to assess if the proposed method is capable of detecting them. In the simulation study, we assume that the outbreaks result in a decrease in the scale parameter. Various surveillance plans are devised to monitor each process. The ATS 0 for the plans associated with the first process (with shape and scale parameters of 0.85 and 0.02, respectively) is set to be equal to 400. The ATS 0 for the plans associated with other processes are equal to 300, 200, and 100, respectively. 100,000 events are employed to estimate the count mean for each process. In addition, 500 events are also considered for the burn-in period of the simulations. The performance of the devised plans is presented in Tables 1-4. The lowest ATS values are colored in the tables below in black bold text to make it easier to see which plans are more efficient in certain situations. Table 1 shows the performance results for the plans employed to monitor counting process where the in-control data are Weibull with scale and shape parameters of 0.035 and 1.25, respectively. ATS 0 (in days) for these plans are set approximately equal to 100. As shown in Table 1 for the early detection of the small outbreaks, smaller values for the smoothing parameter are preferred. Table 2 shows the performance results of the plans when the shape and the scale parameters of the Weibull distribution are equal to 1.15 and 0.03, respectively. The ATS 0 for the plans employed to monitor this counting process is set to be approximately 200. In the results indicated in Table 1, plans with smaller smoothing parameter work better in early detection of small outbreaks than those with larger smoothing parameter. For the larger outbreaks, the detection power of the plans increases as the smoothing parameters increases. Table 3 shows the performance results of the plans when the shape and the scale parameters of the Weibull distribution are equal to 0.95 and 0.025, respectively. The ATS 0 for the plans employed to monitor this counting process is set to be EWMA counts (α), shape = 0.95 approximately 300. Conclusions similar to those regarding Tables 1 and 2 can be drawn from Table 3. It is clearly observed that the larger the magnitude of the outbreak, the larger the smoothing parameter should be. On the other hand, for the detection of small outbreaks, plans with larger values of ATS 0 need even smaller values for smoothing parameters than do the plans with smaller ATS 0 . As shown in Table 1, smoothing parameters equal to or larger than 0.08 work better for outbreak detection. As the ATS 0 increases, analogous to results presented in Tables 2  and 3, even smaller values for α are needed to devise a plan of larger detection power. Last but not the least, similar results can be driven from Table 4, which presents the performance of the monitoring plans applied to a counting process with underlying Weibull distribution for TBEs. The shape and the scale parameters of the aforementioned distribution are equal to 0.85 and 0.028, respectively.

Real-world example
In this section, we apply our proposed method to a real-world example. The counting process to monitor is the number of presentations at Gold Coast University Hospital emergency department for a broad definition of influenza. The events are presentations at Gold Coast University Hospital with flu symptoms. Data is gathered for four consecutive years starting from January 2015. We first check if the TBEs are Weibull distributed. We fitted a Weibull regression model to data, and all the parameters of the conditional distribution of the response variable are modeled using explanatory variables as the hour of the day harmonics, seasonal harmonics, and day of the week. To do so, we employ "gamlss" [13] R package for statistical modeling. This package includes functions for fitting the generalized additive models for location, scale, and shape introduced by Rigby and Stasinopoulos [14]. The R procedure for the model fitting is as follows: gamlss(formula = TBE $ wd + wd=="Monday") + (wd=="Thursday")) * (cos(2 * pi * hr/24) + sin(2 * pi * hr/24)) + cos (2 * pi * nday/365.25) + sin (2 * pi * nday/365.25)sigma.formula = $(cos2 * pi * hr/24)+ sin (2 * pi * hr/24))+cos(2 * pi * nday/365.25)+sin(2 * pi * day/365.25) family=WEI()data=data) where wd, nday, and hr. are week day, number of the day in a year, and the time of day (0-24), respectively. Figure 1 summarizes the analysis of the residuals for the fitted model. As shown in Figure 1, since residuals do not represent any particular pattern in data, and the model describes the response variables quite well, then we conclude that TBEs are Weibull distributed. Details for the analyzing the model adequacy is presented in Appendix B.
As mentioned earlier, the threshold, h c , is determined according to parameters of the underlying Weibull distribution for the TBEs. For any timestamp during the monitoring period, the fitted model is used to predict the parameters of the Weibull distribution. Then these parameters along with the smoothing parameter, α, and the in-control ARL are substituted in models in Appendix A to establish the h c for Weibull-distributed counts. The parameter for the day of the count is taken as the average of the estimated parameters for the day of events. Considering α ¼ 0:1 and ARL 0 for the plan to be 400, the estimated parameters are used to establish the AEWMA plan for monitoring the daily counts of flu presentations at the Gold Coast University Hospital. Figure 2 shows the time series of the monitoring statistic over the study period, and Figure 3 is the time series of the daily flu counts. The proposed plan indicates six potential flu outbreaks during 4 years as the monitoring statistic falls beyond 1. Note that influenza counts in 2015 were bigger than usual, but not unusual in 2016 (apart from early in that year). However, in 2017 this was flagged as a very unusual influenza outbreak and this seems to persist into early 2018.

Concluding remarks
In this chapter, we proposed an adaptive EWMA surveillance plan to monitor a counting process of which the time between its events is Weibull distributed. The proposed method can be applied to both homogeneous and nonhomogeneous processes. To implement the proposed surveillance plan, the scale and shape parameters for the underlying distribution of the TBEs are estimated using a distributional regression approach. Then the threshold for the counts is established using the estimated parameters and the desired ARL. The proposed plan is applied to both simulated and real data. Simulation results indicate that the proposed method is applicable for detecting outbreaks of any magnitude and also signals them in a reasonable time after their incidence. In addition, simulations revealed that for the detection of the large outbreak, plans with larger smoothing parameter are superior. However, for the early detection of small outbreaks, we need to employ smaller smoothing weights. Applying the proposed surveillance method to real data, we conclude that the proposed method is capable of detecting outbreaks in nonhomogeneous counting processes.