A water quality anomaly identification method and system for a water supply plant

Abstract

The application discloses a water quality anomaly identification method and system for a water supply plant, relates to the technical field of water quality anomaly identification, and comprises the following steps: constructing a VARX model containing endogenous variable lag terms and exogenous variable lag terms, and training the VARX model by taking the minimization of a convex loss function as a target, wherein the convex loss function comprises a sample weight balancing term for balancing the fitting priority of normal and abnormal samples and a process output constraint term conforming to the process practice; analyzing historical sample monitoring data based on a box-and-whisker plot to determine the normal fluctuation interval of various variables; inputting real-time monitoring data into the trained VARX model to obtain water quality index predicted values; and comparing the water quality index predicted values with the normal fluctuation interval to determine whether the water quality is abnormal. The application guarantees that model parameter estimation converges to a unique global optimal solution under the premise of improving anomaly identification precision, and improves the stability of the model.

Description

Technical Field

[0001] This invention relates to the field of water quality anomaly identification technology, specifically to a method and system for identifying water quality anomalies in water supply plants. Background Technology

[0002] Existing water quality anomaly detection methods in water supply plants often utilize ordinary least squares (OLS) to characterize univariate linear relationships and achieve correlation analysis of static water quality indicators.

[0003] Alternatively, traditional VAR models based on the MSE loss function can be used. Traditional VAR models capture the dynamic coupling relationships of water quality parameters (such as pH, turbidity, and ammonia nitrogen) through multivariate time series analysis, characterizing the dynamic linkages of multiple endogenous water quality indicators, but they do not incorporate exogenous influencing factors. The VAR model construction and anomaly identification process mainly includes:

[0004] 1) Data preprocessing: Collect historical water quality data from the water supply plant (such as turbidity of filtered water and coagulant dosage), standardize the data to eliminate the influence of dimensions, remove invalid data caused by equipment failure and measurement error, and center the data (eliminate intercept).

[0005] 2) Variable definition and model structure construction: The endogenous water quality indicators are identified as the research objects, lagged terms of endogenous variables are introduced, and regression equations are constructed;

[0006]

[0007] in: Y t Let be the vector of endogenous variables at time t. The parameter matrix to be estimated, The term represents the residual, and p represents the lag order.

[0008] 3) Lag order selection: Determine the optimal lag order p of the endogenous variables based on the AIC criterion to balance model fit and complexity;

[0009] 4) Model training: With the goal of minimizing the MSE loss function, the model parameters are estimated using ordinary least squares, and the parameter matrix is ​​optimized by minimizing the global prediction bias;

[0010] 5) Anomaly identification: After training, calculate the model's predicted value or directly determine whether the actual value of the water quality index exceeds a certain threshold to identify anomalies.

[0011] Traditional VAR models are based on the MSE loss function, which does not consider the effects of exogenous variables, is detached from process practice, and lacks theoretical support for data processing. Furthermore, this loss function has poor adaptability. If the water supply plant has a large proportion of normal samples and a small number of abnormal samples, the MSE loss function will preferentially fit the normal samples and conservatively estimate the abnormal samples. Moreover, this partially improved loss function is non-convex, and the parameter estimation is prone to getting trapped in local optima, resulting in insufficient model stability. OLS is a single-index static analysis that does not incorporate constraints such as national water quality standards and equipment limits, and the output values ​​may be invalid, resulting in weak model analysis capabilities. Summary of the Invention

[0012] The purpose of this invention is to provide a method and system for identifying water quality anomalies in water supply plants. The method incorporates exogenous variables that affect endogenous indicators into the model and improves the model by using the minimization of a convex loss function as the objective of training the model to balance the weights of normal and abnormal samples. At the same time, it incorporates process constraints to improve the accuracy of anomaly identification and ensure that the model parameter estimation converges to a unique global optimum, thereby improving the stability of the model.

[0013] To achieve the above objectives, this application proposes the following solution:

[0014] On the one hand, the present invention provides a method for identifying water quality anomalies in water supply plants, specifically including the following steps:

[0015] S1. Obtain historical sample monitoring data of the water supply plant process link, cluster the historical sample monitoring data into endogenous variables and exogenous variables, construct a VARX model containing the lagged terms of endogenous variables and the lagged terms of exogenous variables, and use the AIC criterion to determine the optimal lag order of the VARX model.

[0016] S2. Use historical sample monitoring data to train the VARX model with the goal of minimizing the convex loss function to obtain a trained VARX model. The convex loss function includes a sample weight balancing term to balance the fitting priority of normal and abnormal samples and a process output constraint term that conforms to process practice.

[0017] S3. Analyze historical sample monitoring data based on box-and-whisker diagrams to determine the normal fluctuation range of various variables;

[0018] S4. Input the real-time monitoring data into the trained VARX model to obtain the predicted values ​​of water quality indicators. Determine whether the predicted values ​​of water quality indicators are within the normal fluctuation range. If they exceed the range, the water quality is judged to be abnormal.

[0019] In some specific implementation schemes, after clustering historical sample monitoring data into endogenous and exogenous variables, the following steps are also included:

[0020] The endogenous and exogenous variables are centered by first-order difference, respectively, to obtain the first-order difference vectors corresponding to the endogenous and exogenous variables.

[0021] In some specific implementation schemes, the core regression equation of the VARX model is:

[0022] ;

[0023] Where p and q represent the optimal lag order, This represents the first-order difference vector corresponding to the endogenous variables at times t-1 to tp. Let be the first-order difference vector of the exogenous variables at times t-1 to tq; Let p be the parameter matrix of the endogenous variables. Let be the parameter matrix of the exogenous variable lagged by order q. , Let Σ denote the vector of random perturbation terms at time t, and let Σ denote the covariance matrix. N (0,Σ) represents a multivariate normal distribution.

