Modification of ARL for detecting changes on the double EWMA chart in time series data with the autoregressive model

Abstract

This research aims to derive the average run length (ARL) evaluation of the double exponentially weighted moving average (double EWMA) control chart for observation data that follows exponential white noise in a time series model with an autoregressive model. Since most real-world data is automatically correlated, autoregressive models are available. Comparisons were made between the ARLs obtained using the explicit formula and the numerical integral equation (NIE) approach. The results showed that the explicit formula's use of the ARL outperformed the NIE approach in terms of computation time. After that, the efficacy of the exponentially weighted moving average (EWMA) and double EWMA charts is then compared using the suggested explicit ARL formula. The ARL of the double EWMA chart was found to perform better than the ARL of the EWMA chart in all situations. It also uses natural gas and diesel prices on stock exchanges around the world as cases studies. The results show that the double EWMA chart has better detection sensitivity than the EWMA chart, and the results are consistent with the experimental results. As a result, the sensitivity of the double EWMA chart in detecting changes makes it a good alternative for monitoring processes with real-world data.

KEYWORDS:

• Average run length
• autoregressive model
• double EWMA control chart
• EWMA control chart
• explicit formulas
• numerical integral equation

1. Introduction

Time-series data analysis has been the subject of extensive research. Currently, it is utilised for predicting a variety of circumstances, including stock trends epidemic diseases, export sales and cost prediction (Cheng et al., 2020; Satrio et al., 2021). There are common uses for these models as an autoregressive (AR), a moving average (MA), and an autoregressive moving average (ARMA). Statistical process control (SPC) has been widely used in the industrial sector to monitor, control, and improve operations. Many researchers have studied SPC using real-world data and case studies. For instance, Tegegne et al. (2022) suggested a new kind of multivariate control chart for the manufacturing process, which was utilised to better decide on quality control during the production process. They were conducted to validate the control chart by using a case study from the cement industry. Yeganeh et al. (2023) employed SPC techniques, particularly control charts, to track the probability of connections between network surveillance system nodes using control charts based on machine learning.

The control chart is widely used in the fields of finance, business, engineering, healthcare, and others to track processes and detects sudden changes in the mean or variance of those processes, especially in the field of industries (Fugen, 2019; Zhang et al., 2022). Using control charts, it is possible to distinguish between common and unique sources of process variation and handle each case separately. The control limits on control charts are determined by common factors. Control charts provide precise instructions on how and when to modify a process. Nowadays, control charts are an essential tool in SPC, and their industrial significance lies in their ability to monitor and control process variability, improve process efficiency, reduce costs, ensure product quality, and support decision-making. There are many types of research that were presented that utilised the control chart in the field of industries. Shrestha (2021) used a control chart to control the sausage manufacturing process. Tegegne et al. (2022) used a control chart to decide on quality control during the production process in their case study of the Ethiopian cement industry, and so on.

The concept of the monitoring process was initially introduced in Shewhart's control chart (1931). The quality of the Shewhart control chart is only effectively achieved by identifying large shifts in the process parameters. There are many viable alternatives to well-known and efficient technology for detecting and monitoring slight changes in processes. The cumulative sum (CUSUM) control chart was created by Page (1954), and Montgomery (2009) provided more information on the chart's specifics. The CUSUM control chart was employed in the current investigation because of its greater capacity to identify minor parameter adjustments and to benefit from the advantages of employing an upper-sided control chart. The exponentially weighted moving average (EWMA) control chart was suggested by Roberts (1959), and it is generally recognised that it may be used to monitor and identify slight changes. After that, several researchers proposed a variety of different control charts utilising the EWMA control chart development to improve the ability to detect slight changes. The modified exponentially weighted moving average (modified EWMA) control chart was devised by Khan et al. (2017) after being initially proposed by Patel and Divecha (2011). The extended exponentially weighted moving average (extended EWMA) control chart was suggested by Naveed et al. (2018). The extended EWMA statistics, which can detect minor changes faster than the EWMA control chart shown in Karoon et al. (2022) and Karoon et al. (2023), were created by adding the EWMA statistics with the exponential smoothing parameters. This approach was shown to be more successful at recognising tiny changes. Additionally, Shamma and Shamma (1992) introduced the double exponentially weighted moving average (double EWMA) control chart before Mahmoud and Woodall (2010) offered it. The double EWMA control chart performs better than the EWMA control chart when the monitor process is utilised with small shift sizes rapidly. It is now generally acknowledged that a reliable control chart is essential for quickly detecting subtle process changes.

The average run length (ARL) is frequently utilised to assess control chart performance. The ARL is separated into two states; the ARL is ARL₀, which is described by the value as it happens from the in-control state and should be large, and ARL₁, which is explained by the value as it occurs from the out-of-control stage and should be as tiny as feasible. A variety of methods can be employed to estimate the ARL of the control chart, including Monte Carlo and Markov Chain. Moreover, there is Numerical Integral Equations (NIE) method (see Karoon et al., 2021). Also, Areepong mentioned the NIE method for SPC in 2009. For quality control, it was used to identify the distribution change point. To evaluate the principal components of EWMA control charts for various distributions, NIE was proposed to derive a numerical algorithm and infer precise analytical formulas. In addition, ARLs were calculated using explicit formulas. There is a large body of literature that cites the explicit ARL formulas for control charts with time series models. Due to reality, the control chart's performance may suffer if the observations are serially associated. For this reason, the ARL has been suggested by several researchers when there is a serial correlation in the observations. Areepong and Sukparungsee (2016) proposed an explicit formula for the ARL on the EWMA control chart with the ARIMA model. After that, Sunthornwat et al. (2017) derived the analytical ARL on the EWMA control chart for a long-memory autoregressive fractionally integrated moving average (ARFIMA) model. On the EWMA control chart, Peerajit et al. (2017) proposed the explicit formulas for the derivation of the ARL using the data are the ARFIMA model on the CUSUM chart. Later, The capacity of the ARL was put out by Supharakonsakun et al. (2020), and it was obtained from an explicit formula for the MA(1) model of the modified EWMA control chart. Supharakonsakun (2021) proposed an explicit formula of the ARL with observations of the MA(p) model. After that, Phanyaem (2022) presented an explicit formula and NIE method of the ARL for the SARX(p,r)_L model based on the CUSUM control chart. Petcharat (2022) proposed the ARL of the CUSUM control chart when the data are SARX(p)_L model. Moreover, Karoon et al. (2022) proposed an explicit formula of ARL on the extended EWMA control chart when the observations are the AR(p) model. Recently, Peerajit and Areepong (2023) proposed the ARL of an autoregressive fractionally integrated process with exponential white noise running on the modified EWMA control chart.

As was already indicated, these results show that the ARL may be used to compare efficiency with control charts and that explicit formulas require considerably less calculation time to assess the ARL than the other approaches. However, the explicit formulas for the double EWMA chart have not yet been released. As a result, the focus of this research was on developing an ARL for using an AR model to monitor change in a time series. This method is based on the double EWMA chart with underlying exponential white noise in different scenarios. The ARL that was obtained by an explicit formula was provided. And then, that is expanded for comparison with the typical EWMA control chart. Besides, the ARL for them is also extended to compare processes with the diesel and natural gas price datasets.

