Assessing Control Sustainability Using L-Moment Ratio Diagrams

Abstract

This paper presents an application of L-moment statistics and the respective L-moment ratio diagrams (LMRD) to assess control performance, in particular, in terms of control system sustainability. L-moment diagrams are common in extreme events analysis and are considered a very powerful tool in this field at the regional level. Control system assessment is a well-established research area that investigates approaches and methodologies for measuring the quality of control systems. Statistical moments can be used to assess the effectiveness of control systems. The same principle applies to L-moments, with a possible further application to the assessment of control system robustness. The incorporation of the time impact into the analysis allows us to examine the evolution of control systems. In this context, measuring sustainability is only one step away. In this research, L-moments and L-moment ratio diagrams are used to assess the quality of PID-based control systems. In addition, the evolution of their performance over time is depicted visually. Moreover, a robust discordance measure is proposed to measure the robustness, evolution, and sustainability of control systems. The proposed approach is successfully validated using real industrial data obtained from PID basic regulatory control within the hierarchical advanced process control (APC) structure of a large ammonia production plant.

1. Introduction

This article presents results that are a creative combination of two separate contexts: control engineering and life sciences. Control engineering introduces an emerging problem, i.e., the measurement of control system sustainability within its hierarchical layout, including basic PID-based regulatory loops supervised by Model Predictive Control (MPC) dynamic optimization. Control system sustainability is a key challenge in modern control engineering, as systems that experience a decline in control quality can limit and degrade overall process performance. This research specifically focuses on this subject.

Sustainability is a concept that is defined as meeting the needs of the present generation without compromising the ability of future generations to meet their own needs [1]. Besides these high-level perspectives, sustainability should be taken into account at all levels, including in the field of control engineering. One target is to design and develop a control system that meets its efficiency and environmental goals. Simultaneously, these goals must be maintained over time, despite all the fluctuations and uncertainties. This is also referred to as sustainability. The performance should remain at least as good as during the initial commissioning phase. Although control system degradation should be detected and corrected as quickly as possible, this does not always occur. There are many reasons for this, for instance, the limited number of experienced personnel on site and the lack of dedicated tools that support plant staff with information. We cannot allow control systems to degrade. This paper proposes a suitable approach that allows for the measurement of control system performance over time and delivers appropriate information through the visualization of L-moment ratio diagrams.

Modern control systems are organized in a hierarchical structure. The main dynamic controls are organized into basic univariate loops that use the PID algorithm, which accounts for more than 95% of the utilized algorithms [2]. Industrial control systems often perform poorly because of inadequate maintenance, uncoupled disturbances, non-stationarity, non-linearity, failures, incorrect tuning, or human activity.

The majority of the CPA approaches are used in a static way and have a limited ability to incorporate time in a continuous, online assessment. Typically, data are collected, indexes are calculated, and decisions are made. The incorporation of time would allow us to observe the evolution of control systems, i.e., how control performance changes over time. One expects and requires that once tuned, a control loop maintains its performance over time. This is referred to as control system sustainability. The present work aims to address this aspect and develop tools for sustainability analysis.

Existing Control Performance Assessment (CPA) indexes are inadequate and L-moments can improve loop monitoring. They are a robust extension of the classical statistical moments [3], and L-moment ratio diagrams offer a modified version of the standard moment ratio diagrams (MRD). These tools have been successfully applied to the frequency analysis of extreme events in [4], seismology [5] and medicine [6]. In contrast, L-moments have not yet been applied in the control engineering context. A review of the control performance assessment (CPA) literature [7,8] does not reveal any relevant mentions of this subject. This work shows that L-moments can be effectively used to express the quality and evolution of control systems.

Moment ratio diagrams are popular and are commonly used in extreme events analysis. They bring simplicity and clarity to the analysis of extreme events, i.e., observations in the tails of distributions. Therefore, control engineers should also be able to make use of them. The proposed approach takes the idea of LMRDs from extreme events analysis and proposes a new interpretation of the results from a control perspective. Control sciences focus on the most frequent events, and extreme observations (tails of distributions) are considered outliers, abnormalities, or errors. However, they are responsible for degradation. The addition of this perspective to the analysis of control systems can only improve the picture and it has been shown that they offer a promising solution.

LMRD reduces the graphical representation of control system performance to a single point on a two-dimensional plane. The change in control system performance changes the statistical properties of the control error signal, which ultimately changes the point position in the LMRD plane. The positions of the respective points in the LMRD plain not only enable us to compare different loops but also allow us to observe system evolution. Classical MRDs could be used but their efficiency is limited due to normalization problems, tail sensitivity, and the non-Gaussian properties of the variables [8]. Even with small samples and tailed distributions, L-moments can be effective. The idea of using LMRDs to measure control system sustainability is thoroughly validated in a real control rehabilitation initiative. The data used during the analysis were obtained through the APC implementation project [9] at the ammonia production plant. Sustainability is assessed and analyzed over a period of 14 months of plant operation.

This paper begins with Section 2, where the sustainability assessment methods are described. Section 3 presents a case study of an ammonia production plant and the data for analysis. Section 4 presents the main results. Finally, Section 5 concludes this paper with some observations, remarks, and suggestions for further research.

2. Methods and Algorithms

The applied approaches and methods are presented below, i.e., control performance assessment, moments, L-moments, and moment ratio diagrams.

2.1. Control Performance Assessment