[0024] In some specific implementation schemes, the optimal lag orders p and q are determined by the AIC criterion. The p and q corresponding to the minimum AIC value are selected as the optimal lag orders. The number of variables calculated in the core regression equation is: the number of selected variables × the corresponding lag order.

[0025] In some specific implementations, the convex loss function L is:

[0026] ;

[0027]

[0028]

[0029] in, For sample weight balancing term, These are process output constraints. w t The weights of the endogenous variables at time t are represented; T represents the total duration of the historical sample monitoring data. The first-order difference vector of the VARX model at time t for the endogenous variables. The predicted value, where λ represents the process constraint penalty coefficient. Let be the predicted value of the i-th endogenous variable at time t; This represents the standard process value of the i-th endogenous variable; Let m be the allowable fluctuation range for the i-th endogenous variable, and m represent the number of endogenous variables.

[0030] In some specific implementation schemes, the weights of endogenous variables at time t w t :

[0031]

[0032]

[0033] Where α∈[1,∞] is the weight amplification factor, S max S represents the upper limit of the weight. max >0; s t This is an indicator of the degree of anomaly at time t; Let be the standard value of the first-order difference vector corresponding to the i-th endogenous variable; Let be the allowable fluctuation range of the first-order difference vector corresponding to the i-th endogenous variable.

[0034] In some specific implementations, the training process also includes the following steps:

[0035] Construct the adjoint matrix M of the VARX model:

[0036] M=

[0037] Among them, I m Let m be the m-order identity matrix, where m is the number of endogenous variables;

[0038] By solving the characteristic equation of the adjoint matrix M , I Represent the identity matrix and obtain all eigenvalues ​​Λ;

[0039] Training of the VARX model begins when the absolute value of all eigenvalues ​​|Λ| is less than 1; otherwise, the VARX model is readjusted.

[0040] In some specific implementation plans, the following steps are also included: when the water quality is determined to be abnormal, the cause of the abnormality is analyzed, and the specific process is as follows:

[0041] Perform Cholesky decomposition on the covariance matrix Σ, and extract the random perturbation term vector. Transform into orthogonalized perturbation terms ;

[0042] The VARX model is transformed into an infinite-order moving average form by hysteresis operator transformation;

[0043] By orthogonalizing the perturbation term The error term coefficients in the VARX model of infinite-order moving average form are decoupled to obtain orthogonalized moving average coefficients. The impulse response of endogenous and exogenous variables is calculated based on the orthogonalized moving average coefficients.

[0044] The endogenous or exogenous variable with the largest impulse response value is marked as an abnormal indicator.

[0045] In some specific implementation schemes, the following steps are also included: For the cause of the anomaly, output the adjustment direction; the specific process is as follows:

[0046] The variance decomposition method is used to quantify the long-term stable contribution of each influencing variable to the abnormal indicator. The influencing variables include endogenous and exogenous variables.

[0047] Identify the variables that contribute the most to long-term performance and implement key controls.

[0048] Secondly, this application provides a water quality anomaly identification system for a water supply plant, comprising:

[0049] The data preprocessing module is used to acquire historical sample monitoring data of the water supply plant's process link and cluster the historical sample monitoring data into endogenous variables and exogenous variables.

[0050] The model building module is used to construct a VARX model that includes lag terms of endogenous and exogenous variables, and to determine the optimal lag order of the VARX model using the AIC criterion.

[0051] The model training module is used to train the VARX model with the goal of minimizing the convex loss function using historical sample monitoring data, and obtain the trained VARX model. The convex loss function includes a sample weight balancing term to balance the fitting priority of normal and abnormal samples and a process output constraint term that conforms to process practice.

[0052] The abnormal data identification module is used to analyze historical sample monitoring data based on box-and-whisker plots to determine the normal fluctuation range of various variables;

[0053] The anomaly detection module is used to obtain real-time monitoring data and input it into the trained VARX model to obtain predicted values ​​of water quality indicators. It then determines whether the predicted values ​​of water quality indicators are within the normal fluctuation range. If they exceed the range, the water quality is judged to be abnormal.

[0054] The advantages of this invention over the prior art are as follows:

[0055] This invention, based on the traditional VAR model, introduces exogenous process influence variables, extending the model into a VARX model. This allows the model to better align with process logic. The classification criteria for endogenous / exogenous variables are clearly defined according to the water treatment process logic (e.g., controllable process parameters are classified as endogenous variables, and external environmental driving parameters as exogenous variables). This enables the VARX model to perform regression solutions according to the actual process operation logic, avoiding the blindness of variable classification, improving model interpretability and practical relevance, and ensuring that the model output accurately matches the water treatment process control requirements. Furthermore, the first-order differencing method is used to eliminate the intercept term, ensuring the simplicity, rationality, and interpretability of the VARX model design.

[0056] To adapt to water quality anomaly analysis, a targeted convex loss function (such as introducing anomaly sample weighting coefficients and time series smoothness constraints) was designed to address the temporal characteristics of water quality data and the specificities of anomaly samples. This improves the fitting accuracy of the VARX model to the fluctuation trends of water quality indicators and enhances the model's sensitivity to anomaly data, providing a fundamental guarantee for the accuracy of subsequent anomaly identification and attribution.

[0057] To address the causes of water quality anomalies, impulse response analysis and variance decomposition techniques are deeply integrated with water management operations and maintenance practices. A practical application framework is established that combines "impulse response analysis to locate short-term instantaneous impact sources" and "variance decomposition to lock in long-term core influencing factors." The specific application procedures and judgment criteria of these two types of analytical techniques in emergency response to water anomalies and long-term steady-state control are clarified, thus solving the problems of disconnect and insufficient targeting of traditional time series analysis techniques in water management practices. Attached Figure Description

[0058] Figure 1 A flowchart of a water quality anomaly identification method for water supply plants provided in an embodiment of the present invention;