2. Process identification on the EWMA and double EWMA control charts

2.1. EWMA control chart

Roberts (1959) created the first suggestion for the EWMA chart. It is a tool that is typically used to monitor and identify slight changes in the process. The equation below can be used to express the statistics of the EWMA control chart: (1) $E_t=\lambda Y_t+(1-\lambda )E_{t-1}$(1) where the EWMA chart parameter $Y_t$ is a sequence of autoregressive AR(p) model and a sequence data at $t=1,2,3,\dots$ with exponential white noise, $\lambda$ is an exponential smoothing parameter $(0<\lambda \le 1)$, $Y_t(t=0)$ is the initial value of the EWMA statistics that is equal to $\mu$ and the variance of $Y_t$ is $\lambda {\sigma }^2/(2-\lambda )$. Therefore, the upper control limit (UCL) and lower control limit (LCL) of the EWMA chart can be established from the mean $\mu$, the standard deviation $\sigma$, and a suitable control width limit $Q_E$ as follows below: (2) $\mu \pm Q_E\sigma \sqrt{\frac{\lambda }{2-\lambda }}.$(2) The stopping time of the EWMA control chart is expressed as: (3) ${\tau }_h=inf\{t\ge 0:E_t>h\},\text{for}h\ge \psi ,$(3) where $h$ is the upper control limit, a constant parameter.

2.2. Double EWMA control chart

The double EWMA control chart was proposed by Shamma and Shamma (1992) and Mahmoud and Woodall (2010) made it available. It was adjusted by two exponential smoothing and then extended from the EWMA chart. The equation below may be used to describe the statistics of the double EWMA chart. (4) $\begin{array}{rl}{D}_t& ={\lambda }_2{E}_t+(1-{\lambda }_2)D_{t-1},\\ E_t& ={\lambda }_1{Y}_t+(1-{\lambda }_1)E_{t-1}.\end{array}$(4) where the double EWMA control chart parameter $Y_t$ is a sequence of autoregressive AR(p) model and a sequence data at $t=1,2,3,\dots$ with exponential white noise, ${\lambda }_1$ and ${\lambda }_2$ is an exponential smoothing parameter $0<{\lambda }_1,{\lambda }_2\le 1$, $Y_t(t=0)$ is the initial value of the double EWMA statistics that is equal to $\mu$ and the variance of $Y_t$ is $\frac{{\lambda }_1^2{\lambda }_2^2}{{({\lambda }_1-{\lambda }_2)}^2}{\sigma }^2\left[\frac{{(1-{\lambda }_2)}^2}{1-{(1-{\lambda }_2)}^2}+\frac{{(1-{\lambda }_1)}^2}{1-{(1-{\lambda }_1)}^2}-2\frac{(1-{\lambda }_1)(1-{\lambda }_2)}{1-(1-{\lambda }_1)(1-{\lambda }_2)}\right]$. Therefore, the upper control limit (UCL) and lower control limit (LCL) of the EWMA control chart can be established from the mean $\mu$, the standard deviation $\sigma$, and a suitable control width limit $Q_D$ as follows below: (5) $\mu \pm Q_D\sigma \sqrt{\frac{{\lambda }_1^2{\lambda }_2^2}{{({\lambda }_1-{\lambda }_2)}^2}\left[\frac{{(1-{\lambda }_2)}^2}{1-{(1-{\lambda }_2)}^2}+\frac{{(1-{\lambda }_1)}^2}{1-{(1-{\lambda }_1)}^2}-2\frac{(1-{\lambda }_1)(1-{\lambda }_2)}{1-(1-{\lambda }_1)(1-{\lambda }_2)}\right]}.$(5) The stopping time of the EWMA control chart is expressed as: (6) ${\tau }_h=inf\{t\ge 0:D_t>h\},\text{for}h\ge \psi ,$(6) where $h$ is the upper control limit, a constant parameter.

However, ${\lambda }_1$ values of the double EWMA statistics were replaced by 1, and the double EWMA statistics is the EWMA statistics.

2.3. Cumulative Sum control chart (CUSUM)

Page (1954) introduced the CUSUM control chart in quality control for identifying minor variations in process mean. The equation below can be used to express the statistics of the CUSUM chart. (7) $C_t=max(0,C_{t-1}+Y_t-a),t=1,2,\dots$(7) where a is non-zero constant, $C_0=u$ is the initial value of CUSUM with $u\in [0,h]$, $h$ is upper control limit (UCL), and the stopping time of the CUSUM chart is defined as (8) ${\tau }_b=inf\{t>0;C_t>h\},h>u.$(8)

2.4. The ARL on double EWMA control chart for time series with AR model

The time series is a collection of chronologically organised data. The time series is observed and studied in order to identify its patterns of change and development and forecast its future trend. Time-series data can be divided into two categories: stationary data and non-stationary data. The time-series data collection with stationary data does not show a trend or a seasonal influence. The data collection solely contains random error as a cause of variation, whereas non-stationary time-series data is a time-series data set that exhibits a trend or a seasonal effect. There are other sources of variation in the data set besides just random error. The Box–Jenkins method is currently the most complete and accurate algorithm for analysing and predicting time series data. Methods for measuring timed data are referred to as times series. Autoregression (AR), Moving Average (MA), Autoregressive Moving Average (ARMA), Autoregressive Integrated Moving Average (ARIMA), and Seasonal Autoregressive Integrated Moving-Average (SARIMA) are examples of common types.