CPA efforts started in 1967 with the univariate PID-based loop assessment implemented for a pulp and paper plant using deviation benchmarking. Thereafter, assessment solutions have evolved in many directions to successfully support the industry with mature approaches, measures, and procedures. Current CPA research covers various areas and approaches. Different classifications are used [8,10,11], such as methods that require a plant experiment [12], methods that require the model of the plant [13,14], model-free methods that directly use only plant data [12] and hybrid methods [15,16].

Statistical measures constitute an important part of the CPA research. Control engineering is especially interested in two factors: the time series location and its variability. Data location may be measured using the mean, the median or by other estimators introduced by the Probabilistic Density Function (PDF). The scale informs about time series variability and can be measured using standard deviation or by other estimators. Standard deviation (or variance) is quite an obvious choice as the CPA measure. They are sometimes followed by the higher order statistics like skewness or Kurtosis. Skewness measures the PDF asymmetry, while the Kurtosis describes the data concentration. Therefore, the analysis should often take into account more than one statistical factor [8].

Apart from the aforementioned normal moments, we may utilize statistical factors of other distributions [17,18], variance band index [19] or the robust statistical estimators [20].

The majority of these methods assumes the system to be linear and Gaussian. The opposite situations are rarely addressed, especially when we consider the fact that industrial systems are complex, nonlinear and they may exhibit strange properties.

2.2. Statistical Moments

Statistical analysis can be performed in many ways. It can follow a theoretical way assuming some distribution function reflecting the underlying stochastic process. Thus, the proper PDF can be used to evaluate its statistical factors and moments, if they exist. The second way is to approach the problem empirically without specific assumptions about the PDF. In such a situation, one has to estimate values of the statistical moments empirically.

Let ${\left\{X_i\right\}}^T$ denote the time series with mean $\mu$ and r-th central moment ${\gamma }_r=E{\left(X-\mu \right)}^r$, where $E\left(·\right)$ denotes the expectation operator. Variance ${\sigma }^2$ represents the second central moment denoted as ${\gamma }_2$ ($\sigma$ standard deviation). The first two moments are often supplemented by the third moment ${\gamma }_3$ called skewness (asymmetry measure) and the fourth one ${\gamma }_4$, the kurtosis (concentration measure).

2.3. L-moments

The theory of L-moments was proposed by Hosking [3] as a linear combination of order statistics. This approach significantly improves conventional methodologies as L-moments introduce new estimates for the statistical moments that may help to better reflect distribution properties and compare different PDF candidates. We define L-moments for any random variable, whose expected value exists. Unlike central moments, L-moments provide almost unbiased statistics, even for a small sample. Moreover, they are less sensitive to the existence of the distribution tails [21]. The above properties have been fully appreciated in life sciences, although observation of the time series of the industrial control systems does not provide grounds for limiting their applicability.

Statistical properties of the data are reflected by L-location (L-shift) $l_1$, L-scale $l_2$, coefficient of L-variation (L-Cv) ${\tau }_2\in \left.〈0,1\right)$, L-skewness ${\tau }_3\in \left(-1,1\right)$ and L-kurtosis ${\tau }_4\in \left.〈-1/4,1\right)$. They are used to fit a distribution to a dataset by the comparison of the empirical L-moments to their theoretical counterparts. L-skewness and L-kurtosis are well suited as the goodness-of-fit measures. They can be analytically evaluated for different PDFs [22].

2.4. Moment Ratio Diagrams

Moment ratio diagrams are thought to have first been introduced by Craig [23] and originated from Karl Pearson’s works in the early XIX century. They graphically present statistical properties of some time series. They are able to perform theoretical PDF fitting, with a comparison of the functions shapes or their classification [24]. The MRD is a graphical representation in a Cartesian coordinate of a pair of standardized moments. Actually, there are two versions [25]. The classical moment-ratio diagram MRD(${\gamma }_3,{\gamma }_4$) shows the third standardized moment ${\gamma }_3$ (or its square ${\gamma }_3^2$) as abscissa and the fourth one ${\gamma }_4$ as ordinate, plotted upside down. There exists theoretical limitation of the accessible area, as ${\gamma }_4-{\gamma }_3^2-1\ge 0$. The locus corresponding to PDF can be a point, a curve or a region. It depends on the number of shape parameters. PDFs lacking shape factor (such as Gauss or Laplace) are represented by a point. Functions with one shape parameter are reflected by a curve. Regions characterize functions with two factors. The second type of the diagram MRD(${\gamma }_2,{\gamma }_3$) proposed by Cox and Oaks [26] represents variance ${\gamma }_2$ as the abscissa and the skewness ${\gamma }_3$ as the ordinate. Though they are location and scale dependent, they gained acceptance.

The MRDs are used to measure the proximity between univariate PDFs, to show distribution versatility, to choose the best fit PDF, to clarify the relations between function families and to identify homogeneous processes.

2.5. L-moment Ratio Diagrams

L-moments were proposed by Hosking [3] in 1990. It has been shown that they outperform other estimation methods. L-moments and LMRDs can be used to identify proper distribution for empirical observations. The plot LMRD(${\tau }_3,{\tau }_4$) shows the L-kurtosis ${\tau }_4$ versus L-skewness ${\tau }_3$. Similarly to MRDs, one can confront the empirical data with the theoretical functions [21]. A sample plot with theoretical shapes is shown in Figure 1. Such curves can be used to compare empirical observations with theoretical functions.

Similarly to MRD(${\gamma }_2,{\gamma }_3$), its L-moment equivalent might be used as well. Apart from PDF fitting, LMRDs can be used to compare samples of different origins. In hydrology these are river discharges and their homogeneity [4], while in automation control errors of the PID loops. The homogeneity testing task is performed on a two-dimensional plane spanned by L-skewness and L-kurtosis. Homogeneity testing may be performed in a way that the visual inspection is supported by the dedicated measures [27].