[0059] Figure 2 This is a schematic diagram of pulse impact between endogenous variables provided in an embodiment of the present invention;

[0060] Figure 3 This is a schematic diagram of pulse impact between exogenous variables provided in an embodiment of the present invention;

[0061] Figure 4 This is a schematic diagram illustrating the contribution percentage of variables to fluctuations, provided in an embodiment of the present invention. Detailed Implementation

[0062] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. The following description of at least one exemplary embodiment is merely illustrative and is in no way intended to limit the present invention or its application or use. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0063] Unless otherwise specifically stated, the relative arrangement, numerical expressions, and values ​​of the components and steps described in these embodiments do not limit the scope of the invention.

[0064] At the same time, it should be understood that, for ease of description, the dimensions of the various parts shown in the accompanying drawings are not drawn according to actual scale.

[0065] Furthermore, for clarity and brevity, descriptions of well-known structures, functions, and configurations may have been omitted. Those skilled in the art will recognize that various changes and modifications can be made to the examples described herein without departing from the spirit and scope of this disclosure.

[0066] Techniques, methods, and equipment known to those skilled in the art may not be discussed in detail, but where appropriate, such techniques, methods, and equipment should be considered part of the specification.

[0067] In all examples shown and discussed herein, any specific values ​​should be interpreted as merely exemplary and not as limitations. Therefore, other examples of exemplary embodiments may have different values.

[0068] Definitions of terms used in this application:

[0069] Cholesky decomposition: A matrix decomposition method that decomposes the covariance matrix into the product of a lower triangular matrix and its transpose (as described below). This can achieve orthogonalization of the disturbance term, eliminate the synchronous correlation between endogenous variables, and provide a basis for the accurate identification of the impact effect in impulse response analysis.

[0070] Mean Squared Error (MSE) Loss Function L MSE Loss metric is used to measure the deviation between model predictions and actual values. It is calculated as the mean of the squared differences between the predicted and actual values ​​and is a commonly used loss metric in model training (such as VAR / VARX model parameter estimation).

[0071]

[0072] in: These are the model's predicted values. The values ​​are the true values. Parameter estimation is achieved by minimizing the loss function, and anomaly identification is finally completed based on the predicted residuals or the actual threshold values ​​of the indicators.

[0073] VAR: Vector Autoregression, a multivariate time series model. Its core is to construct regression relationships through the lagged terms of multiple endogenous variables, characterize the dynamic linkage of multiple endogenous water quality indicators, and capture the dynamic interaction between variables. It is suitable for analyzing the time series fluctuations of multivariate systems (such as the dynamic correlation between water quality indicators), but does not include exogenous influencing factors.

[0074] VARX: Vector Autoregressive Exogenous Model, an extended form of the VAR model. It adds exogenous variables to the endogenous variables, and can simultaneously characterize the dynamic relationship between endogenous variables and the driving effect of exogenous variables (such as the water quality analysis model that incorporates exogenous variables such as influent flow rate in the document).

[0075] AIC: Akaike Information Criterion, a model selection criterion that determines the optimal model by balancing the goodness of fit and complexity (number of parameters). It is often used to select the lag order of VAR / VARX models (such as determining p and q values) to avoid overfitting.

[0076] OLS: Ordinary Least Squares, a parameter estimation method that solves for the parameters of a linear model by minimizing the sum of squares of the model residuals, characterizing univariate linear relationships, and enabling correlation analysis of static water quality indicators; it is the basic method for parameter estimation of linear models such as VAR and VARX.

[0077] Example 1

[0078] like Figure 1 As shown in the figure, this embodiment provides a method for identifying water quality anomalies in a water supply plant, specifically including the following steps:

[0079] S1. Obtain historical sample monitoring data of the water supply plant process link, cluster the historical sample monitoring data into endogenous variables and exogenous variables, construct a VARX model containing the lagged terms of endogenous variables and the lagged terms of exogenous variables, and use the AIC criterion to determine the optimal lag order of the VARX model.

[0080] Endogenous variables include PAC dosage, pre-chlorination dosage, and pre-filtration turbidity, which are used to characterize water quality indicators; exogenous variables include process operation and environmental variables that affect endogenous variables, including influent turbidity and influent flow rate.

[0081] After clustering historical sample monitoring data into endogenous and exogenous variables, the following steps are also included:

[0082] The endogenous and exogenous variables are centered by first-order difference, respectively, to obtain the first-order difference vectors corresponding to the endogenous and exogenous variables.

[0083] The optimal lag orders p and q are determined by the AIC criterion. The p and q corresponding to the smallest AIC value are selected as the optimal lag orders. The number of variables calculated in the core regression equation is: the number of selected variables × the corresponding lag order.

[0084] S2. Use historical sample monitoring data to train the VARX model with the goal of minimizing the convex loss function to obtain a trained VARX model. The convex loss function includes a sample weight balancing term to balance the fitting priority of normal and abnormal samples and a process output constraint term that conforms to process practice.

[0085] The core regression equation of the VARX model is:

[0086]

[0087] Where p and q represent the optimal lag order, This represents the first-order difference vector corresponding to the endogenous variables at times t-1 to tp. Let be the first-order difference vector of the exogenous variables at times t-1 to tq; Let p be the parameter matrix of the endogenous variables. Let be the parameter matrix of the exogenous variable lagged by order q. , Let Σ denote the vector of random perturbation terms at time t, and let Σ denote the covariance matrix. N (0,Σ) represents a multivariate normal distribution.

[0088] The convex loss function L is:

[0089] ;

[0090]

[0091]