One of the most commonly used models is the autoregressive model, which was also used in this study. Regression models, known as autoregressive (AR) models, have a dependent or response variable that is a linear function of its previous values. The autoregressive model has an important parameter r, p, which is the order of an autoregressive term and is often used for actual data. The AR(p) model is defined by the following equation: (9) $Y_t=\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p}+{\zeta }_t$(9) where $\omega$ is the constant of the model, ${\varphi }_i$, $i=1,2,\dots ,p$ are coefficients of autoregressive $(0\le {\varphi }_i\le 1)$. It is assumed that the error term (${\zeta }_t$) is an exponential white noise (${\zeta }_t\sim \text{Exp}(\beta )$). The probability density function of ${\zeta }_t$ is given as $f(y,\beta )=\frac{1}{\beta }{e}^{-\frac{y}{\beta }}$; $\beta >0$, and then initial values of the autoregressive AR(p) model; $Y_{t-1},Y_{t-2},\dots ,Y_{t-p}$ are equal to 1. In earlier studies, ARL formulae for monitoring change processes using data that were characterised as the AR model were developed utilising the equation of the AR model in Equation (9). For instance, Phanthuna et al. (2021), and Karoon et al. used the AR model in Equation (9) to fit with ARL for the EWMA type control chart. For the double EWMA control chart running an autoregressive model under exponential white noise in this study, let $L(\psi )$ represent the ARL conditional on the initial value. The initial value for monitoring the double EWMA statistic $D_0=\psi$ is established at $\psi \in [0,h]$. Let $P_{\theta }$ denote the induced expectation and probability measure associated with $\psi$ at change-point ($\theta$) under the density function $f(y,\beta )$. The ARL for AR(p) model on the double EWMA control chart defined with $L(\psi )$ can be described as follows: $L(\psi )=\text{ARL}=\{\begin{array}{lll}\text{AR}{\text{L}}_0=P_{\mathrm{\infty }}({\tau }_h),& \theta =\mathrm{\infty },& (\text{no}\text{change})\\ \text{AR}{\text{L}}_1=P_1({\tau }_h),& \theta =1,& (\text{change})\end{array},$ where ARL₀ refers to the in-control ARL with $\theta =\mathrm{\infty }$, which denotes no change in the statistical process. On the other hand, ARL₁ refers to the out-of-control ARL with $\theta =1$, which designates the first-time point in the statistical process at which a change occurs place from ${\beta }_0$ to ${\beta }_1$.

3. The methodology and the performance measurement

3.1. Explicit formulas of ARL on Double EWMA control chart for AR model

In this section, the explicit formula of the ARL is solved on the double EWMA chart with the AR model. Let’s start by substituting Equation (9) into Equation (4) as follows: $D_t={\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p}+{\zeta }_t)+{\lambda }_2(1-{\lambda }_1)E_{t-1}+(1-{\lambda }_2)D_{t-1}$ when the first time $t=1$ such that $D_0=\psi$ is determined, then the initial values $E_0=\gamma$ and $Y_{t-1},Y_{t-2},\dots ,Y_{t-p}$ equal to 1. The double EWMA statistics with AR(p) can be described as: $D_1={\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p}+{\zeta }_1)+{\lambda }_2(1-{\lambda }_1)\gamma +(1-{\lambda }_2)\psi$ In control process, the interval of $D_1$ between the lower and upper bound control limits are expressed to be $l$ and $h$ can be written as follows below. $l<{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})+{\lambda }_2(1-{\lambda }_1)\gamma +(1-{\lambda }_2)\psi +{\lambda }_1{\lambda }_2{\zeta }_1<h.$ On the variable ${\zeta }_1$, this interval can be rewritten as: $\begin{array}{rl}& \left[\begin{array}{c}{\displaystyle \frac{l-(1-{\lambda }_2)\psi }{{\lambda }_1{\lambda }_2}}-{\displaystyle \frac{{\lambda }_2(1-{\lambda }_1)\gamma }{{\lambda }_1{\lambda }_2}}-{\displaystyle \frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})}{{\lambda }_1{\lambda }_2}}\end{array}\right]<{\zeta }_1\\ & <\left[\begin{array}{c}{\displaystyle \frac{h-(1-{\lambda }_2)\psi }{{\lambda }_1{\lambda }_2}}-{\displaystyle \frac{{\lambda }_2(1-{\lambda }_1)\gamma }{{\lambda }_1{\lambda }_2}}-{\displaystyle \frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})}{{\lambda }_1{\lambda }_2}}\end{array}\right].\end{array}$ The Fredholm integral equation is used to express the integral equation of the ARL on the double EWMA control chart for the AR(p) model with an initial value $D_0=\psi$. The equation represented is as follows: $\begin{array}{rlr}L(\psi )& =\underset{\left[\frac{l-(1-{\lambda }_2)\psi }{{\lambda }_1{\lambda }_2}-\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{{\lambda }_1{\lambda }_2}-\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{{\lambda }_1{\lambda }_2}\right]}{\overset{\left[\frac{h-(1-{\lambda }_2)\psi }{{\lambda }_1{\lambda }_2}-\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{{\lambda }_1{\lambda }_2}-\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{{\lambda }_1{\lambda }_2}\right]}{\int }}& \times L\left(\begin{array}{c}{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})\\ +{\lambda }_2(1-{\lambda }_1)\gamma +(1-{\lambda }_2)\psi +{\lambda }_1{\lambda }_2{\zeta }_1\end{array}\right)\cdot f(\eta )d\eta .\end{array}$ If $\eta ={\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})+{\lambda }_2(1-{\lambda }_1)\gamma +(1-{\lambda }_2)\psi +{\lambda }_1{\lambda }_2{\zeta }_1$ is obtained by substituting the integration variable, $L(\psi )$ is restructured as follows: (10) $\begin{array}{rl}L(\psi )& =1+\frac{1}{{\lambda }_1{\lambda }_2}\underset{l}{\overset{h}{\int }}L(\eta )f\left(\frac{\begin{array}{c}\eta -(1-{\lambda }_2)\psi -{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}\\ +{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})-{\lambda }_2(1-{\lambda }_1)\gamma \end{array}}{{\lambda }_1{\lambda }_2}\right)d\eta .\end{array}$(10) As the function of ${\zeta }_1$ has an exponential distribution, $L(\psi )$ can be expressed as follows: (11) $L(\psi )=1+\frac{e^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}}{\beta {\lambda }_1{\lambda }_2}\underset{l}{\overset{h}{\int }}L(\eta )\cdot e^{\frac{-\eta }{\beta {\lambda }_1{\lambda }_2}}d\eta .$(11) The ARL solution is verified by Banach's fixed-point theorem. This is considered an ARL solution in terms of its existence and uniqueness (see Sofonea et al., 2005), from Equation (11), suppose that $\ell (\psi )=\frac{e^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}}{\beta {\lambda }_1{\lambda }_2}\text{and}\rho =\underset{l}{\overset{h}{\int }}L(\eta )\cdot e^{\frac{-\eta }{\beta {\lambda }_1{\lambda }_2}}d\eta .$ Consequently, the ARL solution that is provided in Equation (11) can be reformulated by establishing additional variables as follows: (12) $L(\psi )=1+\ell (\psi )\cdot \rho$(12) After that, the integral equation $\rho$, which can be written as: (13) $\begin{array}{rl}\rho & =\underset{l}{\overset{h}{\int }}{e}^{\frac{-\eta }{\beta {\lambda }_1{\lambda }_2}}\cdot \left(1+\frac{e^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}}{\beta {\lambda }_1{\lambda }_2}\cdot \rho \right)d\eta \\ \rho & =-\beta {\lambda }_1{\lambda }_2\left[e^{\frac{-h}{\beta {\lambda }_1{\lambda }_2}}-e^{\frac{-l}{\beta {\lambda }_1{\lambda }_2}}\right]-\frac{e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})+{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot \rho }{{\lambda }_2}\left[e^{\frac{-h}{\beta {\lambda }_1}}-e^{\frac{-l}{\beta {\lambda }_1}}\right]\\ \rho & =\frac{-\beta {\lambda }_1{\lambda }_2\left[e^{\frac{-h}{\beta {\lambda }_1{\lambda }_2}}-e^{\frac{-l}{\beta {\lambda }_1{\lambda }_2}}\right]}{1+\frac{e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})+{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}}{{\lambda }_2}\cdot \left[e^{\frac{-h}{\beta {\lambda }_1}}-e^{\frac{-l}{\beta {\lambda }_1}}\right]}.\end{array}$(13) Finally, Equation (12) is substituted into the solution of $\rho$ as shown in Equation (13), which is obtained as: (14) $\begin{array}{rl}L(\psi )& =1-\frac{e^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}}{\beta {\lambda }_1{\lambda }_2}\\ & \cdot \frac{\beta {\lambda }_1{\lambda }_2\left[e^{\frac{-h}{\beta {\lambda }_1{\lambda }_2}}-e^{\frac{-l}{\beta {\lambda }_1{\lambda }_2}}\right]}{1+\frac{e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})+{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}}{{\lambda }_2}\cdot \left[e^{\frac{-h}{\beta {\lambda }_1}}-e^{\frac{-l}{\beta {\lambda }_1}}\right]}\\ L(\psi )& =1-\frac{{\lambda }_2{e}^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot \left[e^{\frac{-h}{\beta {\lambda }_1{\lambda }_2}}-e^{\frac{-l}{\beta {\lambda }_1{\lambda }_2}}\right]}{{\lambda }_2{e}^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})+{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}+\left[e^{\frac{-h}{\beta {\lambda }_1}}-e^{\frac{-l}{\beta {\lambda }_1}}\right]}.\end{array}$(14) Since the AR(p) model is used with the double EWMA chart, the explicit ARL formula is given in Equation (14). Moreover, the in-control process is replace $\beta$ with ${\beta }_0$, whereas the out-of-control process is replace $\beta$ with ${\beta }_1$, as well as ${\beta }_1=(1+\delta ){\beta }_0$, and $\delta$ denotes the shift size in the monitoring process.