Let ${\mu }_V$ and ${\sigma }_V$ denote the mean and the standard deviation of the test statistics V shown in Equations (1)–(3)

$$V_1=\sqrt{\frac{{\sum }_{i=1}^N{n}_i{\left({\tau }_2^i-{\tau }_2^R\right)}^2}{{\sum }_{i=1}^N{n}_i}},$$

(1)

$$V_2=\frac{{\sum }_{i=1}^N{n}_i\sqrt{{\left({\tau }_2^i-{\tau }_2^R\right)}^2+{\left({\tau }_3^i-{\tau }_3^R\right)}^2}}{{\sum }_{i=1}^N{n}_i},$$

(2)

$$V_3=\frac{{\sum }_{i=1}^N{n}_i\sqrt{{\left({\tau }_3^i-{\tau }_3^R\right)}^2+{\left({\tau }_4^i-{\tau }_4^R\right)}^2}}{{\sum }_{i=1}^N{n}_i},$$

(3)

where $n_i$ is the site i record length, ${\tau }_2^i,{\tau }_3^i,{\tau }_4^i$ are the sample L-moments and ${\tau }_2^R,{\tau }_3^R,{\tau }_4^R$ are regional mean L-moments.

They are derived from a large number of simulated replications from that region of data, having a four-parameter Kappa distribution as the frequency distribution.

The values of the $H_i$ statistics (4) indicate the region to be acceptably homogeneous, when $H_i<1$, possibly, when $1\le H_i<2$, and definitely heterogeneous in case of $H_i\ge 2$

$$H_i=\frac{V_i-{\mu }_V}{{\sigma }_V}.$$

(4)

2.6. Discordance Measure

The homogeneity $H_i$ statistics do not indicate the observation, which causes heterogeneity. The discordance measure $D_i$ [27] in terms of L-moments ${\tau }_2^i,{\tau }_3^i,{\tau }_4^i$ has been proposed to address this issue. It identifies unusual distributions that are generally discordant with the whole set of sources. In its original definition, distributions reflect river discharge sites. PDFs of the variables reflecting loop performance are the data sources in our case. We compute the $D_i$ measure for each data source separately, keeping in mind that the data are represented by three dimensional vectors ${\mathit{u}}_i=\left[{\tau }_2^i,{\tau }_3^i,{\tau }_4^i\right]$, $i=1,\dots ,N$:

$$\begin{array}{cc}\hfill D_i^2\left({\mathit{u}}_i\right)& ={\left({\mathit{u}}_i-\overline{\mathit{u}}\right)}^T{\mathit{S}}^{-1}\left({\mathit{u}}_i-\overline{\mathit{u}}\right),\hfill \end{array}$$

(5)

$$\begin{array}{cc}\hfill \overline{\mathit{u}}& =\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$N$}\right.\sum _{i=1}^N{\mathit{u}}_i\hfill \end{array}$$

(6)

$$\begin{array}{cc}\hfill \mathit{S}& =\frac{1}{(N-1)}\sum _{i=1}^N\left({\mathit{u}}_i-\overline{\mathit{u}}\right){\left({\mathit{u}}_i-\overline{\mathit{u}}\right)}^T\hfill \end{array}$$

(7)

where $\overline{\mathit{u}}$ denotes the sample mean and $\mathit{S}$ the group (region) covariance matrix. Discordance measure $D_i$ represents how far a certain vector ${\mathit{u}}_i$ is from the region center. Large values mean that the the data source does not follow the major pattern. This task is similar to outlier identification [28], although in a multivariate space using the Mahalanobis distance [29]. The only thing left is to define the threshold, which will delimit the discordant observations. Assuming that vectors ${\mathit{u}}_i$ are drawn from the multivariate Gaussian population, the $D_i^2$ follows approximately 3 degrees of freedom ${\chi }^2$ distribution. In such a case, the square root of the 0.975 quantile of the ${\chi }^2$ PDF threshold is used, i.e., $d_0=\sqrt{{\chi }_{3,0.975}^2}\approx 3.06$.

We have to keep in mind that mean and covariance estimators are sensitive to the outliers [30] and it may be significantly biased. Few outliers may attract the mean $\overline{\mathit{u}}$ and inflate the covariance in their direction, thereby protecting their discovery (masking effect). The solution is to use robust estimators with a higher breakdown point, optimally equal to $50\%$ [31]. In our case, the FAST-MCD (Minimum Covariance Determinant) method is used. The FAST-MCD algorithm [32] looks for n observations with the smallest scatter and is characterized by a breakdown point of almost 50%, i.e., the best that can be achieved. The estimator is defined for a subset of $M<N$ observations for m-variate data, whose covariance matrix has the smallest determinant out of all possible subsets of length M. The MCD location T and scale S are:

$$\begin{array}{cc}\hfill \mathit{T}& =\frac{1}{M}\sum _{j=1}^M{\mathit{u}}_{i_j}\hfill \end{array}$$

(8)

$$\begin{array}{cc}\hfill \mathit{S}& =c_m\frac{1}{M}\sum _{j=1}^M\left({\mathit{u}}_{i_j}-\mathit{T}\right){\left({\mathit{u}}_{i_j}-\mathit{T}\right)}^T.\hfill \end{array}$$

(9)

The subset size M can be chosen in many ways, however, $M=⌊N+m+1⌋$ is suggested. Matlab algorithm implementation is used for the FAST-MCD computation.