[0092] in, For sample weight balancing term, These are process output constraints. w t The weights of the endogenous variables at time t are represented; T represents the total duration of the historical sample monitoring data. The first-order difference vector of the VARX model at time t for the endogenous variables. The predicted value, where λ represents the process constraint penalty coefficient. Let be the predicted value of the i-th endogenous variable at time t; This represents the standard process value of the i-th endogenous variable; Let m be the allowable fluctuation range for the i-th endogenous variable, and m represent the number of endogenous variables.

[0093] Weights of endogenous variables at time t w t :

[0094]

[0095]

[0096] Where α∈[1,∞] is the weight amplification factor, S max S represents the upper limit of the weight. max >0; s t This is an indicator of the degree of anomaly at time t; Let be the standard value of the first-order difference vector corresponding to the i-th endogenous variable; Let be the allowable fluctuation range of the first-order difference vector corresponding to the i-th endogenous variable.

[0097] The training process also includes the following steps:

[0098] Construct the adjoint matrix M of the VARX model:

[0099] M=

[0100] Among them, I m Let m be the m-order identity matrix, where m is the number of endogenous variables;

[0101] By solving the characteristic equation of the adjoint matrix M , I Represent the identity matrix and obtain all eigenvalues ​​Λ;

[0102] Training of the VARX model begins when the absolute value of all eigenvalues ​​Λ, |Λ|, is less than 1; otherwise, the VARX model is readjusted.

[0103] S3. Analyze historical sample monitoring data based on box-and-whisker diagrams to determine the normal fluctuation range of various variables;

[0104] S4. Input the real-time monitoring data into the trained VARX model to obtain the predicted values ​​of water quality indicators. Determine whether the predicted values ​​of water quality indicators are within the normal fluctuation range. If they exceed the range, the water quality is judged to be abnormal.

[0105] When water quality is determined to be abnormal, the cause of the abnormality is analyzed. The specific process is as follows:

[0106] Perform Cholesky decomposition on the covariance matrix Σ, and extract the random perturbation term vector. Transform into orthogonalized perturbation terms ;

[0107] The VARX model is transformed into an infinite-order moving average form by hysteresis operator transformation;

[0108] By orthogonalizing the perturbation term The error term coefficients in the VARX model of infinite-order moving average form are decoupled to obtain orthogonalized moving average coefficients. The impulse response of endogenous and exogenous variables is calculated based on the orthogonalized moving average coefficients.

[0109] The endogenous or exogenous variable with the largest impulse response value is marked as an abnormal indicator.

[0110] Based on the cause of the anomaly, the output adjustment direction is determined, and the specific process is as follows:

[0111] The variance decomposition method is used to quantify the long-term stable contribution of each influencing variable to the abnormal indicator. The influencing variables include endogenous and exogenous variables.

[0112] Identify the variables that contribute the most to long-term performance and implement key controls.

[0113] It is understandable that, in response to the current technical limitations of existing OLS and traditional VAR models based on the MSE loss function in water quality anomaly identification in water supply plants, such as "poor adaptability of loss function, detachment from process practice, insufficient model stability, and weak system analysis capability," this paper proposes a method that applies the VARX analysis method to the field of water quality anomaly identification in water supply plants. This method overcomes the limitations of traditional OLS single-index static analysis and traditional VAR neglecting exogenous process variables. It constructs a model based on the multi-parameter coupled process logic of water supply plants, aiming to solve problems such as insufficient fitting of anomaly samples, output values ​​detached from process constraints, unstable parameter estimation, and disconnect between system analysis and process theory. This embodiment employs first-order difference preprocessing to eliminate the intercept term in the regression equation. This not only ensures the rationality of the model setting but also has practical technological significance—effectively eliminating data trend terms and promoting the analysis system to reach a stable state (verified by the absolute value of the eigenvalues ​​of the adjoint matrix being less than 1), laying the foundation for subsequent accurate modeling. A customized convex loss function is designed, incorporating sample balancing into the loss function: by designing differentiated weights to balance the fitting priority of normal samples and scarce abnormal samples, while incorporating process constraints, the model parameter estimation converges to a unique global optimum solution by leveraging the characteristics of the convex loss function, while improving the accuracy of anomaly identification. This solves the problem of the traditional MSE loss function's conservative estimation of abnormal samples and insufficient stability. Based on the VARX model combined with Cholesky decomposition, impulse response analysis and variance decomposition are carried out to accurately locate the key impulse influencing factors and their long-term stable contribution ratios affecting the target water quality indicators, enabling rapid discovery, in-depth analysis, and precise control of anomaly sources, thus compensating for the weakness of traditional model systems in terms of analytical capabilities.

[0114] The method provided in this embodiment is widely used in water quality anomaly identification and analysis scenarios, such as water quality anomaly analysis scenarios in water supply plants. Its potential application areas include:

[0115] 1) Large-scale integrated water supply plants: Suitable for large-scale water supply plants along rivers and coasts where the quality of incoming water fluctuates frequently (with complex exogenous variables). It can accurately identify water quality anomalies caused by exogenous factors such as sudden increases in the concentration of pollutants in the incoming water and tidal effects, and provide support for the collaborative operation and maintenance of multiple process units.

[0116] 2) Small-scale rural water supply stations: Adapted to the characteristics of small-scale water supply stations with simple processes and limited operation and maintenance resources, relying on the stability of convex loss functions and the adaptability of process constraints, it can output anomaly handling suggestions that can be directly implemented without complicated debugging;

[0117] 3) Cross-regional water affairs operation and maintenance platform: Supports the water affairs group to centrally manage and control multiple water plants under its jurisdiction. Based on the rationality and interpretability of the model settings, it realizes standardized analysis and source tracing of abnormal data from different water plants, and improves cross-regional operation and maintenance efficiency.

[0118] To better illustrate the solution concept of this embodiment, the following uses a water supply plant's supply chain as an example to explain the overall solution from data preprocessing, variable setting, model building and training, to cause identification and attribution. Specifically:

[0119] Step 1: Data Preprocessing