3.2. Numerical integral equation of ARL on double EWMA control chart for AR model

The numerical integral equation method, or NIE method, of the ARL is used to compare with the explicit formula in terms of the computational time spent. Let $N(\psi )$ denote the ARL of the double EWMA control chart on the AR(p) model with an exponential white noise, with the ARL solution being computed according to the midpoint quadrature rule. Specifically, after selecting a quadrature rule, the interval $[l,h]$ is divided into $l\le v_1\le v_2\le \dots \le v_m\le h$ using a set of constant weights $w_j=(h-l)/m$.

In accordance with the quadrature rule, the approximation for an integral is assessed as follows: (15) $\underset{l}{\overset{h}{\int }}L(\eta )f(\eta )d\eta \approx \sum _{j=1}^m{w}_{j}f(v_j)$(15) Let $N(v_i)$ represent a numerical approximation to the integral equation. It is the result of solving the following linear equations: (16) $\begin{array}{rl}N(v_i)& =1+\frac{1}{{\lambda }_1{\lambda }_2}\sum _{j=1}^{m}N(v_j)\cdot f\left(\frac{\begin{array}{c}{\eta }_j-(1-{\lambda }_2)\psi -{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}\\ +{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})-{\lambda }_2(1-{\lambda }_1)\gamma \end{array}}{{\lambda }_1{\lambda }_2}\right)\\ & i=1,2,\dots ,m.\end{array}$(16) The m linear equation system is shown as $N_{m\times 1}+R_{m\times m}{N}_{m\times 1}$. If the inverse $(N_m-R_{m\times m}{)}^{-1}$ exists, the unique solution is shown as: $N_{m\times 1}=(I_m-R_{m\times m}{)}^{-1}\cdot 1_{m\times 1},$ where $N_{m\times 1}=[N(v_1),N(v_2),\dots ,N(v_m){]}^{\mathrm{\prime }}$, $diag(1,1,\dots ,1)$ is the unit matrix order m, $1_{m\times 1}=[1,1,\dots ,1{]}^{\mathrm{\prime }}$ is a column vector of $N(v_i)$, and $R_{m\times m}$ is a matrix, the definition of (m,mth) elements of R matrix is expressed as $[R_{ij}]\approx \frac{1}{{\lambda }_1{\lambda }_2}{w}_{j}f\left(\frac{\begin{array}{c}{\eta }_j-(1-{\lambda }_2)\psi -{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})\\ -{\lambda }_2(1-{\lambda }_1)\gamma \end{array}}{{\lambda }_1{\lambda }_2}\right).$ Finally, $v_i$ is instead of $\psi$ into $N(v_i)$, the approximation of the NIE method of the ARL is represented as follows: (17) $\begin{array}{rl}N(\psi )& \approx 1+\frac{1}{{\lambda }_1{\lambda }_2}\sum _{j=1}^{m}N({\psi }_j)\cdot f\left(\frac{\begin{array}{c}{\eta }_j-(1-{\lambda }_2)\psi -{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}\\ +{\varphi }_2{Y}_{t-2}+\cdots +{\varphi }_p{Y}_{t-p})-{\lambda }_2(1-{\lambda }_1)\gamma \end{array}}{{\lambda }_1{\lambda }_2}\right)\end{array}$(17) where ${\eta }_j$ is a set of the division point within the interval ${\eta }_j=(j-0.5)w_j+l$ for $j=1,2,\dots ,m$ and $w_j$ is a weight of the composite midpoint formula; $w_j=(h-l)/m$.

3.3. The existence and uniqueness of the explicit formula for the ARL based on the double EWMA control chart

In this section, the existence and uniqueness of the ARL solution are confirmed using Banach's fixed-point theorem. Banach’s Fixed-point Theorem, sometimes referred to the contraction theorem, is concern with certain mappings of a complete metric space into itself. As we shall see, a fixed point is a point that is mapped to itself, and this statement establishes requirements adequate for the existence and uniqueness of a fixed point.

Theorem 3.1.

Banach’s Fixed-point Theorem:

Let $(Y,d)$ be a complete metric space and let $T:Y\to Y$ be a contraction mapping on Y. Then T has a unique fixed point $y\in Y$(such that $T(y)=y$) with contraction constant $r\in [0,1)$ such that $||T(L{(\psi )}_1)-T(L{(\psi )}_2)||\le r||L{(\psi )}_1-L{(\psi )}_2||$, for all $L(\psi {)}_1,L(\psi {)}_2\in Y$. Then there exists a unique $L(\psi )\in Y$ such that $T(L(\psi ))=L(\psi )$, i.e. a unique fixed-point in $Y$ (Sofonea et al., 2005).

Proof:

Let T defined in Equation (11) is a contraction mapping for $L(\psi {)}_1,L(\psi {)}_2\in u[l,h]$.