2.7. LMRD in Control Sustainability Task

The above methods and algorithms are successfully used in life sciences. An analogous application in control engineering does not exist (at least, one has not been found). The possibility to represent control loop statistical properties as a single point on a plane spanned by moments seems highly attractive. Moments, such as the mean value, standard deviation, kurtosis or skewness are used as common measures in the CPA [8]. Robust L-moments offer an even more attractive alternative, as they are robust against outliers and are effective for short datasets.

The ability to present loop properties in the form of a single point in the LMRD plan provides further opportunity. Among other possibilities, we may analyze and fit underlying empirical statistical process to the known theoretical ones. In case of many loops, we may compare them, check their homogeneity or identify the outlying ones. These tasks may be supported by various algorithms, but even visual inspection is informative for practitioners. The proposed use of the LMRD diagrams goes further, as it attempts to incorporate time into the study. Time introduces a new degree of freedom. Since we are only able to see the variation of statistical properties over time, we may attempt to measure system evolution. Once the system does not evolve, its statistical properties remain constant and it sustains its performance. Thus, we obtain a tool to measure control system sustainability.

The travel of a point in the LMRD plane can help to visualize the phenomenon, but it is not enough. One needs a quantitative measure, a number saying that the travel is small enough to be considered as sustaining its properties. The homogeneity test might solve the issue. However, the result of such a test is not practical enough. It would say that the system is homogeneous (which is good) or heterogeneous (which is bad), but in the second case we just do not know why is it so. Discordance analysis answers both questions. When no discordance is found, the system is considered homogeneous. In the opposite case, we have indication, which observation is responsible and to what extent. As a result, a root-cause analysis can be performed. A custom measure is proposed to better reflect the concentration of the data. A robust data center is evaluated as a two dimensional geometric median (GeoMed) ${\overline{x}}_{med}$ defined as a value of the argument $x_0$ from where the sum of all Euclidean distances to points $x_i$ is minimized

$${\overline{x}}_{med}=\underset{x_0\in {\mathbb{R}}^m}{arg\; min}\sum _{i=1}^N{\parallel x_i-x_0\parallel }_2,$$

(10)

where N is the number of points. It is evaluated with Weiszfeld’s algorithm [33]. Once the center is known, the range of surrounding points can be obtained. We use two-dimensional modification of the median absolute deviation around median $r_{cc}$ and maximum distance $r_{max}$ to compare different datasets

$$\begin{array}{cc}\hfill r_{cc}& =1.5·\underset{i}{median}\left|x_i-{\overline{x}}_{med}\right|,\hfill \end{array}$$

(11)

$$\begin{array}{cc}\hfill r_{max}& =\underset{i}{max}\left|x_i-{\overline{x}}_{med}\right|.\hfill \end{array}$$

(12)

The above methods can be used in the following sustainability assessment methodology (Figure 2) that can help to assess any control system’s sustainability:

• Select periods of the consistent system operation, i.e., of the comparable plant load and performance. This choice will justify that we will be evaluating comparable and appropriate plant modes. One can also imagine a situation where a separate analysis is conducted for different operating regimes, such as different loads.
• Collect the respective time series of the loop control errors, as the basis for evaluations.
• Calculate L-moments for the considered data.
• Plot LMRD for the selected loops.
• Calculate discordance measures.
• Identify and label discordant observations.

It is worth noting that the ability to measure and analyze the sustainability of a control system enables its improvement, which entails facilitating the work of the control engineer. Correct information on the evolution of control system quality helps to to extend its life, as well as to maintain the highest possible efficiency of the controlled plant.

3. Ammonia Plant APC Implementation

3.1. Plant Description

Hydrogen is indispensable for ammonia synthesis. Grupa Azoty Zakłady Azotowe “Puławy” SA produces ammonia through the autothermal reforming of methane $C{H}_4$ with the utilization of oxygen: pure and from the air (see Figure 3).

The hydrogen preparation for further processing follows the subsequent sub-processes:

• Methane conversion (raw materials heating: natural gas, process air, oxygen and 3.2 MPa steam) in the preheaters; next autothermal methane reforming,
• Carbon oxide conversion in shift reaction,
• CO${}_2$ removal in Benfield unit by the absorption in hot potassium carbonate and activator solution,
• CO and CO${}_2$ residuals removal from process gas: Copper-Ammonia Cleaning.

The produced hydrogen catalytically reacts with nitrogen (obtained from the air) to form anhydrous liquid ammonia. Such a synthesis is a balanced process, where gases leaving the reactor consists of an approximately $17\%$ volume of ammonia. The process runs in the synthesis loop. Circulating synthesis gas mixes with the fresh one. An ammonia synthesis reaction occurs in the ferric catalyst, while the generated heat heats up the steam and preheats the entering gas. The process runs at a pressure of approximately 28–30 MPa. Produced ammonia is condensed in heat exchangers and chillers after the reaction. Liquefied ammonia is gathered in separators, decompressed and sent to further installations.

3.2. Control System Layout