[0120] Collect historical sample monitoring data of the water supply plant's process chain, and remove invalid data caused by equipment failure and measurement errors. Taking the coagulation process section as an example, historical sample monitoring data can be collected, such as: influent flow rate Q_in, influent turbidity T_in, PAC dosage PAC_in, pre-chlorination dosage Cl_in, and pre-filtration turbidity T_out. The data granularity can be sampled once every 20 minutes, and the sampling time range is set to a continuous year (covering different seasonal operating conditions and peak water consumption periods to ensure data representativeness). The total data volume is 525,600 sets.

[0121] Based on the interconnectedness of the coagulation process units in the water supply plant, historical sample monitoring data are categorized according to variable type:

[0122] ① Endogenous variables: core water quality indicators that need to be monitored (such as PAC dosage, pre-chlorination dosage, and turbidity of pre-filtration water).

[0123] ② Exogenous variables: These are process operation and environmental variables that affect endogenous indicators, such as influent turbidity and influent flow rate. These variables are explicitly not controlled by the current process unit. Specifically, this can be understood as follows: the operation of the flocculation tank unit will not affect the value of exogenous variables such as influent turbidity. The setting of variables is related to the process logic of the water plant, ensuring that the variable settings closely match the actual process logic and avoiding interference from variables unrelated to the process to improve model accuracy. Simple VAR models lack the concept of exogenous variables, treating all variables as endogenous. This can produce regression equations that do not conform to process logic. For example, influent turbidity being affected by the historical lag term of PAC dosage is clearly unreasonable.

[0124] To ensure the regression equations do not require an intercept term under reasonable process logic, thereby promoting a stable analytical system and laying the foundation for subsequent model specification and stability verification, first-order differencing is used to center all the variables (endogenous and exogenous variables), yielding first-order differencing vectors for the endogenous and exogenous variables, respectively.

[0125] Step 2: VARX model construction, loss function setting, and parameter training

[0126] Model Construction: Based on the first-order difference vectors preprocessed in step 1, a VARX(p,q) model containing lagged terms of endogenous and exogenous variables is constructed, where (p,q) represents the optimal lag order. Its core regression equation is:

[0127] (1)

[0128] in: This represents the first-order difference vector corresponding to the endogenous variables at times t-1 to tp. Let be the first-order difference vector of the endogenous variables at time t, such as , Let be the first-order difference vector corresponding to the turbidity of the water before filtration. Let be the first-order difference vector corresponding to the previous chlorination amount. This is the first-order difference vector corresponding to the PAC dosage;

[0129] Let be the first-order difference vector of the exogenous variables at times t-1 to tq. Let be the first-order difference vector of the exogenous variables at time t, such as , This is the first-order difference vector corresponding to the influent flow rate. This is the first-order difference vector corresponding to the turbidity of the influent;

[0130] Let p be the parameter matrix of the endogenous variables. The expression is as follows:

[0131]