To show that $||T(L{(\psi )}_1)-T(L{(\psi )}_2)||\le r||L{(\psi )}_1-(L{(\psi )}_2||$ for all $L(\psi {)}_1,L(\psi {)}_2\in Y$ with $r\in [0,1)$.

Consider: $||T(L{(\psi )}_1)-T(L{(\psi )}_2)|{|}_{\mathrm{\infty }}$ $\begin{array}{rl}& =\underset{u\in [l,h]}{sup}|(L{(\psi )}_1-(L{(\psi )}_2|\\ & =\underset{u\in [l,h]}{sup}\left|\frac{e^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}}{\beta {\lambda }_1{\lambda }_2}\underset{l}{\overset{h}{\int }}(L_1(g)-L_2(g))\cdot e^{\frac{-g}{\beta {\lambda }_1{\lambda }_2}}dg\right|\\ & \le \underset{u\in [l,h]}{sup}|{||L{(\psi )}_1-L{(\psi )}_2||}_{\mathrm{\infty }}{e}^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}\\ & \times \left(e^{\frac{-l}{\beta {\lambda }_1{\lambda }_2}}-e^{\frac{-h}{\beta {\lambda }_1{\lambda }_2}}\right)|\\ & =||L{(\psi )}_1-L{(\psi )}_2|{|}_{\mathrm{\infty }}\underset{u\in [l,h]}{sup}\\ & \times \left|e^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}\right|\left|e^{\frac{-l}{\beta {\lambda }_1{\lambda }_2}}-e^{\frac{-h}{\beta {\lambda }_1{\lambda }_2}}\right|\\ & \le r||L{(\psi )}_1-L{(\psi )}_2|{|}_{\mathrm{\infty }}\end{array}$ where $r=\underset{u\in [l,h]}{sup}\left|e^{\frac{(1-{\lambda }_2)\psi }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_2(1-{\lambda }_1)\gamma }{\beta {\lambda }_1{\lambda }_2}}\cdot e^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})}{\beta }}\right|\left|e^{\frac{-l}{\beta {\lambda }_1{\lambda }_2}}-e^{\frac{-h}{\beta {\lambda }_1{\lambda }_2}}\right|;r\in [0,1).$

3.4. The performance measurement of control chart

The performance of control charts is commonly evaluated in terms of the average run length (ARL). It is used to detect changes in the mean of the monitoring process. As a result, the performances of the ARL, which are produced by explicit equations, and the NIE technique, which is used to monitor shifts with an autoregressive model, or AR model, on the double EWMA chart, were compared.

The absolute percentage relative change (%APRC) indicates the relative performance of two methods of the ARL (Phanthuna & Areepong, 2022), which is defined as follows: (18) $\mathrm{\%}APRC=\left|\frac{L(\psi )-N(\psi )}{L(\psi )}\right|\times 100\mathrm{\%}$(18) After that, the performance of the ARL is then evaluated with different parameter values based on the double EWMA chart. And then it is extended for comparison with the EWMA chart. The ARL measurement is appropriate to evaluate the effectiveness of control charts in the process parameters when shift sizes are different. The performance of a control chart throughout a range of shifts (${\delta }_{min}\le \delta \le {\delta }_{max}$) can be evaluated using overall performance measures, which are recommended by several researchers. Some of them have metrics that are used to assess their performance, such as the average extra quadratic loss (AEQL) and the performance comparison index (PCI) (see Alevizakos et al., 2021).

The AEQL is described mathematically as (19) $\text{AEQL}=\frac{1}{\mathrm{\Delta }}\sum _{\delta ={\delta }_{min}}^{{\delta }_{max}}({\delta }^2\times 1-\frac{{\lambda }_2{e}^{\frac{(1-{\lambda }_2)\psi }{(1+\delta )\beta \cdot {\lambda }_1{\lambda }_2}}\cdot \left[e^{\frac{-h}{(1+\delta )\beta \cdot {\lambda }_1{\lambda }_2}}-e^{\frac{-l}{(1+\delta )\beta \cdot {\lambda }_1{\lambda }_2}}\right]}{\begin{array}{c}{\lambda }_2{e}^{\frac{{\lambda }_1{\lambda }_2(\omega +{\varphi }_1{Y}_{t-1}+{\varphi }_2{Y}_{t-2}+\dots +{\varphi }_p{Y}_{t-p})+{\lambda }_2(1-{\lambda }_1)\gamma }{(1+\delta )\beta \cdot {\lambda }_1{\lambda }_2}}\\ +\left[e^{\frac{-h}{(1+\delta )\beta \cdot {\lambda }_1}}-e^{\frac{-l}{(1+\delta )\beta \cdot {\lambda }_1}}\right]\end{array}})$(19) where the specific of shift in the process is represented by $\delta$, and then $\mathrm{\Delta }$ is the sum of number of divisions from ${\delta }_{min}$ to ${\delta }_{max}$. In the study, $\mathrm{\Delta }=11$ is determined from ${\delta }_{min}=0$ to ${\delta }_{max}=1$. The most effective control chart is those with the lowest AEQL values.

The PCI measurement is the ratio of the AEQL of the control chart and the AEQL of the most effective control chart (shown as the AEQL_lowest). It is described mathematically as (20) $\text{PCI}=\frac{\text{AEQL}}{\text{AEQ}{\text{L}}_{\text{lowest}}}.$(20) The most effective control chart's PCI value is 1, whereas the values of the other control charts in competition are more than 1.

4. The results

A common criterion for assessing the efficacy of a control chart is the ARL for identifying alterations in the process. That is an evaluation of its practicality and sensitivity. In this study, the efficiency of the explicit formula and the NIE method for determining the ARL for monitoring shifts were compared. These were evaluated utilising the double EWMA chart and the autoregressive (AR) model. The AR(1), AR(2), and AR(3) autoregressive models were expressed. The number of division points m = 500 is determined for the approximation of the ARL by the NIE method. The solutions have been split into two categories: simulated and real-world data. For an initial ARL, In the in-control process scenario, ARL₀ is followed by exponential white noise with a mean of $\beta ={\beta }_0=1$. In the out-of-control process scenario, ARL₁ was computed with parameter values as exponential white noise with a mean of $\beta ={\beta }_1$, as ${\beta }_1={\beta }_0(1+\delta )$, and then $\delta$ equals 0.001, 0.002, 0.003, 0.005, 0.01, 0.03, 0.05, 0.1, 0.5, 1. As previously stated, “exponential white noise” is defined as the residual of an exponential distribution with uncorrelated data. As a result, $\delta =0$ indicates that the process is in-control, whereas $\delta >0$ indicates that the process is out-of-control. ARL₁ could have a low value on the efficacy control chart. Also, ARL₀ was determined at 370, and the Mathematica, spec Intel(R) Xeon(R) CPU X5680 @ 3.33 GHz 3.33 GHz (2 processors) RAM 32.0 GB was used to evaluate through study.

The process can be succinctly stated as follows:

1. Step 1: Input parameters for example, the coefficients of autoregressive; ${\varphi }_i$, the initial values of the autoregressive; $Y_{t-1},Y_{t-2},\dots ,Y_{t-p}$, and the control chart parameters; ${\lambda }_1,{\lambda }_2.$

2. Step 2: Impose the parameter of exponential white noise;${\beta }_0$ for the in-control process and assign the parameter of exponential white noise;${\beta }_1=(1+\delta ){\beta }_0$ for the out-of-control process.