The implemented APC system follows an industry standard hierarchical control layout with a PID-based basic regulatory layer supervised by the MPC dynamic optimization [34]. Existing PID base regulatory control loops take responsibility for the dynamic responses, i.e., for the realization of setpoint demands. APC is implemented above and it operates with longer time sampling and horizons. Supervisory APC focuses on the evaluation of setpoint signals and eventually of the controller output biases, considering overall process dynamics defined by longer time constants than basic PID loops. In the considered application, the Honeywell APC supervisory system uses a single MPC controller decomposed into two subcontrollers responsible for the process of gas preparation and ammonia synthesis, respectively. It is integrated with the Emerson deltaV Distributed Control System (DCS) via the OPC (OLE for Process Control) link. APC on/off switching is performed using a bumpless procedure specifically designed for the system safety. Therefore, the regulatory PID loops participating in the APC coordinated work gain a new mode of operation. Apart from standard manual (MAN) and automatic (AUTO) modes, there is the Remote CAScaded mode (RCAS), which uses optimal remote setpoint values from the upper-level MPC (Manipulated Variables—MVs).

The APC system is implemented using a proven methodology. The implementation starts and concludes with the plant and control system assessment to measure and confirm the benefits [9]. A positive outcome of APC commissioning is required. However, it is only the starting point for the system operation. It should work and maintain its performance for many years. Thus, positive results must be sustained. This work focuses on that subject and proposes the methodology to assess control system sustainability. Exemplary industrial application demonstrates the effectiveness of the approach.

4. System Sustainability Analysis

The analysis follows the procedure proposed in Section 2.7.

4.1. Data Selection

Data selection is crucial for further analysis. Improperly selected and, worst of all, incomparable data may result in the collection of erroneous results. Much attention must be paid during the data-review process. The procedure should follow one main assumption: the data must be comparable. They should reflect the same normal operating regimes. Though this process is crucial for the credibility and correctness of further analysis, authors often do not pay attention to it. In the considered case, the data from 14 months of operation (07/2020–08/2021) of the installation are collected with a 1 min sampling interval. Much attention is paid to the selection of the comparable time periods, characterized by the similar and constant throughput. A high plant load is selected. Natural gas consumption varies on a monthly basis $\pm 1.5\%$. The load is constant during each period and its variations are only process related, as shown in Figure 4. In our analysis, we had to keep in mind the fact that we are working on real processes with real data. Therefore, there was a level of unpredictability. To mitigate this, (or at least to minimize it significantly) we have also excluded periods during which hardware problems were detected to ensure that all possible external effects do not appear.

As is visible in the sample time trend, the data are subject to outlying observations, which justifies the selection of robust statistical estimators.

4.2. Data Collection

The data are collected from DCS and SCADA (Supervisory Control And Data Acquisition) systems. The DCS stores process data, while weather data are kept in the SCADA. There are 12 loops for flow control (denoted as FLOW_1 … FLOW_12), 2 for level control (LVL_1 and LVL_12), 1 for pressure (PRESS_1) and 7 for temperature (TEMP_1 … TEMP_7), amounting to 22 control loops overall. A more detailed loops description is limited by the Company’s data security rules. Weather data include ambient temperature, pressure and humidity. Air density is evaluated using the weather data.

4.3. L-moments Evaluation

At first, L-moments are evaluated for each loop and for each time period. As we have 14 datasets and 22 loops, there are 308 time series that need to be evaluated. The following statistics are calculated: median, mean ${\gamma }_1$, standard deviation ${\gamma }_2$, skewness ${\gamma }_3$, kurtosis ${\gamma }_4$, L-shift $l_1$, factor $l_2$, L-scale ${\tau }_2$, L-skewness ${\tau }_3$ and L-kurtosis ${\tau }_4$. Showing all values is not crucial. The sample set for the flow control loop during 07/2021 is shown in

Table 1. Statistics for sample control loop during 07/2021.

normal	${\mathit{\gamma }}_1$	median	${\mathit{\gamma }}_2$	${\mathit{\gamma }}_3$	${\mathit{\gamma }}_4$
	0.1334	−0.0322	5.2014	4.6643	97.0951
L-moments	$l_1$	$l_2$	${\tau }_2$	${\tau }_3$	${\tau }_4$
	0.1334	2.3418	17.5593	0.0593	0.1957

. Data contain outliers: the median is close to zero in contrast to the mean.

4.4. LMRD Plots

Presentation of pure numbers for statistical factors does not help much. Therefore, moment ratio diagrams are evaluated. The analysis is derived for a single sample loop, for the same flow controller as in the previous chapters. At first, normal moment ratio diagrams are presented in Figure 5.