[0132] Elements in the parameter matrix This represents the response coefficient of the j-th endogenous variable lagged p times to the i-th endogenous variable (e.g., (This represents the response coefficient of the m-th endogenous variable after lag p times to the m-th endogenous variable), where m is the number of endogenous variables.

[0133] The parameter matrix is ​​the q-order lagged exogenous variable; The expression is as follows:

[0134]

[0135] In the parameter matrix, This represents the response coefficient of the j-th exogenous variable lagged q times to the i-th endogenous variable (e.g., (This represents the response coefficient of the m-th exogenous variable lagged by q to the n-th endogenous variable), where n is the number of exogenous variables.

[0136] is a vector of random perturbation terms (mean is 0, covariance matrix is ​​Σ). N (0,Σ) represents a multivariate normal distribution; the covariance matrix Σ is expressed as follows:

[0137]

[0138] In the covariance matrix (For example Let be the variance of the disturbance term of the i-th endogenous variable, and the off-diagonal elements are... Let be the covariance of the disturbance terms of the i-th and j-th endogenous variables. Since the variables will not actually be completely collinear, the covariance matrix can usually be considered positive definite, which provides the basis for the subsequent Cholesky decomposition.

[0139] The optimal lag order (p,q) is determined by the AIC criterion, the formula for which is: ( k This represents the total number of parameters to be estimated in the model. L For the model likelihood function value, select (p,q) corresponding to the minimum AIC value as the optimal lag order. It should be noted that, due to the limited actual sample size, the number of variables included in the regression equation is equivalent to "number of selected variables × corresponding lag order", which will additionally occupy the model's degrees of freedom. If the lag order p and q are too high or the number of selected variables is too large, it is easy to cause the model to overfit and the estimation accuracy to decrease. Therefore, the lag order p and q and the number of selected variables should not be too large.

[0140] ① Loss Function Setting: Construct a convex loss function that includes a "sample weight balancing term" and a "process output constraint term". The formula for the total loss function L is:

[0141]

[0142] in:

[0143] 1) Sample weight balancing term T represents the total duration of historical sample monitoring data, used to balance the fitting priority between normal and abnormal samples; the variable weights at time t. α∈[1,∞] is the weight amplification factor; S max S represents the upper limit of the weight. max >0, to avoid excessive interference from extreme samples (can be customized according to process requirements); This is an indicator of the degree of anomaly at time t; The standard value of the first-order difference variable corresponding to the i-th endogenous variable (can be customized according to process requirements); Its allowable fluctuation range (can be customized according to process requirements); For the model pair The predicted value, It is the square of the L2 norm.

[0144] 2) Process output constraints Ensure that the model output conforms to the actual process; This is the penalty coefficient for process constraints (which can be customized according to process requirements). The predicted value of the i-th endogenous variable can be obtained through the differenced predicted value of the i-th endogenous variable. After inverse difference operation get; The process annotation value for the i-th endogenous variable (can be customized according to process requirements); Its allowable fluctuation range.

[0145] 3) The loss function is a convex function, as proven below:

[0146] Based on publicly available information, the following convex function form is known:

[0147] absolute value function

[0148]

[0149] Note: It concerns the parameter matrix to be estimated. The affine function, since the composition of a convex function and an affine function is still a convex function, is further defined as follows: ,therefore, It is a convex function. The essence is also about the parameter matrix. Since the composition of a convex function and an affine function is still a convex function, and the composite function of convex functions is also a convex function, therefore... It is a convex function. Since the non-negative linear combination of convex functions is still a strictly convex function, the total loss function L is a convex function.

[0150] ② Parameter Training: Before parameter training, it is necessary to ensure the stationarity of the VARX model; otherwise, the regression interpretation will be meaningless or even erroneous. It is necessary to check whether the absolute values ​​of all eigenvalues ​​of the adjoint matrix are less than 1. The formula for constructing the adjoint matrix is:

[0151]

[0152] in:

[0153] I m Given an m-order identity matrix (where m is the number of endogenous variables); solve the characteristic equation. If the absolute value |Λ| is less than 1, and this holds true for all eigenvalues, then the model is stable; otherwise, the lag order or variable settings need to be readjusted.

[0154] After passing the model stability test, the parameters to be estimated in the model are solved with the objective of minimizing the convex loss function designed above. The optimization objective formula is:

[0155]

[0156] The gradient descent algorithm can be used to solve the above optimization problem. Combined with the characteristic of convex functions that "there is no local optimum and only a unique global optimum", it can be guaranteed that the parameter estimation converges to the unique global optimum. The iteration termination condition is that the parameter update step size is less than a preset threshold or the change in the loss function value is less than a preset threshold.

[0157] Step 3: Water Quality Anomaly Identification

[0158] ① Historical anomaly identification: Based on the preprocessed historical data from box-and-whisker plot analysis, the normal fluctuation range (interquartile range) of each water quality indicator is determined, and historical anomaly samples are identified;

[0159] ② Future anomaly prediction: After the real-time monitoring data is preprocessed in step 1, it is input into the trained VARX model to obtain the predicted values ​​of water quality indicators; the normal range determined by the box plot is input to determine whether the future data exceeds the normal fluctuation range. If it does, an anomaly warning is triggered.

[0160] The anomaly detection results are shown in Table 1 below:

[0161] Table 1 Evaluation Table of Water Quality Anomaly Identification Effectiveness

[0162]

[0163] Note: The data in Table 1 comes from 9120 data points of particle size over a 20-minute period at a water supply plant, and all are predictions for the future. Identification accuracy = (Number of identified abnormal samples / Actual number of abnormal samples) × 100%; False alarm rate = (Number of false alarm samples / Predicted number of samples) × 100%, where false alarm samples refer to normal samples that are incorrectly identified as abnormal.

[0164] Step 4: Attribution of Water Quality Anomalies and Recommendations

[0165] 1) Anomaly Attribution - Impulse Response: By analyzing impulse analysis coefficients, endogenous / exogenous variables that have a significant impact on the target variable are analyzed to achieve anomaly attribution analysis. The steps are as follows:

[0166] a) Perform Cholesky decomposition on the covariance matrix Σ, such that: The original disturbance term Transform into orthogonalized perturbation terms The relationship between the two is as follows:

[0167] or (2)

[0168] At this point, the orthogonalization perturbation term satisfy I is the identity matrix, meaning that the structural impact terms are independent of each other and there is no synchronous correlation, thus achieving accurate identification of the impact effect.

[0169] b) Transform the VARX(p,q) model into an infinite-order moving average MA(∞) form, and rewrite formula (1) using the lag operator transformation:

[0170]

[0171] in:

[0172] For endogenous variable lag operator polynomials;

[0173] For exogenous variable lag operator polynomials;

[0174] L For lag operators (satisfying) ).

[0175] c) To Inverse equation of the model yields the MA(∞) expansion:

[0176] (3)

[0177] For endogenous variables, since there may be contemporaneous correlations between error terms (for example, the PAC dosage and the previous chlorination dosage may be subjected to the same shock in the same period, or the machine may suddenly malfunction and the dosage may increase, causing the error terms of the two to be correlated in the same period), it is not possible to analyze the effect of the change of a single error term by keeping other error terms constant when interpreting the model, because when this error term changes, other error terms may change at the same time; this requires decoupling the error term coefficients in formula (3) based on Cholesky decomposition (formula (2)). Substituting formula (2) into formula (3), we get:

[0178] (4)

[0179] in:

[0180] The orthogonalized moving average coefficient matrix. for m×m The coefficient matrix. For exogenous variables, assuming they are independent, the coefficient matrix in formula (4) can be used. Calculate the impact of the shock (i.e., the traditional analysis method, which keeps other variables constant, and increases the target variable by a coefficient when the target variable increases by 1 unit).

[0181] Based on the above expansion, the impulse response function Defined as: when the j-th structural impact term When a positive shock of one unit occurs, while other shock terms and exogenous variables remain unchanged, the i-th endogenous variable... In the t+s Expect( s The response value is the result of the impact hysteresis order. Its mathematical expression is:

[0182]

[0183] in:

[0184] Coefficient matrix The element in the i-th row and j-th column.

[0185] d) Examples of impulse shocks between endogenous / exogenous variables, such as... Figure 2 , Figure 3 As shown: Among them, Figures 2-3 The horizontal axis represents the response period (lag period), and the vertical axis represents the impulse response value. Figure 2 This represents the impulse response of the lag period of PAC dosage to the turbidity of the pre-filtration water, which is an endogenous variable. Figure 3This represents the impulse response of the lag period of the influent turbidity in the exogenous variable to the PAC dosage. Figure 2 The point (1, -0.128) means that when the change in PAC dosage with a lag of 1 period is impacted by 1 standard deviation, the change in turbidity of the pre-filtration water in the current period will change by -0.128 units. When attributing water quality anomalies, endogenous / exogenous variables with the largest impulse response values ​​can be prioritized for anomaly investigation.