3. Step 3: Determine the initial value of known parameter; $\omega$ and the initial value of the double EWMA statistic; $\psi$.

4. Step 4: Specify the lower control limit ($l$) equals 0 and fixed the initial value ARL₀ for computing compute the upper control limit ($h$).

5. Step 5: Compare the ARL of the explicit formula and the NIE method.

6. Step 6: Compute the solution of ARL₁ for shift sizes in the monitoring process where ${\beta }_1=(1+\delta ){\beta }_0$, and then the upper control limit ($h$) is the result of Step 4.

4.1. Simulated data

The aforementioned algorithm was used to determine ARL₀ and ARL₁, as well as to evaluate the efficiency of the proposed ARL by the explicit formula and the ARL by the NIE method in situations where the in-control and out-of-control processes were present. The results computed based on specified parameters such as ${\varphi }_i$ and $\omega$ for ARL₀ = 370 using the ARL on the double EWMA chart with the AR(p) model as AR(1), AR(2), and AR(3) models obtained from two techniques are reported in Table 1. The coefficient parameters for the AR(1), AR(2), and AR(3) models were set as ${\varphi }_1=\pm 0.1$, ${\varphi }_1={\varphi }_2=\pm 0.1$, ${\varphi }_1={\varphi }_2={\varphi }_3=\pm 0.1$, respectively. The ARL values of the explicit formula technique denote as $L(\psi )$, and it is computed by using Equation (14). And then, the NIE method denotes as $N(\psi )$, it is calculated by Equation (17). The %APRC of all situations is extremely low, nearly equal to 0. However, the computational times for the ARL values calculated using the NIE technique are about 3.7–4.3 s, while the explicit formula appears practically instantly for all situations. Therefore, using these explicit formulas as the next step is a good idea.

Table 1. Comparing the ARL values of the explicit formula against the NIE method for AR(p) models on the double EWMA control chart with the initial ARL₀ = 370 and $\omega =0$ based on different situations.

	Model	AR (1)		AR (2)		AR (3)
	$({\lambda }_{1,},{\lambda }_2)$	(0.05,0.10)		(0.10,0.10)		(0.20,0.10)
	${\varphi }_i$	${\varphi }_1=0.1$	${\varphi }_1=-0.1$	${\varphi }_1={\varphi }_2=0.1$	${\varphi }_1={\varphi }_2=-0.1$	${\varphi }_1={\varphi }_2={\varphi }_3=0.1$	${\varphi }_1={\varphi }_2={\varphi }_3=-0.1$
$\delta$	$h$	0.000261692	0.000319821	0.00214	0.0032099	0.00828217	0.01536251
0.000	$L(\psi )$	370.4220733	370.2648847	370.3991671	370.1362546	370.2329561	370.1003211
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	370.4220725	370.2648835	370.3991529	370.1362213	370.2328989	370.1001001
	(sec.)	(3.719)	(3.719)	(3.999)	(3.954)	(4.125)	(4.125)
	%APRC	0.000000217	0.000000326	0.000003852	0.00000899	0.00001544	0.00005973
0.001	$L(\psi )$	148.8651477	153.0333562	182.6354932	194.4213028	202.5517767	225.6187725
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	148.8651475	153.0333560	182.6354890	194.4212921	202.5517571	225.6186822
	(sec.)	(3.703)	(3.766)	(3.891)	(3.937)	(4.048)	(4.235)
	%APRC	0.000000114	0.000000173	0.000002270	0.00000550	0.00000968	0.00004003
0.002	$L(\psi )$	93.35627913	96.65270512	121.3940906	132.0288133	139.6031530	162.4531770
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	93.35627905	96.65270499	121.3940885	132.0288077	139.6031426	162.4531260
	(sec.)	(3.686)	(3.751)	(4.031)	(3.907)	(4.140)	(4.109)
	%APRC	0.000000088	0.000000133	0.000001751	0.00000426	0.00000751	0.00003140
0.003	$L(\psi )$	68.10924767	70.73971718	91.01980482	100.0631702	106.6152042	127.0307882
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	68.10924765	70.73971710	91.01980346	100.0631666	106.6151974	127.0307545
	(sec.)	(3.735)	(3.734)	(3.921)	(3.969)	(4.078)	(4.125)
	%APRC	0.000000076	0.000000115	0.000001493	0.00000362	0.00000637	0.00002654
0.005	$L(\psi )$	44.33856516	46.18590815	60.80811349	67.56761750	72.55240610	88.61615789
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	44.33856513	46.18590810	60.80811274	67.56761550	72.55240235	88.61613906
	(sec.)	(3.718)	(3.734)	(3.859)	(3.938)	(4.126)	(4.126)
	%APRC	0.000000064	0.000000097	0.000001232	0.00000296	0.00000518	0.00002124
0.010	$L(\psi )$	23.90454235	24.95647627	33.48319690	37.55530006	40.60863865	50.76516199
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	23.90454233	24.95647625	33.48319657	37.55529918	40.60863702	50.76515391
	(sec.)	(3.844)	(3.720)	(3.890)	(3.875)	(4.123)	(4.110)
	%APRC	0.000000053	0.000000080	0.000000985	0.00000233	0.00000403	0.00001592
0.030	$L(\psi )$	8.801532543	9.191897196	12.42105636	14.00582314	15.21062636	19.31715963
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	8.801532540	9.191897190	12.42105626	14.00582290	15.21062591	19.31715751
	(sec.)	(3.733)	(3.875)	(3.953)	(3.890)	(4.156)	(4.188)
	%APRC	0.000000041	0.000000062	0.000000744	0.00000173	0.00000294	0.00001099
0.050	$L(\psi )$	5.625763186	5.868546426	7.890248987	8.890661342	9.654033441	12.26764061
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	5.625763184	5.868546423	7.890248935	8.890661208	9.654033193	12.26763945
	(sec.)	(3.766)	(3.751)	(3.937)	(3.969)	(4.093)	(4.187)
	%APRC	0.000000036	0.000000054	0.000000650	0.00000151	0.00000256	0.00000945
0.100	$L(\psi )$	3.215706744	3.343586631	4.419555676	4.958260657	5.371177151	6.788110226
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	3.215706743	3.343586630	4.419555654	4.958260598	5.371177042	6.788109719
	(sec.)	(3.766)	(3.766)	(3.906)	(4.015)	(4.156)	(4.093)
	%APRC	0.000000027	0.000000040	0.000000506	0.00000118	0.00000202	0.00000747
0.500	$L(\psi )$	1.326293806	1.355810005	1.620815219	1.761932128	1.872436279	2.254169620
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	1.326293806	1.355810005	1.620815217	1.761932122	1.872436268	2.254169564
	(sec.)	(3.766)	(3.843)	(3.891)	(3.985)	(4.156)	(4.063)
	%APRC	0.000000005	0.000000008	0.000000131	0.00000033	0.00000060	0.00000247