Note that interpreting the individual graphs for each loop requires caution, as shown below. Next, the same loop is analyzed, but using L-moment ratio diagrams. Two standard LMRDs are plotted, i.e., standard kurtosis versus skewness LMRD(${\tau }_3,{\tau }_4$ vs. L-Cv LMRD(${\tau }_2,{\tau }_3$) diagram. The LMRD($l_2,{\tau }_3$) L-skewness versus $l_2$ diagram is compared in Figure 6b.

The review of the above plots highlights several observations. The LMRD(${\tau }_3,{\tau }_4$) plot highlights interesting information about the loop behavior (see Figure 7). It shows that the loop operates in the same way during all time periods except the first one and the 12th (except July-2020 and June-2021). This unusual behavior is caused by unsteady process parameters that failed to return to normal after an emergency shutdown just prior to the data periods under consideration. All other points are well grouped around a theoretical point depicting normal distribution, i.e., $(0,0)$. As there is a common opinion that a good loop should be characterized by normal distribution (stationary with all disturbances are properly rejected), it is a positive observation. The other two plots show similar behavior, though they are scale dependent. Additionally, the loop maintains variability (information represented by scaling), except in month #12, i.e., June-2021.

This example reveals interesting observations about loop properties and performance sustainability. Further analysis should aim to determine the root-cause for such a behavior. It is suggested to add the third degree of freedom to LMRD diagrams. During discussions with plant technology staff, potential causes of the observed phenomenon were proposed. First off all, no loop tuning happened during these periods, so the human activity reason is rejected. The other reasons can be technology dependent. The process uses atmospheric air to obtain nitrogen and oxygen. It must be remembered that the amount of oxygen or nitrogen in the volume of air depends on its temperature, humidity and pressure, i.e., actual air density. Finally, there might be slight changes in the installation load (natural gas consumption). A specific change to the LMRD plot is proposed. The color of the circles depicting a single month may depend on an external variable. Five selections of such variables are investigated: natural gas consumption $F_{\mathrm{gas}}$, ambient air temperature $T_{\mathrm{air}}$, pressure $p_{\mathrm{air}}$, humidity $h_{\mathrm{air}}$ and density ${\rho }_{\mathrm{air}}$. Figure 8 presents the respective plots. We notice that two outlying months are characterized by low humidity.

There is a period of month #13 (July-2021) that is characterized by low humidity. We notice that both the load $F_{\mathrm{gas}}$ and ambient temperature $T_{\mathrm{air}}$ can be used to differentiate these situations. Thus, low humidity with a high load and high temperature can cause performance deterioration in the considered case. A similar analysis has been performed for other loops, thereby confirming the initial observations. However, the visual interpretation of the LMRD plots is not sufficient. One would expect to have a more automatic tool determining whether the points in the LMRD are homogeneous. For that, the discordance measure analysis is performed with its results described below.

4.5. Discordance Measure Evaluation

Discordance measures have been evaluated for all considered loops and time periods. Sample plots for loop #1 are shown in Figure 9. The observation of the disclosed behavior confirms and extends observations derived from LMRD plots. What is even more important, the observations considered to be outlying (months #1 and #12) are not considered as outlying observations by the discordance measure D. As one can see, relying only on visual inspection might be misleading.

Table 2. Discordance analysis results for all the loops.

		Month Number
		1	2	3	4	5	6	7	8	9	10	11	12	13	14
ambient	temp	17.2	19.8	12.7	10.7	7.5	1.8	4.0	4.9	5.4	7.6	8.7	23.6	22.8	17.6
	press	1001.8	992.9	993.7	997.7	1014.8	1003.9	981.0	1015.4	1005.1	1000.9	990.9	999.7	997.3	1001.8
	humi	76.8	90.0	100.0	97.9	99.1	99.7	95.6	95.1	82.7	95.0	85.7	81.1	83.3	89.3
	dens	5.16	2.47	5.34	6.41	7.76	9.41	8.63	8.81	9.06	7.76	7.75	1.25	1.41	3.87
Loop Number	1	2.35	1.23	1.31	1.05	1.16	0.95	1.15	0.99	0.95	1.13	0.95	2.22	1.12	0.80
	2	0.81	2.01	1.02	0.77	1.16	0.54	2.25	1.05	1.19	3.39	2.23	3.30	0.70	1.11
	3	0.81	1.22	2.66	2.70	1.25	1.08	1.02	0.96	2.30	3.01	0.91	0.84	0.99	0.90
	4	1.23	1.01	0.67	2.21	0.90	1.09	1.29	1.37	0.85	3.97	1.21	0.71	1.01	1.05
	5	1.24	0.96	3.46	0.84	2.19	2.67	0.74	1.06	1.26	0.83	0.96	10.05	1.06	0.89
	6	1.09	0.90	1.28	0.91	0.93	3.04	1.33	1.22	1.17	0.90	2.17	0.96	1.09	0.90
	7	0.89	1.20	1.08	3.17	1.29	2.35	1.03	0.78	1.04	0.74	2.82	0.81	1.01	1.26
	8	0.97	4.84	1.06	1.16	5.73	3.68	0.99	0.99	0.62	2.14	2.67	0.76	0.98	1.10
	9	1.16	2.20	0.75	1.14	3.16	2.43	0.86	1.21	0.93	1.02	0.81	2.68	0.82	1.17
	10	0.80	0.82	1.20	3.44	1.72	1.29	0.66	2.13	0.98	1.07	1.16	1.10	1.07	1.02
	11	1.36	0.96	1.17	0.89	0.92	0.90	1.38	0.71	1.81	3.80	0.54	1.01	1.16	1.19
	12	1.89	2.43	4.72	1.11	2.52	0.63	0.98	1.08	1.05	1.19	5.00	1.01	0.72	0.71
	13	1.69	0.86	1.06	2.12	1.25	1.00	1.13	1.13	0.69	1.15	5.49	0.85	2.27	0.50
	14	1.21	1.09	1.00	0.85	3.56	1.22	0.90	4.01	0.77	0.92	0.94	4.05	1.97	1.05
	15	0.90	1.31	1.34	0.91	1.35	1.06	0.96	1.20	0.95	6.10	0.99	1.15	1.16	0.93
	16	1.23	0.89	1.58	0.75	1.02	2.22	0.84	1.19	1.25	2.99	1.07	1.02	0.85	1.12
	17	3.37	1.18	2.40	1.19	1.22	0.79	0.84	1.03	1.06	3.18	1.13	0.81	0.81	1.17
	18	2.44	3.76	2.11	1.27	0.99	1.01	1.20	1.15	0.79	2.56	0.58	0.99	0.81	0.95
	19	1.13	2.14	1.21	0.90	0.99	1.04	1.07	4.69	1.00	1.15	0.98	0.99	3.31	1.07
	20	1.18	1.92	1.00	0.62	2.28	0.96	2.30	0.84	0.82	7.72	3.30	1.16	0.77	1.14
	21	1.21	1.04	0.93	1.01	1.00	1.11	1.23	0.97	0.82	0.85	4.48	1.61	1.19	1.91
	22	0.72	1.02	1.05	1.20	1.20	0.95	1.13	0.78	2.07	2.55	2.77	0.91	0.94	14.48

shows all the results, while

Table 3. The results of outlier detection.

Loop	Outliers	Month	Loop	Outliers	Month
FLOW_1	0		FLOW_12	1	3
FLOW_2	2	10, 12	LVL_1	0
FLOW_3	0		LVL_2	3	5, 8, 12
FLOW_4	1	10	PRESS_1	0
FLOW_5	1	3	TEMP_1	0
FLOW_6	1	6	TEMP_2	2	1, 10
FLOW_7	1	4	TEMP_3	1	2
FLOW_8	2	2, 6	TEMP_4	2	8, 13
FLOW_9	1	5	TEMP_5	1	11
FLOW_10	1	4	TEMP_6	1	11
FLOW_11	1	10	TEMP_7	0

summarizes detection using the discordance measure. Each loop has different characteristics. It seems that the derivation of conclusions requires greater attention and comparisons between LMRDs and discordance diagrams.

4.6. Concentration Analysis

The analysis cannot be performed in isolation from reality, i.e., the properties of the plant and the control system. One option is a mathematical concentration analysis, while the technological process information helps to place these observations in a real-world context. This is clearly observable once we compare three loops that exemplify the issue: loop #7 (FLOW_7), #14 (LVL_2) and #15 (PRESS_1).

Figure 10 shows the LMRD(${\tau }_3,{\tau }_4$) plot and discordance measures diagram for loop #7. We notice that points are very well concentrated and the cluster is very close to the point depicting normal distribution $(0,0.1226)$. Visual inspection of the LMRD plot can help to draw a conclusion that this loop sustains its features. However, the discordance measure points out one outlying month. The concentration measures can be taken into account. Actually, the evaluated concentration radius $r_{cc}=0.025$ and maximum distance from the geometric median $r_{max}=0.043$ are significantly low. Given the above, one can safely determine that this loop sustains its performance.

Loop #15 (PRESS_1), shown in Figure 11, seems to be totally opposite. Points in the LMRD(${\tau }_3,{\tau }_4$) plot show very high dispersion, but the discordance measure shows safe values and none observation is considered as an outlier. The explanation for this is quite simple. All the points are equally spread; therefore, they exhibit similar homogeneous behavior. None of them can be pointed out as being discordant. With these observations in mind, let us now look at the results of the concentration analysis. The concentration radius is equal to $r_{cc}=0.383$ and the maximum distance from the geometric median is $r_{max}=0.566$. As one can see, both values are two orders of magnitude higher than for the previously discussed loop. It is clear that this loop does not sustain its performance.

Loop #14 (LVL_2), which is sketched in Figure 12, is considered as the third example. Its behavior can be positioned somehow between the above two examples, but it shows another interesting feature.

Observed points are characterized by a fairly significant dispersion reflected in $r_{cc}=0.179$ and $r_{max}=0.536$. The discordance analysis determined three outlying observations. These discordant loops raise concerns about their sustainability. The above investigation focuses on a purely mathematical aspect. It is necessary to confront these facts with the actual loop properties and the process conditions. One should be aware that the mode of loop operation, i.e., MAN, AUTO, CAS or RCAS should impact the observed results. The majority of loops is in one of the types of automatic modes, i.e., AUTO, CAS or RCAS. In opposition, one loop (PRESS_1) is all the time kept in the manual mode. This fact clearly explains the results of the analysis. Large points scattering, expressing changing statistical properties of the loop, occurs due to the open loop control (the output of the controller is manually set by the human operator). One can be tempted to conclude that L-moments coincide with process conditions of the system.

Loop #7 operates continuously in a closed-loop fashion, similarly to Loop #1 (FLOW_1), which is fully visualized in previous chapters. The applied well-tuned controller keeps the process variables in check, and as can be seen from the small scattering of points on the LMRD plane, the elapsed time does not degrade the results. As one can see, the extreme properties are represented quite well by the proposed approach. As always, the observations that fall between the extremes are the most difficult to interpret. There is a group of loops characterized by a switching between the following two modes: RCAS and AUTO. As a result, operators decided to exclude them from supervisory APC control (RCAS). Moreover, these loops were subject to the MAN mode operations, which confirm the drawn conclusions. It is interesting to notice that even short MAN mode operation (3.4% or 8.7% of time) is noticeable in the analysis. There are three more aspects that should be taken into account, namely the loop type, its goodness and the time.

Table 4. Concentration analysis for the loops (green color denotes the closest loops, red color the most discordant ones).

Loop	Type	Mode	${\overline{\mathit{x}}}_{\mathit{m}\mathit{e}\mathit{d}}$		${\mathit{r}}_{\mathit{c}\mathit{c}}$	${\mathit{r}}_{\mathit{m}\mathit{a}\mathit{x}}$	${\mathit{r}}_{\mathit{G}}$
			${\mathbf{\tau }}_{\mathbf{3}}$	${\mathbf{\tau }}_{\mathbf{4}}$
1	FLOW	CAS	0.006	0.168	0.060	0.131	0.051
2		CAS	0.010	0.203	0.045	0.107	0.091
3		AUTO	0.009	0.181	0.038	0.047	0.067
4		CAS	−0.006	0.160	0.073	0.119	0.044
5		CAS	0.005	0.089	0.037	0.084	0.038
6		RCAS	−0.022	0.135	0.054	0.079	0.035
7		CAS	−0.014	0.125	0.025	0.043	0.017
8		RCAS/AUTO	0.000	0.145	0.137	0.367	0.022
9		RCAS	−0.040	0.111	0.064	0.073	0.052
10		RCAS	−0.030	0.098	0.054	0.109	0.055
11		AUTO	0.017	0.130	0.018	0.042	0.025
12		AUTO	0.026	0.146	0.035	0.224	0.049
13	LVL	RCAS/AUTO	0.018	0.239	0.187	0.331	0.134
14		RCAS/AUTO	−0.020	0.163	0.179	0.539	0.061
15	PRESS	MAN	0.157	0.159	0.383	0.566	0.194
16	TEMP	RCAS	0.004	0.243	0.123	0.181	0.124
17		RCAS	0.025	0.145	0.041	0.383	0.046
18		RCAS	−0.007	0.155	0.038	0.067	0.040
19		RCAS	−0.011	0.204	0.125	0.173	0.093
20		AUTO	−0.012	0.156	0.091	0.202	0.046
21		RCAS	0.021	0.254	0.313	0.352	0.153
22		RCAS	−0.080	0.153	0.073	0.120	0.110

summarizes all the results from the concentration analysis. It shows the type of control, the dominating control mode, the geometric median ${\overline{x}}_{med}$ for the LMRD(${\tau }_3,{\tau }_4$) plot, the concentration radius $r_{cc}$, the maximum outlying distance to the geometric median $r_{max}$ and the $L_1$ distance $r_G$ of the ${\overline{x}}_{med}$ to the point representing normal PDF.

Experience shows that the type of controlled variable (flow, temperature, pressure or level in the tank) matters. We see this in the data. The observed level of control is especially challenging. It is a completely new aspect requiring further research. A small number of loops cannot be used to draw reliable opinions, but the problem is worth further consideration.

Loop goodness is the second issue. There is a common agreement that the well controlled loop exhibits normal properties. It means that all the disturbances are decoupled, the process is stationary, linear and is not impeded by any external impact. The LMRD(${\tau }_3,{\tau }_4$) plot is therefore useful. A normal PDF is represented by the point $(0.0,0.1226)$. Measuring the distance of the geometric median to that point is useful for determining the loop goodness. Table 4 shows that the loop goodness varies. Green depicts the three shortest distances, while red indicates the worst three. The loop operating in the MAN mode is furthest from the normal PDF. The scattered loops show large distance as well. The loops that are well clustered are the closest. The $L_1$ distance $r_G$ from the geometric median ${\overline{x}}_{med}$ to the normal PDF point is proposed for the new CPA measure, which is worth further investigation.

The results show that there is no connection between the discordance $D_i$ and concentration measures: ${\overline{x}}_{med}$, $r_{cc}$ and $r_G$. We do not observe performance degradation in time. The variations are dependent on the process itself or on external weather conditions. The control system sustains its properties.

4.7. Discussion of the Results

The analyses are performed for the real plant data during the 14 months of its operation. Two subjects are addressed: the loop performance change, i.e., the loop control sustainability and the loop homogeneity/discordance monitoring. It is shown that the LMRD(${\tau }_3,{\tau }_4$) diagram can be used to demonstrate how the loop control quality changes in time. Once we assume that the loop is initially properly tuned we may observe whether this quality is maintained. Once the point representing the loop starts to travel along the diagram, we may determine that something wrongs happens to that loop.

A homogeneity analysis can be used to compare the loops. We may compare loops and detect the discordant one. Moreover, this analysis may be performed in time. It can even be used to detect the impact of the loop operational mode when it switches to MAN. Therefore, the LMRD diagrams can be used in two separate instances: to observe loop performance fluctuation in time and to compare different loops.

5. Conclusions and Further Research

The paper addresses a very important aspect of the control performance assessment—its eventual degradation in time. This phenomenon is called control system sustainability. The subject is of a very high practical significance, but it rarely appears in the research. There are no common measures for control system sustainability. A novel approach to the subject is proposed using the measures derived from the statistical L-moment ratio diagrams. These tools are common in the life sciences extreme analysis, though unprecedented in control engineering.

It is shown that the proposed methodology using the LMRD(${\tau }_3,{\tau }_4$) plot accompanied with the discordance and concentration measures addresses real industrial issues. This valuable observation reveals new research opportunities. The paper shows that we may observe control performance change in a single LMRD plot. If the system works properly, we may observe whether its control quality changes or not. Assuming that the initial tuning is proper, the observed change is towards worse performance. Thus, we may monitor and indicate that the loop performance changes (degrades). However, two important issues associated with the proposed approach immediately appear: the degradation pattern showing, which parts in the LMRD are wrong or bad and the root cause analysis.

These both subjects represent quite different area of interest and require specific tools. The LMRD pattern requires more research on the control systems using simulations and industrial data. Moreover, it requires theoretical investigations using probabilistic density functions of distributions common in control engineering. Causality analysis is a separate subject, which uses dedicated tools, like for instance transfer entropy or Granger causality. Both issues are open tasks for further research.

The L-moments and their further variants such as trimmed TL-moments, higher-order LH-moments, partial PL-moments or linear interpolation quantile LQ-moments can help in control engineering. Moment ratio diagrams and homogeneity and discordance analyses may also be helpful. The authors consider this research as a promising opening for further investigations, showing that the extreme statistics research is worth considering in control engineering.