[0186] 2) Anomaly Attribution - Stability Control: By quantifying the long-term stable contribution of each influencing variable to the anomaly indicator through variance decomposition, targeted operation and maintenance adjustment suggestions are output to support precise anomaly control. The specific steps are as follows:

[0187] a) Form of prediction error: Combining formula (4), the i-th target endogenous variable s-period prediction error It can be represented in the form of a hysteresis operator matrix, that is:

[0188]

[0189] in:

[0190] The contribution of endogenous shocks to the prediction error of the i-th target variable in period s. For the orthogonalized moving average operator matrix, The orthogonalized impact vector is defined by the superscript s, which indicates that the lag contribution from the previous s periods is truncated.

[0191] The contribution of exogenous shocks to the prediction error of the i-th target variable in period s. For the exogenous shock propagation operator matrix, This is the difference vector of exogenous variables, and the superscript s has the same meaning as above.

[0192] b) Total variance form: due to endogenous shocks (orthogonalized shock vector) Exogenous shocks They are mutually independent (VARX model assumption, see the derivation process of formula (4)), the total variance of the prediction error Represented as:

[0193]

[0194] c) Form of the proportion of variance contribution of endogenous / exogenous variables: Based on the matrix expression of the total variance, the proportion of contribution of endogenous and exogenous variables can be uniformly expressed as the ratio of the variance of the corresponding shock contribution to the total variance:

[0195] The contribution percentage of the j-th endogenous variable:

[0196]

[0197] in: The submatrix formed by the column vectors corresponding to the j-th endogenous shock;

[0198] The contribution percentage of the r-th exogenous variable:

[0199]

[0200] in:

[0201] It is a submatrix formed by the column vectors corresponding to the r-th exogenous variable in the exogenous transitive operator matrix.

[0202] when When, operator matrix As the system stabilizes, the contribution ratio of each variable also stabilizes, meaning that the contribution ratio of each variable to the system stabilizes. The long-term stable contribution of fluctuations is the core basis for formulating control recommendations. An example of variance decomposition is... Figure 4 As shown, the variance decomposition of the total effluent turbidity from the flocculation tank, including exogenous variables, is presented, along with the contribution of each variable to the fluctuation. Figure 4 The horizontal axis represents the future period 's', and the vertical axis represents the percentage of variance contribution. It can be seen that as 's' increases, the variance of the total turbidity in the flocculation tank effluent is mainly contributed by the PAC dosage.

[0203] 3) Process Integration and Operation & Maintenance Recommendations: Variance decomposition and impulse response analysis are not used in isolation; they form a complementary and synergistic analytical system. Variance decomposition focuses on the long-term stable contribution weight of variables to the target indicator, clarifying the priority of impact; impulse response analysis focuses on the short-term dynamic transmission process of the impact, capturing the intensity and duration of the instantaneous impact. If the core operation and maintenance objective is to ensure the long-term stable operation of the water quality system (e.g., controlling the long-term fluctuation of turbidity in pre-filtration water within the acceptable range), then the variance decomposition results should be used as the core basis to prioritize the identification of variables with the highest long-term contribution and implement key control measures. If it is necessary to deal with sudden water quality anomalies (e.g., a sudden surge in turbidity), then the impulse response diagram should be relied upon to quickly identify the key variables that cause the instantaneous impact and formulate an emergency response plan.

[0204] Example 2

[0205] This embodiment provides a water quality anomaly identification system for a water supply plant, used to implement the water quality anomaly identification method proposed in Embodiment 1, including:

[0206] The data preprocessing module is used to acquire historical sample monitoring data of the water supply plant's process link and cluster the historical sample monitoring data into endogenous variables and exogenous variables.

[0207] The model building module is used to construct a VARX model that includes lag terms of endogenous and exogenous variables, and to determine the optimal lag order of the VARX model using the AIC criterion.

[0208] The model training module is used to train the VARX model with the goal of minimizing the convex loss function using historical sample monitoring data, and obtain the trained VARX model. The convex loss function includes a sample weight balancing term to balance the fitting priority of normal and abnormal samples and a process output constraint term that conforms to process practice.

[0209] The abnormal data identification module is used to analyze historical sample monitoring data based on box-and-whisker plots to determine the normal fluctuation range of various variables;

[0210] The anomaly detection module is used to obtain real-time monitoring data and input it into the trained VARX model to obtain predicted values ​​of water quality indicators. It then determines whether the predicted values ​​of water quality indicators are within the normal fluctuation range. If they exceed the range, the water quality is judged to be abnormal.

[0211] The above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention in any way. Based on the technical essence of the present invention, any simple modifications, equivalent substitutions, and improvements made to the above embodiments within the spirit and principles of the present invention shall still fall within the protection scope of the present invention.

Claims

1. A method for identifying water quality anomalies in a water supply plant, characterized in that, Specifically, the following steps are included: S1. Obtain historical sample monitoring data of the water supply plant process link, cluster the historical sample monitoring data into endogenous variables and exogenous variables, construct a VARX model containing the lagged terms of endogenous variables and the lagged terms of exogenous variables, and use the AIC criterion to determine the optimal lag order of the VARX model. Endogenous variables are used to characterize water quality indicators; exogenous variables include process operations and environmental variables that affect endogenous variables. S2. Use historical sample monitoring data to train the VARX model with the goal of minimizing the convex loss function to obtain a trained VARX model. The convex loss function includes a sample weight balancing term to balance the fitting priority of normal and abnormal samples and a process output constraint term that conforms to process practice. The convex loss function L is: ； in, For sample weight balancing term, Indicates process output constraints; w t The weights of the endogenous variables at time t are represented; T represents the total duration of the historical sample monitoring data. The first-order difference vector of the VARX model at time t for the endogenous variables. The predicted value, where λ represents the process constraint penalty coefficient. Let be the predicted value of the i-th endogenous variable at time t; This represents the standard process value of the i-th endogenous variable; Let m be the allowable fluctuation range for the i-th endogenous variable, and m represent the number of endogenous variables. S3. Analyze historical sample monitoring data based on box-and-whisker diagrams to determine the normal fluctuation range of various variables; S4. Input the real-time monitoring data into the trained VARX model to obtain the predicted values ​​of water quality indicators. Determine whether the predicted values ​​of water quality indicators are within the normal fluctuation range. If they exceed the range, the water quality is judged to be abnormal.