1.000	$L(\psi )$	1.127463864	1.142484646	1.286428785	1.367487603	1.432223859	1.659007178
	(sec.)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)	(<0.001)
	$N(\psi )$	1.127463864	1.142484646	1.286428784	1.367487602	1.432223856	1.659007162
	(sec.)	(3.751)	(3.703)	(3.953)	(3.938)	(4.156)	(4.141)
	%APRC	0.000000001	0.000000002	0.000000043	0.00000012	0.00000022	0.00000098

Next, the effectiveness of the explicit ARL based on the double EWMA chart running on the AR(p) model is compared with the EWMA and CUSUM charts and then investigated using various $\lambda$. Moreover, the double EWMA chart serves as the EWMA chart if instead of ${\lambda }_1$ equals 1. ARL₀ value was fixed as 370 for in-control situation. ARL₁ values were used to compare the efficiency of the EWMA, double EWMA, and CUSUM charts and were set to values mentioned above. The results in Tables 2–4 show that the lower ${\lambda }_1$ shows a decrease in ARL₁, and the double EWMA chart (denoted as DEWMA) showed the effectiveness exceeding the EWMA and CUSUM charts for all situations. Additionally, the AEQL and PCI values are used to confirm their effectiveness. The AEQL values of double EWMA chart have lower than AEQL values of EWMA and CUSUM charts and the PCI values of double EWMA chart with ${\lambda }_1=0.05$ are equals 1, same for AR(1), AR(2), and AR(3).

Table 2. Comparing the ARL values on the double EWMA and EWMA control charts for the AR(1) model with different ${\lambda }_1$ and determined ${\lambda }_2=0.1$ and $\omega =0$.

${\varphi }_1=0.1$					${\varphi }_1=-0.1$
H	0.000261692	0.00236781	0.0946601	4.887	0.000319821	0.00289988	0.1169051	4.343
$\delta$	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	EWMA $({\lambda }_1=1)$	CUSUM $(a=2)$	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	EWMA $({\lambda }_1=1)$	CUSUM $(a=2)$
0.000	370	370	370	370	370	370	370	370
0.001	148.865	185.521	235.416	367.173	153.033	191.464	247.433	367.500
0.002	93.356	123.934	172.794	364.365	96.653	129.269	185.973	364.938
0.003	68.109	93.155	136.596	361.585	70.740	97.684	149.074	362.399
0.005	44.339	62.387	96.425	356.104	46.186	65.767	106.881	357.389
0.010	23.905	34.426	55.878	342.859	24.957	36.459	62.905	345.251
0.030	8.802	12.785	21.436	295.768	9.192	13.576	24.401	301.720
0.050	5.626	8.120	13.619	256.675	5.869	8.618	15.510	265.076
0.100	3.216	4.543	7.517	184.532	3.344	4.811	8.530	265.076
0.500	1.326	1.653	2.445	32.563	1.356	1.723	2.696	38.071
1.000	1.127	1.305	1.771	12.151	1.142	1.345	1.915	13.790
AEQL	0.138	0.164	0.229	2.099	0.140	0.169	0.250	2.449
PCI	1.000	1.187	1.661	15.212	1.000	1.209	1.781	17.463

Table 3. Comparing the ARL values on the double EWMA and EWMA control charts for the AR(2) model with different ${\lambda }_1$ and determined ${\lambda }_2=0.1$ and $\omega =0$.

${\varphi }_1={\varphi }_2=0.1$					${\varphi }_1={\varphi }_2=-0.1$
h	0.000236728	0.00214	0.0852505	5.441	0.000353579	0.0032099	0.130044	4.145
$\delta$	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	EWMA $({\lambda }_1=1)$	CUSUM $(a=2)$	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	EWMA $({\lambda }_1=1)$	CUSUM $(a=2)$
0.000	370	370	370	370	370	370	370	370
0.001	146.814	182.635	230.436	366.813	155.213	194.421	254.613	367.614
0.002	91.766	121.394	167.423	363.590	98.386	132.029	194.181	365.117
0.003	66.850	91.020	131.582	360.403	72.127	100.063	157.031	362.643
0.005	43.460	60.808	92.296	354.130	47.163	67.568	113.727	357.758
0.010	23.407	33.483	53.154	339.028	25.515	37.555	67.627	345.916
0.030	8.618	12.421	20.301	286.056	9.400	14.006	26.430	303.344
0.050	5.512	7.890	12.895	243.033	5.998	8.891	16.804	267.373
0.100	3.156	4.420	7.127	166.402	3.412	4.958	9.217	199.090
0.500	1.313	1.621	2.344	25.057	1.372	1.762	2.857	39.764
1.000	1.121	1.286	1.712	10.225	1.151	1.367	2.005	14.346
AEQL	0.137	0.161	0.221	1.733	0.141	0.173	0.263	2.479
PCI	1.000	1.177	1.613	12.655	1.000	1.221	1.857	17.523

Table 4. Comparing the ARL values on the double EWMA and EWMA control charts for the AR(3) model with different ${\lambda }_1$ and determined ${\lambda }_2=0.1$ and $\omega =0$.

${\varphi }_1={\varphi }_2={\varphi }_3=0.1$					${\varphi }_1={\varphi }_2={\varphi }_3=-0.1$
H	0.000214151	0.001934325	0.076812	3.392	0.000390914	0.00355368	0.14477	2.681
$\delta$	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	EWMA $({\lambda }_1=1)$	CUSUM $(a=2)$	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	EWMA $({\lambda }_1=1)$	CUSUM $(a=2)$
0.000	370	370	370	370	370	370	370	370
0.001	144.844	179.831	225.662	367.859	157.409	197.626	262.971	368.040
0.002	90.235	118.945	162.476	365.536	100.159	134.975	204.032	365.805
0.003	65.638	88.969	127.045	363.232	73.555	102.596	166.765	363.588
0.005	42.616	59.298	88.623	358.681	48.175	69.483	122.290	359.206
0.010	22.930	32.586	50.767	347.626	26.095	38.722	73.671	348.554
0.030	8.442	12.076	19.318	307.639	9.616	14.464	29.071	309.908
0.050	5.402	7.673	12.268	273.517	6.133	9.181	18.488	276.771
0.100	3.098	4.303	6.788	207.714	3.483	5.115	10.105	212.370
0.500	1.300	1.591	2.254	44.948	1.389	1.804	3.054	48.208
1.000	1.114	1.270	1.659	16.236	1.159	1.392	2.111	17.595
AEQL	0.136	0.159	0.213	2.778	0.143	0.176	0.278	2.981
PCI	1.000	1.168	1.570	20.432	1.000	1.234	1.949	20.880

4.2. Real-world data

The prices of diesel and natural gas are shown in Figure 1, with the time series brought up for study in this part using the AR(1) and AR(2) models, respectively. Two datasets, which comprise quarterly price data, were selected as the models by SPSS. The results are in Table 5 shows that two datasets are suitable parameters for the AR model. The Kolmogorov–Smirnov test was then used to determine that white noise significantly fits the exponential mean, as shown in Table 6. For the diesel price, it was fitted to the AR(1), which expressed as $Y_t=3.772+0.943Y_{t-1}+{\zeta }_t\begin{array}{cc},& {\zeta }_t\sim Exp(0.2194)\end{array}.$For the natural gas price, it was fitted to the AR(2), which expressed as $Y_t=16.473+0.373Y_{t-1}-0.256Y_{t-2}+{\zeta }_t\$,\${\zeta }_t\sim Exp(2.8641).$

Figure 1. The time series of datasets (a) Diesel price with AR(1) and (b) Natural gas price with AR(2).

Table 5. The AR(p) coefficients for two datasets from the first quarter 2000 to the second quarter 2022.

Variable	Coefficient	SE	t-Statistic	Sig.
Diesel_price	3.772	0.607	6.214	0.000
AR(1)	0.943	0.039	24.120	0.000
Natural gas_price	16.473	0.424	38.844	0.000
AR(1)	0.373	0.104	3.589	0.001
AR(2)	−0.256	0.104	−2.463	0.016

Table 6. Testing the datasets for the exponential white noise of the exponential distribution.

The datasets	Diesel price	Natural gas price
One-sample Kolmogorov-Smirnov test	0.528	1.181
Exponential parameter	0.2194	2.8641
Asymp. Sig. (2-tailed)	0.943	0.123

Tables 7–8 show that the results of the comparison of control charts between the double EWMA and EWMA charts and the results show the double EWMA chart with the lower ${\lambda }_1$ had a lower ARL₁ and was more effective than the EWMA chart in all situations. However, the performance of CUSUM charts is not compared in this section because CUSUM charts give ARL1 values that are much higher than EWMA and EWMA charts. Moreover, the results show the outcomes of two datasets in a way similar to simulated data, as illustrated in Figure 2. The AEQL and PCI values supported the control chart's effectiveness. The results verify that the double EWMA chart with ${\lambda }_1=0.05$ had the lowest AEQL and PCI equals 1 outperformed the double EWMA chart with higher ${\lambda }_1$ and the EWMA chart, all of which have higher AEQL and PCI > 1, as illustrated in Figure 3.

Figure 2. ARL values of double EWMA and EWMA charts for (a) the AR(1) and (b) the AR(2) models.

Figure 3. The AEQL and PCI values of two models: (a) and (b) for the AR(1), (c) and (d) for the AR(2).

Table 7. Comparing the ARL values on the double EWMA and EWMA control charts for the AR(1) model of the diesel price ($/gallon) with different ${\lambda }_1$ and given ${\lambda }_2=0.1$ and $\omega =0$.

	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	DEWMA $({\lambda }_1=0.2)$	EWMA $({\lambda }_1=1)$
$\delta$	h = 0.0000144156	h = 0.000283249	h = 0.001795748	h = 0.02307312
0.000	370	370	370	370
0.001	132.285	177.741	210.530	248.701
0.002	80.741	117.146	147.296	187.447
0.003	51.112	77.664	101.705	137.031
0.005	34.449	53.721	72.141	101.016
0.010	18.373	29.306	40.358	58.939
0.030	6.779	10.825	15.113	22.722
0.050	4.374	6.889	9.593	14.440
0.100	2.570	3.898	5.358	7.989
0.500	1.186	1.488	1.867	2.561
1.000	1.061	1.212	1.429	1.838
AEQL	0.128	0.151	0.181	0.239
PCI	1.000	1.180	1.421	1.871

Table 8. Comparing the ARL values on the double EWMA and EWMA control charts for the AR(2) model of the natural gas price ($/mcf) with different ${\lambda }_1$ and given ${\lambda }_2=0.1$ and $\omega =0$.

	DEWMA $({\lambda }_1=0.05)$	DEWMA $({\lambda }_1=0.1)$	DEWMA $({\lambda }_1=0.2)$	EWMA $({\lambda }_1=1)$
$\delta$	h = 0.000929664	h = 0.00759138	h = 0.0309506	h = 0.2763760
0.000	370	370	370	370
0.001	152.201	187.520	210.591	235.377
0.002	97.242	127.211	148.887	174.331
0.003	70.883	95.536	114.336	137.527
0.005	46.382	64.291	78.677	97.356
0.010	25.069	35.572	44.425	56.493
0.030	9.224	13.217	16.719	21.670
0.050	5.887	8.390	10.610	13.765
0.100	3.354	4.688	5.890	7.596
0.500	1.358	1.691	2.012	2.465
1.000	1.144	1.326	1.515	1.783
AEQL	0.140	0.167	0.193	0.231
PCI	1.000	1.188	1.377	1.644

Additionally, Figures 4 and 5 also show how well the control chart performs to detect shift changes in the monitoring process. For the AR(1) model in Figure 4, the double EWMA control chart (${\lambda }_1=0.05$) recognises shifts as being out of control at the first observation, whereas the EWMA control chart (${\lambda }_1=1$) does so at the 15th observation. For the AR(2) model in Figure 5, the double EWMA control chart (${\lambda }_1=0.05$) recognises shifts as being out of control at the first observation, whereas the EWMA control chart (${\lambda }_1=1$) does so at the 34th observation. The results show that in the monitoring process, the double EWMA control chart may detect shift changes faster than the EWMA control chart. Besides, An excellent option for spotting change is the double EWMA control chart, which performs better than the EWMA control chart. And then, the exponential smoothing parameter ${\lambda }_1=0.05$ is recommended for the double EWMA control chart.

Figure 4. The effectiveness of detecting shift change in the monitoring process on the (a) EWMA and (b) double EWMA charts for the AR(1) model.

Figure 5. The effectiveness of detecting shift change in the monitoring process on the (a) EWMA and (b) double EWMA charts for the AR(2) model.

5. Discussions and conclusions

The explicit formula for the ARL on the double EWMA control chart was put out in this study. It is faster than the NIE method in computing the ARL for the AR(p) model with exponential white noise. The proposed ARL is compared with the ARL on the EWMA control chart based on in-control and out-of-control scenarios. It was found that the double EWMA control chart's ARL consistently outperformed the traditional control chart's ARL. The double EWMA chart with the lower ${\lambda }_1$ had a lower ARL₁ and was superior to the EWMA chart in terms of effectiveness for all cases. The exponential smoothing parameter ${\lambda }_1=0.05$ is recommended. The efficacy of the control chart was validated by the AEQL and PCI values. Moreover, the performance of the ARL was expanded to include use with real-world datasets. The prices of diesel and natural gas were used as real-world datasets. The outcomes of two datasets are displayed in the findings in a style equivalent to simulated data. However, the double EWMA control chart can be expanded for future study in order to integrate with more models and cover more distributions.

Acknowledgements

The authors would like to express our gratitude to National Research Council of Thailand (NRCT) and King Mongkut’s University of Technology North Bangkok for supporting the research fund with contact no. N42A650318.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

This project was supported by National Research Council of Thailand (NRCT) and King Mongkut’s University of Technology North Bangkok [grant number N42A650318].