2. The method for identifying abnormal water quality in a water supply plant according to claim 1, characterized in that, After clustering historical sample monitoring data into endogenous and exogenous variables, the following steps are also included: The endogenous and exogenous variables are centered by first-order difference, respectively, to obtain the first-order difference vectors corresponding to the endogenous and exogenous variables.

3. The method for identifying water quality anomalies in a water supply plant according to claim 2, characterized in that, The core regression equation of the VARX model is: ； Where p and q represent the optimal lag order, This represents the first-order difference vector corresponding to the endogenous variables at times t-1 to tp. Let be the first-order difference vector of the exogenous variables at times t-1 to tq; Let p be the parameter matrix of the endogenous variables. Let be the parameter matrix of the exogenous variable lagged by order q. , Let Σ denote the vector of random perturbation terms at time t, and let Σ denote the covariance matrix. N (0,Σ) represents a multivariate normal distribution.

4. The method for identifying water quality anomalies in a water supply plant according to claim 3, characterized in that, The optimal lag orders p and q are determined by the AIC criterion. The p and q corresponding to the smallest AIC value are selected as the optimal lag orders. The number of variables calculated in the core regression equation is: the number of selected variables × the corresponding lag order.

5. The method for identifying water quality anomalies in a water supply plant according to claim 1, characterized in that, Weights of endogenous variables at time t w t : Where α∈[1,∞] is the weight amplification factor, S max S represents the upper limit of the weight. max >0; s t This is an indicator of the degree of anomaly at time t; Let be the standard value of the first-order difference vector corresponding to the i-th endogenous variable; Let be the allowable fluctuation range of the first-order difference vector corresponding to the i-th endogenous variable.

6. The method for identifying water quality anomalies in a water supply plant according to claim 3, characterized in that, The training process also includes the following steps: Construct the adjoint matrix M of the VARX model: M= Among them, I m Let m be the m-order identity matrix, where m is the number of endogenous variables; By solving the characteristic equation of the adjoint matrix M , I Represent the identity matrix and obtain all eigenvalues ​​Λ; Training of the VARX model begins when the absolute value of all eigenvalues ​​|Λ| is less than 1; otherwise, the VARX model is readjusted.

7. The method for identifying water quality anomalies in a water supply plant according to claim 3, characterized in that, It also includes the following steps: When the water quality is determined to be abnormal, analyze the cause of the abnormality. The specific process is as follows: Perform Cholesky decomposition on the covariance matrix Σ, and extract the random perturbation term vector. Transform into orthogonalized perturbation terms ; The VARX model is transformed into an infinite-order moving average form by hysteresis operator transformation; By orthogonalizing the perturbation term The error term coefficients in the VARX model of infinite-order moving average form are decoupled to obtain orthogonalized moving average coefficients. The impulse response of endogenous and exogenous variables is calculated based on the orthogonalized moving average coefficients. The endogenous or exogenous variable with the largest impulse response value is marked as an abnormal indicator.

8. The method for identifying abnormal water quality in a water supply plant according to claim 7, characterized in that, It also includes the following steps: Based on the cause of the anomaly, adjust the output direction; the specific process is as follows: The variance decomposition method is used to quantify the long-term stable contribution of each influencing variable to the abnormal indicator. The influencing variables include endogenous and exogenous variables. Identify the variables that contribute the most to long-term performance and implement key controls.

9. A water quality anomaly identification system for a water supply plant, characterized in that, include: The data preprocessing module is used to acquire historical sample monitoring data of the water supply plant's process link and cluster the historical sample monitoring data into endogenous variables and exogenous variables. Endogenous variables are used to characterize water quality indicators; exogenous variables include process operations and environmental variables that affect endogenous variables. The model building module is used to construct a VARX model that includes lag terms of endogenous and exogenous variables, and to determine the optimal lag order of the VARX model using the AIC criterion. The model training module is used to train the VARX model with the goal of minimizing the convex loss function using historical sample monitoring data, and obtain the trained VARX model. The convex loss function includes a sample weight balancing term to balance the fitting priority of normal and abnormal samples and a process output constraint term that conforms to process practice. The convex loss function L is: ； in, For sample weight balancing term, Indicates process output constraints; w t The weights of the endogenous variables at time t are represented; T represents the total duration of the historical sample monitoring data. The first-order difference vector of the VARX model at time t for the endogenous variables. The predicted value, where λ represents the process constraint penalty coefficient. Let be the predicted value of the i-th endogenous variable at time t; This represents the standard process value of the i-th endogenous variable; Let m be the allowable fluctuation range for the i-th endogenous variable, and m represent the number of endogenous variables. The abnormal data identification module is used to analyze historical sample monitoring data based on box-and-whisker plots to determine the normal fluctuation range of various variables; The anomaly detection module is used to obtain real-time monitoring data and input it into the trained VARX model to obtain predicted values ​​of water quality indicators. It then determines whether the predicted values ​​of water quality indicators are within the normal fluctuation range. If they exceed the range, the water quality is judged to be abnormal.

Images