A method for judging accuracy of water pollution monitoring data based on multivariate statistics

Abstract

The application discloses a kind of water pollution monitoring data accuracy determination method based on multivariate statistics.The method starts from the source of monitoring online monitoring instrument operating principle, real-time acquisition monitoring data and source parameters such as measurement original value, calibration curve slope, calibration curve intercept;The monitoring data and source parameters at the same time are organized as multidimensional time sequence feature vector;Using multivariate statistics to establish the regression model of monitoring data, the consistency matching of fitting parameter and online monitoring instrument reported parameter is carried out, to determine the source end calculation rule actually followed by online monitoring instrument;Real-time determination, the measurement original value is substituted into source end calculation rule to calculate theoretical concentration, and the output monitoring data of online monitoring instrument is dynamically compared, and the suspicion of artificial tampering is identified according to the deviation degree.The application verifies the physical logic consistency between the data in online monitoring instrument and source parameters, and improves the accuracy and reliability of tampering determination from the root of data generation.

Description

Technical Field

[0001] This invention relates to the field of water pollution monitoring data accuracy assessment technology, and more specifically, to a method for assessing the accuracy of water pollution monitoring data based on multivariate statistics. Background Technology

[0002] The authenticity, accuracy, and completeness of online pollution source monitoring data are the lifeline and cornerstone of environmental management. Their importance permeates all levels, including government decision-making, environmental law enforcement, corporate operations, market fairness, and public health. Furthermore, online monitoring data serves as the direct legal basis for determining whether a company is complying with pollution regulations. High-quality data helps law enforcement agencies accurately target illegal companies under "non-site supervision" models while protecting the legitimate rights and interests of compliant companies. Conversely, low-quality data can lead to disputes and even the awkward situation where "companies violate the law, and environmental protection departments are held responsible."

[0003] In response to the phenomenon of human tampering with monitoring data, early identification primarily relied on manual experience analysis. Currently, the industry is gradually shifting towards intelligent identification systems, with the main technological approaches including: Data trend and process correlation analysis: By analyzing the matching between pollutant concentration and flow rate, production conditions, and purification facility operating parameters, scenarios such as sudden drops in concentration but unchanged flow rate, and data mutations under abnormal operating conditions are identified; Data regularity detection: By identifying whether the data is too "perfect" (such as being stable near a certain low value for a long time), it is inferred whether there is a manipulated "safe value" for fraud; Multi-indicator cross-validation: By analyzing the inherent logical relationship between various indicators (such as the correlation between COD and BOD, ammonia nitrogen and total nitrogen), situations where single indicators are abnormal or the overall data is inconsistent are identified; Time series anomaly detection: By comparing the fluctuation patterns of historical data in the same period, data mutations during non-working periods such as early morning are identified.

[0004] However, the above methods share a common limitation: they are essentially all about performing correlation analysis and pattern mining on "result data" (i.e., the generated pollutant concentration values) from different dimensions. They rely on specific field scenario assumptions, involve a wide variety of data, require a large amount of computation, and the accuracy of the identification results is often questionable—especially in distinguishing between real extreme operating conditions and human data tampering. More importantly, these methods fail to address the "source" of the data, that is, they cannot verify whether the monitoring instruments that generate the data are in a real, tamper-free operating state.

[0005] Currently, the monitoring of water pollution sources such as COD, ammonia nitrogen, total phosphorus, and total nitrogen mainly employs the principle of spectrophotometry. Spectrophotometry is a method for qualitative and quantitative analysis of substances by measuring their absorbance or luminescence intensity at specific wavelengths or within a certain wavelength range. Its quantitative analysis is based on the Lambert-Beer Law: when a beam of monochromatic light passes through a solution containing an absorbing substance, the absorbance of the solution is directly proportional to the concentration of the absorbing substance and the thickness of the absorption layer. A commonly used quantitative analysis method is the standard curve method, which involves preparing a series of standard solutions of known concentrations, measuring their absorbance to plot a standard curve, and then calculating the sample concentration based on the sample's absorbance using this standard curve.

[0006] Based on the correlation between pollutant concentration and key parameters inside the monitoring instrument, if the monitoring data is maliciously modified, it will directly cause the correlation between the data to fail, and such tampering is extremely difficult to detect in routine data audits. Summary of the Invention

[0007] The purpose of this invention is to provide a method for determining the accuracy of water pollution monitoring data based on multivariate statistics, in order to solve the above-mentioned problems existing in the prior art.

[0008] The application is as follows: A method for determining the accuracy of water pollution monitoring data based on multivariate statistics, the method comprising: S1. Real-time acquisition of result data and source parameters output by the online water pollution source monitoring instrument in measurement mode; the result data is the final monitoring data output by the online monitoring instrument; the source parameters include the original measurement value, the slope of the calibration curve, and the intercept of the calibration curve; S2. After data cleaning and denoising, the monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept at the same measurement time are organized according to the time window to construct a multi-dimensional time series feature vector. S3. Based on multivariate statistical analysis, a regression model is established using the multidimensional time-series feature vector to connect monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept. The regression model matches the fitted parameters with the source parameters reported by the online monitoring instrument. Based on the matching results, the source-end calculation rules actually followed by the online monitoring instrument are determined. The source-end calculation rules characterize the true mathematical logic of calculating monitoring data from raw measurement values ​​within the online monitoring instrument. S4. During the real-time determination process, according to the determined source-end calculation rules, the real-time collected raw measurement values ​​are substituted into the source-end calculation rules to calculate the theoretical concentration value. The theoretical concentration value is then dynamically compared with the monitoring data output in real time by the online monitoring instrument to calculate the degree of deviation. When the degree of deviation exceeds a preset threshold, it is automatically determined that the current monitoring data is suspected of being tampered with.

[0009] Furthermore, S1 collects real-time data from the online water pollution source monitoring instrument in measurement mode, including the output data and source parameters, such as: Using IoT monitoring technology, at the end of each complete measurement cycle of the online water pollution source monitoring instrument, three sets of data streams are read simultaneously: The first data stream is the monitoring data finally output by the online monitoring instrument, which includes at least one of the concentration values ​​of chemical oxygen demand, ammonia nitrogen, total phosphorus, and total nitrogen; the second data stream is the raw measurement value directly collected by the online monitoring instrument during the measurement process, which is the raw absorbance value or luminescence intensity value reflecting the concentration of the analyte during spectrophotometric determination. The raw measurement value exists internally in the online monitoring instrument in the form of an electrical signal or a digital quantized value and is not directly displayed externally; the third data stream is the slope and intercept of the calibration curve used by the online monitoring instrument in the current measurement cycle, where the slope and intercept of the calibration curve are parameter values ​​stored in the internal register of the online monitoring instrument after the most recent calibration operation. All three data streams are accompanied by a precise timestamp and an identifier indicating the current working status of the online monitoring instrument. The working status identifier includes at least measurement status, calibration status, cleaning status, and standby status. Only when the working status identifier is measurement status is the data set marked as valid and stored in the local cache queue. When the working status identifier is not measurement status, the data set is marked as invalid and discarded.

[0010] Furthermore, S2 constructs a multi-dimensional temporal feature vector, including: The data preprocessing module reads raw records marked as valid data in batches from the local cache queue. The raw records at least include timestamps, monitoring data, raw measurement values, calibration curve slopes, and calibration curve intercepts. The data preprocessing module sequentially performs a first cleaning operation, a second cleaning operation, and a third cleaning operation on the raw records. The first cleaning operation involves traversing all raw records and removing records where the monitoring data is null or less than or equal to zero. The second cleaning operation involves removing records where the measured raw values ​​exceed the preset effective measurement range of the online monitoring instrument, where the preset effective measurement range is defined by the lower and upper limits of the online monitoring instrument's range. The third cleaning operation involves removing records where the slope of the calibration curve is less than or equal to zero, and records where the intercept of the calibration curve is null. After the first, second, and third cleaning operations, the data preprocessing module arranges the original records retained after cleaning in ascending order of timestamps to form a continuous time-series dataset. The feature vector construction module extracts the cleaned monitoring data, original measurement values, calibration curve slope, and calibration curve intercept corresponding to each timestamp from the time series dataset, and uses a dynamic adaptive time window algorithm to organize the sliding window. The dynamic adaptive time window algorithm automatically adjusts the window size according to the degree of data fluctuation, including: calculating the coefficient of variation (CV) of the monitoring data of a preset number of data points before the current time. When the CV is less than a first fluctuation threshold, it is determined to be a stable period, and a first window length is used; when the CV is greater than or equal to the first fluctuation threshold and less than a second fluctuation threshold, it is determined to be a normal fluctuation period, and a second window length is used; when the CV is greater than or equal to the second fluctuation threshold, it is determined to be a period of severe fluctuation, and a third window length is used; wherein the first window length is less than the second window length, and the second window length is less than the third window length; the feature vector construction module takes as input all data points within the current dynamic window and outputs as a set of multi-dimensional time-series feature vectors.

[0011] Furthermore, in S3, a regression model is established using the aforementioned multidimensional time-series feature vector to connect the monitoring data, the original measurement values, the slope of the calibration curve, and the intercept of the calibration curve, including: The monitoring data and original measurement values ​​of the multidimensional time series feature vectors under all valid measurement states within a preset time window are extracted to form an initial fitting dataset; a robust fitting operation based on outlier points based on confidence ellipse is performed on the initial fitting dataset. The robust fitting operation for outlier points based on confidence ellipses includes: A two-dimensional feature space is constructed, and the mean vector and covariance matrix of the two-dimensional feature space are calculated. The Mahalanobis distance of each data point is calculated based on the mean vector and the covariance matrix. According to a preset quantile threshold of 2 degrees of freedom for the chi-square distribution, data points whose Mahalanobis distance is greater than the preset quantile threshold are identified as outliers and removed from the initial fitted dataset to obtain a robust fitted dataset. The mean vector, covariance matrix, and preset quantile threshold used for outlier removal are recorded as reference benchmarks for subsequent determinations. Using the robust fitted dataset, a univariate linear regression is performed using the least squares method to obtain the fitting slope and fitting intercept. The fitting slope and fitting intercept are compared with the data sample... The module compares the average slope of the original calibration curve and the average intercept of the original calibration curve reported by the online monitoring instrument within the same preset time window. If the deviation between the fitted parameter and the reported parameter is less than a preset first deviation threshold, it is determined that the regression model accurately reflects the internal calculation logic of the online monitoring instrument, and the regression model is identified as the source calculation rule of the online monitoring instrument. At the same time, the mean vector, the covariance matrix, and the preset quantile critical value used for outlier removal are stored as auxiliary parameters of the source calculation rule. If the deviation between the fitted parameter and the reported parameter is greater than or equal to the preset first deviation threshold, the back-analysis process is triggered.

[0012] Furthermore, the reverse analysis process includes: Using the mean vector, the covariance matrix, and the preset quantile threshold, outlier removal is performed on the dataset within the current time window to obtain a robust dataset; for the robust dataset obtained after outlier removal, a theoretical inverse value sequence is constructed: , Where L is the length of the time window, C u For monitoring data reported by online monitoring instruments, u represents the index corresponding to the time window length, and k u and b u The slope and intercept of the calibration curve reported by the online monitoring instrument; the theoretical back-calculated value sequence A calc The rounding process yields the rounded theoretical back-calculated value sequence, while the original measurement value sequence reported by the online monitoring instrument is also obtained. The empirical distribution consistency test is used to compare the empirical cumulative distribution functions of the rounded theoretical inverse value sequence and the measured original value sequence, and the maximum vertical distance D between the two is calculated. ks : , Among them, F calc(x) is the empirical cumulative distribution function of the theoretical inverse value sequence, Fi raw (x) is the empirical cumulative distribution function of the values ​​reported by the online monitoring instrument, where x represents the original measurement value variable; Based on the critical value of the empirical distribution consistency test, if the maximum vertical distance is less than or equal to the critical value, it is determined that the theoretical back-calculated value sequence and the value sequence reported by the online monitoring instrument follow the same distribution, confirming that the online monitoring instrument follows the slope k of its built-in calibration curve. u and calibration curve intercept b u The determined linear relationship is used as the source calculation rule of the online monitoring instrument, and the recorded mean vector and covariance matrix are stored as verification parameters of the source calculation rule. If D KS If the value exceeds the threshold, a biased distribution estimation based on resampling is employed: multiple resampled samples are drawn with replacement from the original data, each sample being the same size as the original sample, and D is recalculated for each resampled sample. KS Statistic, to obtain D KS The resampling distribution; if the original D KS If the value is above the preset high percentile of the resampling distribution, a systematic bias is determined to exist, triggering the formula library matching process.

[0013] Furthermore, the real-time determination process in S4 includes: After the source-end calculation rules are determined, the regression model enters the daily real-time judgment process; The daily real-time judgment process adopts a cumulative deviation sliding detection algorithm based on cumulative summation; The cumulative deviation sliding detection algorithm includes: defining the residual e at time t. t =C real,t -C calc,t , where C real,t For the real-time monitoring data output by online monitoring instruments, C calc,t The theoretical concentration value calculated according to the source-end calculation rules; cumulative concentration is defined. and statistics : , , Where μ e This is the mean estimate of the residuals during the historical normal period, where K is the allowable drift amplitude parameter, and t-1 represents the time of the previous measurement period; The cumulative deviation sliding detection algorithm maintains a sliding window with a window length set to a preset number of measurement cycles; for each newly arrived data point within the window, the cumulative sum statistic is updated. Define upper control limit H upper =β·σ e and lower control limit H lower =-β·σ e Where β represents the preset control limit coefficient, σ e This represents the standard deviation estimate of the residuals during historical normal periods; when >H upper When a positive cumulative bias is detected, it indicates that the monitoring data is consistently higher than the theoretical value; when <H lower When a negative cumulative deviation is detected, it indicates that the monitoring data is consistently lower than the theoretical value; When the cumulative deviation exceeds the control limit, the deviation confirmation mechanism is triggered: the cumulative and statistical values ​​must remain outside the control limit for a subsequent preset number of measurement periods before the data is finally confirmed as abnormal; if the cumulative and statistical values ​​fall back to within the control limit within a preset number of periods, the current exceedance is marked as a random disturbance and no abnormal alarm is output. When the cumulative deviation exceeds the control limit, after confirmation for a preset number of consecutive cycles, the real-time comparison and judgment module outputs a multi-dimensional anomaly alarm signal. The multi-dimensional anomaly alarm signal includes cumulative and statistical quantities, average deviation, duration, and deviation direction.

[0014] Furthermore, the real-time comparison and judgment module outputs a multi-dimensional anomaly alarm signal. This multi-dimensional anomaly alarm signal is scored using a dynamic scoring mechanism based on Bayesian inference tampering confidence. This Bayesian tampering confidence scoring mechanism includes: The cumulative sum statistic and upper control limit are extracted from the multidimensional anomaly alarm signal, and the cumulative deviation intensity is calculated using the cumulative sum statistic and upper control limit; the Mahalanobis distance of the current data point relative to the established benchmark is calculated from the mean vector and the covariance matrix stored in S3 to obtain the outlier deviation degree; the residual time series autocorrelation is calculated from the deviation direction sequence of the multidimensional anomaly alarm signal. Construct a Bayesian probability model, treating data anomaly detection as a posterior probability inference problem. Let H0 represent the hypothesis that the data is real and has not been tampered with, and H1 represent the hypothesis that the data has been tampered with. The Bayesian probability model integrates multiple sources of evidence to calculate the likelihood ratio and posterior probability, whereby the multiple sources of evidence include the cumulative bias intensity. m The outlier deviation, the residual time series autocorrelation, and the consistency of cross-validation of multiple online monitoring instruments; Through the cumulative deviation intensity E cusu Outlier deviation E mahal Residual time series autocorrelation E autocorr Cross-validation of consistency among multiple online monitoring instrumentscross Calculate the likelihood ratio: , Wherein, P(·|H1) and P(·|H0) are the probability density functions observed under the tampering hypothesis and the truth hypothesis, respectively, and LR represents the likelihood ratio. The probability density function uses the mean vector and the covariance matrix as historical benchmarks to estimate the kernel density. The posterior probability is calculated as follows: , Wherein, P(H1) is the prior tampering probability, which is updated online using the Bayesian update formula based on the result of each judgment; The tampering confidence score is output based on the posterior probability P(H1|E). When the tampering confidence score exceeds the first confidence threshold, a high-confidence tampering alarm is output. When the tampering confidence score is between the second confidence threshold and the first confidence threshold, a medium-confidence tampering warning is output. When the tampering confidence score is lower than the second confidence threshold, the current judgment result is temporarily stored in the pending confirmation queue, and a comprehensive judgment is made after the accumulation of evidence for a preset number of data periods.

[0015] Furthermore, the cross-validation consistency of the multiple online monitoring instruments is calculated according to the following process: When the posterior probability of the first online monitoring instrument A reaches the preset trigger threshold, the multi-online monitoring instrument cross-validation module obtains the original monitoring data of the second online monitoring instrument B and the third online monitoring instrument C with the same monitoring index on the same cross section within the same time window; and independently calculates the posterior probability of the second online monitoring instrument B and the posterior probability of the third online monitoring instrument C according to the cumulative deviation sliding detection algorithm and the Bayesian tampering confidence scoring mechanism described in S4, respectively. Based on the posterior probability I of the first online monitoring instrument A The posterior probability I of the second online monitoring instrument B The posterior probability I of the third online monitoring instrument C Calculate the anomaly consistency index I of the same cross section. consensus =(I A +I B +I C ) / 3; The same cross-sectional anomaly consistency index I consensus As the consensus of the cross-validation of the multiple online monitoring instruments E cross The output value, i.e., E cross =I consensus ; The consistency of the cross-validation of the multiple online monitoring instruments is fed back into the Bayesian probability model to update the likelihood ratio of the first online monitoring instrument A. The deviation of the first online monitoring instrument A relative to the second online monitoring instrument B and the third online monitoring instrument C is calculated. The updated posterior probability is corrected based on the deviation. The corrected posterior probability is used as the basis for the final tamper confidence score. The results of this cross-validation are stored in the knowledge graph associated with the online monitoring instruments.

[0016] Furthermore, when the tampering confidence score reaches the first confidence threshold and the cross-validation consistency of the multiple online monitoring instruments indicates that at least two online monitoring instruments are simultaneously abnormal, the data points marked as normal in the most recent maintenance cycle of the current online monitoring instrument are extracted from the knowledge graph associated with the online monitoring instruments, and the source-end calculation rules are self-optimized and updated. The specific process of the self-optimization update is as follows: using the incremental Bayesian linear regression method, the current slope coefficient and intercept coefficient of the source calculation rule are used as the mean of the prior distribution, the normal data points extracted from the knowledge graph associated with the online monitoring instrument are used as observation samples, the posterior distribution is calculated, and the updated slope coefficient and intercept coefficient are obtained based on the posterior distribution. After the update is completed, the updated slope coefficient and intercept coefficient are stored as new source-end calculation rules, and the baseline parameters of the corresponding online monitoring instruments in the online monitoring instrument association knowledge graph are updated synchronously.

[0017] Compared with the prior art, the present invention achieves the following beneficial effects: This invention starts from the source of the operating principle of online monitoring instruments, collecting monitoring data in real time, as well as source parameters such as raw measurement values, calibration curve slope, and calibration curve intercept. The monitoring data and source parameters at the same moment are organized into a multi-dimensional time-series feature vector. A regression model is established using multivariate statistics to establish the relationship between the monitoring data and the raw measurement values, calibration curve slope, and calibration curve intercept. The fitted parameters are matched with the parameters reported by the online monitoring instrument to determine the source-end calculation rules actually followed by the online monitoring instrument. In real-time judgment, the raw measurement values ​​are substituted into the source-end calculation rules to calculate the theoretical concentration, which is dynamically compared with the monitoring data output by the online monitoring instrument. Based on the degree of deviation, suspected human tampering is identified. This invention improves the accuracy and reliability of tampering detection by verifying the physical and logical consistency between the internal data of the online monitoring instrument and the source parameters, addressing the root cause of data generation. Attached Figure Description

[0018] Figure 1 This is a flowchart illustrating a method for determining the accuracy of water pollution monitoring data based on multivariate statistics, provided in an embodiment of the present invention. Figure 2This is a flowchart of the regression model initialization process for a method for determining the accuracy of water pollution monitoring data based on multivariate statistics, provided in an embodiment of the present invention. Figure 3 This is a daily real-time judgment flowchart of a method for determining the accuracy of water pollution monitoring data based on multivariate statistics, provided in an embodiment of the present invention. Detailed Implementation

[0019] The present invention will now be described in detail with reference to the accompanying drawings.

[0020] Example 1

[0021] This invention provides a method for determining the accuracy of water pollution monitoring data based on multivariate statistics, such as... Figure 1 As shown, the method includes: S1. Real-time acquisition of result data and source parameters output by the online water pollution source monitoring instrument in measurement mode; the result data is the final monitoring data output by the online monitoring instrument; the source parameters include the original measurement value, the slope of the calibration curve, and the intercept of the calibration curve; wherein, the original measurement value is the absorbance value or luminescence intensity value directly measured during the spectrophotometric determination process without concentration conversion; the slope of the calibration curve and the intercept of the calibration curve are key conversion parameters built into the online monitoring instrument to convert the original measurement value into monitoring data; S2. After data cleaning and denoising, the monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept at the same measurement time are organized according to the time window to construct a multi-dimensional time series feature vector. S3. Based on multivariate statistical analysis, a regression model is established using the multidimensional time-series feature vector to connect monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept. The regression model matches the fitted parameters with the source parameters reported by the online monitoring instrument. Based on the matching results, the source-end calculation rules actually followed by the online monitoring instrument are determined. The source-end calculation rules characterize the true mathematical logic of calculating monitoring data from raw measurement values ​​within the online monitoring instrument. It should be noted that the source-end calculation rule refers to the mathematical logic actually used inside the online water pollution source monitoring instrument to convert the measured raw values ​​(absorbance or luminescence intensity) into the final monitoring data (concentration value). This rule is the core embodiment of the instrument's operating principle and is usually expressed as a univariate linear formula or a polynomial formula; The following explanation uses a total phosphorus online monitoring instrument (quadratic polynomial rule): Instrument model: TP-100 online monitoring instrument (high concentration range nonlinear response); Calibration parameters (obtained through multi-point calibration): Quadratic term coefficient a = 0.003. The coefficient of the linear term is b = 0.52, and the coefficient of the constant term is c = 0.10; the source-end calculation rule is: C = 0.003 × A 2 +0.52×A+0.10, a complete measurement process: molybdenum-antimony spectrophotometric determination, obtaining absorbance A=120.0, total phosphorus online monitor calculation: C = 0.003 × 120 2 +0.52×120+0.10=43.2+62.4+0.10=105.7, Instrument output monitoring data: 105.7 mg / L; S4. During the real-time determination process, according to the determined source-end calculation rules, the real-time collected raw measurement values ​​are substituted into the source-end calculation rules to calculate the theoretical concentration value. The theoretical concentration value is then dynamically compared with the monitoring data output in real time by the online monitoring instrument to calculate the degree of deviation. When the degree of deviation exceeds a preset threshold, it is automatically determined that the current monitoring data is suspected of being tampered with.

[0022] Specifically, this embodiment starts from the source of the online monitoring instrument's operating principle, collecting monitoring data in real time, as well as source parameters such as raw measurement values, calibration curve slope, and calibration curve intercept. The monitoring data and source parameters at the same time are organized into a multi-dimensional time-series feature vector. A regression model is established using multivariate statistics between the monitoring data and the raw measurement values, calibration curve slope, and calibration curve intercept. The fitted parameters are matched for consistency with the parameters reported by the online monitoring instrument to determine the source-end calculation rules actually followed by the online monitoring instrument. During real-time judgment, the raw measurement values ​​are substituted into the source-end calculation rules to calculate the theoretical concentration, which is dynamically compared with the monitoring data output by the online monitoring instrument. Based on the degree of deviation, suspected human tampering is identified. This embodiment improves the accuracy and reliability of tampering judgment by verifying the physical and logical consistency between the data internal to the online monitoring instrument and the source parameters, addressing the root cause of data generation.

[0023] It should be noted that in online monitoring instruments for water pollution sources based on spectrophotometry (such as COD, ammonia nitrogen, total phosphorus, and total nitrogen analyzers), the slope and intercept of the calibration curve are two core conversion parameters built into the instrument. Together, they constitute a linear function that converts the original measured value (absorbance or luminescence intensity) into the final monitoring data (concentration value). The slope k of the calibration curve represents the concentration change corresponding to a unit change in the original measured value, reflecting the instrument's sensitivity to the analyte. Taking a COD monitoring instrument as an example, after calibration, the slope of the calibration curve of a certain online COD monitor is k=0.85. This means that when the original measured value (absorbance) changes by 1 unit, the COD concentration changes by 0.85 mg / L. The physical meaning is that the instrument's sensitivity to COD is 0.85. The intercept b0 of the calibration curve represents the concentration value corresponding to the original measured value being zero, reflecting the zero-point drift or background interference of the instrument system. Taking an ammonia nitrogen monitor as an example, the intercept of the calibration curve of a certain ammonia nitrogen monitor is b0 = -1.5 (negative value). This means that when the absorbance is zero, the concentration is displayed as a negative number. The judgment is that the instrument's zero-point drift is serious and it needs to be recalibrated immediately.

[0024] In the above embodiments, specifically, the real-time acquisition of result data and source parameters output by the online water pollution source monitoring instrument in measurement state during S1 includes: Using IoT monitoring technology, at the end of each complete measurement cycle of the online water pollution source monitoring instrument, three sets of data streams are read simultaneously: The first data stream is the monitoring data finally output by the online monitoring instrument, which includes at least one of the concentration values ​​of chemical oxygen demand, ammonia nitrogen, total phosphorus, and total nitrogen; the second data stream is the raw measurement value directly collected by the online monitoring instrument during the measurement process, which is the raw absorbance value or luminescence intensity value reflecting the concentration of the analyte during spectrophotometric determination. The raw measurement value exists internally in the online monitoring instrument in the form of an electrical signal or a digital quantized value and is not directly displayed externally; the third data stream is the slope and intercept of the calibration curve used by the online monitoring instrument in the current measurement cycle, where the slope and intercept of the calibration curve are parameter values ​​stored in the internal register of the online monitoring instrument after the most recent calibration operation. All three data streams are accompanied by a precise timestamp and an identifier indicating the current working status of the online monitoring instrument. The working status identifier includes at least measurement status, calibration status, cleaning status, and standby status. Only when the working status identifier is measurement status is the data set marked as valid and stored in the local cache queue. When the working status identifier is not measurement status, the data set is marked as invalid and discarded.

[0025] In the above embodiments, specifically, the construction of a multi-dimensional temporal feature vector in S2 includes: The data preprocessing module reads raw records marked as valid data in batches from the local cache queue. The raw records at least include timestamps, monitoring data, raw measurement values, calibration curve slopes, and calibration curve intercepts. The data preprocessing module sequentially performs a first cleaning operation, a second cleaning operation, and a third cleaning operation on the raw records. The first cleaning operation involves traversing all raw records and removing records where the monitoring data is null or less than or equal to zero. The second cleaning operation involves removing records where the measured raw values ​​exceed the preset effective measurement range of the online monitoring instrument, where the preset effective measurement range is defined by the lower and upper limits of the online monitoring instrument's range. The third cleaning operation involves removing records where the slope of the calibration curve is less than or equal to zero, and records where the intercept of the calibration curve is null. After the first, second, and third cleaning operations, the data preprocessing module arranges the original records retained after cleaning in ascending order of timestamps to form a continuous time-series dataset. The feature vector construction module extracts the cleaned monitoring data, original measurement values, calibration curve slope, and calibration curve intercept corresponding to each timestamp from the time series dataset, and uses a dynamic adaptive time window algorithm to organize the sliding window. The dynamic adaptive time window algorithm automatically adjusts the window size according to the degree of data fluctuation, including: calculating the coefficient of variation (CV) of the monitoring data of a preset number of data points before the current time. When the CV is less than a first fluctuation threshold, it is determined to be a stable period, and a first window length is used; when the CV is greater than or equal to the first fluctuation threshold and less than a second fluctuation threshold, it is determined to be a normal fluctuation period, and a second window length is used; when the CV is greater than or equal to the second fluctuation threshold, it is determined to be a period of severe fluctuation, and a third window length is used; wherein the first window length is less than the second window length, and the second window length is less than the third window length; the feature vector construction module takes as input all data points within the current dynamic window and outputs as a set of multi-dimensional time-series feature vectors.

[0026] It should be noted that this embodiment utilizes multivariate statistical analysis to establish a regression model between the monitoring data and the original measured values, the slope of the calibration curve, and the intercept of the calibration curve. By comparing the consistency between the fitted parameters and the original parameters reported by the instrument, it is determined whether the monitoring data has been tampered with. The tampering detection principle is as follows: if the monitoring data is manually modified (e.g., changed from 40.25 to 30.00), but the original measured value A and the calibration parameters k and b0 are not modified simultaneously, then the fitted k will be tampered with. fit and b fit will be compared with the k reported by the instrument. org and b0 org Significant deviations occur; this embodiment identifies human tampering by detecting the failure of this inherent logical relationship.

[0027] In the above embodiments, specifically, S3 utilizes the multidimensional time-series feature vector to establish a regression model between monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept, including: The monitoring data and original measurement values ​​of the multidimensional time series feature vectors under all valid measurement states within a preset time window are extracted to form an initial fitting dataset; a robust fitting operation based on outlier points based on confidence ellipse is performed on the initial fitting dataset. The robust fitting operation for outlier points based on confidence ellipses includes: A two-dimensional feature space is constructed, and the mean vector and covariance matrix of the two-dimensional feature space are calculated. The Mahalanobis distance of each data point is calculated based on the mean vector and the covariance matrix. According to a preset quantile threshold of 2 degrees of freedom for the chi-square distribution, data points whose Mahalanobis distance is greater than the preset quantile threshold are identified as outliers and removed from the initial fitted dataset to obtain a robust fitted dataset. The mean vector, covariance matrix, and preset quantile threshold used for outlier removal are recorded as reference benchmarks for subsequent determinations. Using the robust fitted dataset, a univariate linear regression is performed using the least squares method to obtain the fitting slope and fitting intercept. The fitting slope and fitting intercept are compared with the data sample... The module compares the average slope of the original calibration curve and the average intercept of the original calibration curve reported by the online monitoring instrument within the same preset time window. If the deviation between the fitted parameter and the reported parameter is less than a preset first deviation threshold, it is determined that the regression model accurately reflects the internal calculation logic of the online monitoring instrument, and the regression model is identified as the source calculation rule of the online monitoring instrument. At the same time, the mean vector, the covariance matrix, and the preset quantile critical value used for outlier removal are stored as auxiliary parameters of the source calculation rule. If the deviation between the fitted parameter and the reported parameter is greater than or equal to the preset first deviation threshold, the back-analysis process is triggered.

[0028] In the above embodiments, specifically, the back-analysis process includes: Using the mean vector, the covariance matrix, and the preset quantile threshold, perform the same outlier removal operation as in step S3 on the dataset within the current time window to obtain a robust dataset; for the robust dataset obtained after outlier removal, construct a theoretical inverse value sequence: , Where L is the length of the time window, C u For monitoring data reported by online monitoring instruments, u represents the index corresponding to the time window length, and k u and b u The slope and intercept of the calibration curve reported by the online monitoring instrument; the theoretical back-calculated value sequence A calcThe rounding process yields the rounded theoretical back-calculated value sequence, while the original measurement value sequence reported by the online monitoring instrument is also obtained. The empirical distribution consistency test is used to compare the empirical cumulative distribution functions of the rounded theoretical inverse value sequence and the measured original value sequence, and the maximum vertical distance D between the two is calculated. ks : , Among them, F calc (x) is the empirical cumulative distribution function of the theoretical inverse value sequence, Fi raw (x) is the empirical cumulative distribution function of the values ​​reported by the online monitoring instrument, where x represents the original measurement value variable; Based on the critical value of the empirical distribution consistency test, if the maximum vertical distance is less than or equal to the critical value, it is determined that the theoretical back-calculated value sequence and the value sequence reported by the online monitoring instrument follow the same distribution, confirming that the online monitoring instrument follows the slope k of its built-in calibration curve. u and calibration curve intercept b u The determined linear relationship is used as the source calculation rule of the online monitoring instrument, and the recorded mean vector and covariance matrix are stored as verification parameters of the source calculation rule. If D KS If the value exceeds the threshold, a biased distribution estimation based on resampling is employed: multiple resampled samples are drawn with replacement from the original data, each sample being the same size as the original sample, and D is recalculated for each resampled sample. KS Statistic, to obtain D KS The resampling distribution; if the original D KS If the value is above the preset high percentile of the resampling distribution, a systematic bias is determined to exist, triggering the formula library matching process.

[0029] In the above embodiments, specifically, the real-time determination process in S4 includes: After completing the regression model initialization and determining the source-end calculation rules in step S3, the regression model enters the daily real-time judgment process; The daily real-time judgment process adopts a cumulative deviation sliding detection algorithm based on cumulative sum to replace single-point deviation comparison, so as to enhance the detection sensitivity of continuous small-amplitude tampering. The cumulative deviation sliding detection algorithm includes: defining the residual e at time t. t =C real,t -C calc,t , where C real,t For the real-time monitoring data output by online monitoring instruments, C calc,t The theoretical concentration value calculated according to the source-end calculation rules; cumulative concentration is defined. and statistics : , , Where μ e This is an estimate of the mean of the residuals during historical normal periods, where K is the allowable drift magnitude parameter; The cumulative deviation sliding detection algorithm maintains a sliding window with a window length set to a preset number of measurement cycles; for each newly arrived data point within the window, the cumulative sum statistic is updated. Define upper control limit H upper =β·σ e and lower control limit H lower =-β·σ e Where β represents the preset control limit coefficient, σ e This represents the standard deviation estimate of the residuals during historical normal periods; when >H upper When a positive cumulative bias is detected, it indicates that the monitoring data is consistently higher than the theoretical value; when <H lower When a negative cumulative deviation is detected, it indicates that the monitoring data is consistently lower than the theoretical value; When the cumulative deviation exceeds the control limit, the deviation confirmation mechanism is triggered: the cumulative and statistical values ​​must remain outside the control limit for a subsequent preset number of measurement periods before the data is finally confirmed as abnormal; if the cumulative and statistical values ​​fall back to within the control limit within a preset number of periods, the current exceedance is marked as a random disturbance and no abnormal alarm is output. When the cumulative deviation exceeds the control limit, after confirmation for a preset number of consecutive cycles, the real-time comparison and judgment module outputs a multi-dimensional abnormal alarm signal that includes the cumulative sum and statistics, average deviation, duration, and deviation direction.

[0030] It should be noted that the regression model includes a regression model initialization process and a daily real-time judgment process; Step 1: Regression Model Initialization Process: 1. Data cleaning: Extract monitoring data from online water pollution source monitoring instruments under measurement conditions for a period of time, including raw measurement values, calibration curve slope, calibration curve intercept, or non-empty data such as quadratic coefficient a, linear coefficient b, and constant coefficient c. 2. Data Acquisition: Acquire monitoring data, raw measurement values, calibration curve slope, calibration curve intercept or quadratic coefficient a, one-time coefficient b, and constant coefficient c at the same time. 3. Site Data Pattern Judgment: Based on the s built-in calculation formulas in the model, the original measurement value (x1) is calculated simultaneously, and the deviation of the original measurement value (diff) is calculated. If n% of diff ≤ m%, then this formula is the calculation formula for the monitoring data of the online monitoring instrument of the water pollution source at the online monitoring site. If the deviation of the original measurement value calculated by the s built-in calculation formulas does not meet the condition that n% of diff ≤ m1%, then the output pattern needs to be investigated and manually analyzed to improve the linear regression model. Step 2: Daily Real-Time Judgment Process: 1. Real-time data change judgment: Determine whether the real-time data has changed in real time. Only if the real-time data has changed should the following process continue. 2. Status Judgment: Only data generated under measurement status is judged. If it is under measurement status, the process continues; if it is not under measurement status, the process ends. 3. Data Acquisition: Acquire monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept (for quadratic functions, this includes the quadratic coefficient a, the linear coefficient b, and the constant coefficient c) at the same time. Determine if any parameters are missing. If so, the process ends; otherwise, continue. 4. Data accuracy judgment: Based on the calculation formula of the monitoring data of the online monitoring instrument of the water pollution source determined in the initialization process, calculate the original measurement value (x1), and then calculate the deviation of the original measurement value (diff). If diff≤m2%, the process ends; if diff>m2%, the output data is abnormal. Here, m1 and m2 both represent thresholds.

[0031] In the above embodiments, specifically, the real-time comparison and judgment module outputs a multi-dimensional anomaly alarm signal, and scores the multi-dimensional anomaly alarm signal by introducing a dynamic scoring mechanism for tamper confidence based on Bayesian inference. The Bayesian tamper confidence scoring mechanism includes: The cumulative sum statistic and upper control limit are extracted from the multidimensional anomaly alarm signal, and the cumulative deviation intensity is calculated using the cumulative sum statistic and upper control limit; the Mahalanobis distance of the current data point relative to the established benchmark is calculated from the mean vector and the covariance matrix stored in step S3 to obtain the outlier deviation degree; the residual time series autocorrelation is calculated from the deviation direction sequence of the multidimensional anomaly alarm signal. Construct a Bayesian probability model, treating data anomaly detection as a posterior probability inference problem. Let H0 represent the hypothesis that the data is real and has not been tampered with, and H1 represent the hypothesis that the data has been tampered with. The Bayesian probability model integrates multiple sources of evidence to calculate the likelihood ratio and posterior probability, whereby the multiple sources of evidence include the cumulative bias intensity. m The outlier deviation, the residual time series autocorrelation, and the consistency of cross-validation of multiple online monitoring instruments; Through the cumulative deviation intensity E cusu Outlier deviation E mahal Residual time series autocorrelation E autocorr Cross-validation of consistency among multiple online monitoring instruments cross Calculate the likelihood ratio: , Wherein, P(·|H1) and P(·|H0) are the probability density functions observed under the tampering hypothesis and the truth hypothesis, respectively, and LR represents the likelihood ratio. The probability density function uses the mean vector and the covariance matrix as historical benchmarks to estimate the kernel density. The posterior probability is calculated as follows: , Where P(H1) is the prior tampering probability, which is updated online using the Bayesian update formula based on the result of each judgment, and E represents the multidimensional evidence source; The tampering confidence score is output based on the posterior probability P(H1|E). When the tampering confidence score exceeds the first confidence threshold, a high-confidence tampering alarm is output. When the tampering confidence score is between the second confidence threshold and the first confidence threshold, a medium-confidence tampering warning is output. When the tampering confidence score is lower than the second confidence threshold, the current judgment result is temporarily stored in the pending confirmation queue, and a comprehensive judgment is made after the accumulation of evidence for a preset number of data periods.

[0032] In the above embodiments, specifically, the cross-validation consistency of the multiple online monitoring instruments is calculated according to the following process: When the posterior probability of the first online monitoring instrument A reaches the preset trigger threshold, the multi-online monitoring instrument cross-validation module obtains the original monitoring data of the second online monitoring instrument B and the third online monitoring instrument C with the same monitoring index on the same cross section within the same time window; and independently calculates the posterior probability of the second online monitoring instrument B and the posterior probability of the third online monitoring instrument C according to the cumulative deviation sliding detection algorithm and the Bayesian tampering confidence scoring mechanism described in S4, respectively. Based on the posterior probability I of the first online monitoring instrument A The posterior probability I of the second online monitoring instrument B The posterior probability I of the third online monitoring instrument C Calculate the anomaly consistency index I of the same cross section. consensus =(I A +I B +I C ) / 3; The same cross-sectional anomaly consistency index I consensus As the consensus of the cross-validation of the multiple online monitoring instruments Ecross The output value, i.e., E cross =I consensus ; The consistency of the cross-validation of the multiple online monitoring instruments is fed back into the Bayesian probability model to update the likelihood ratio of the first online monitoring instrument A. The deviation of the first online monitoring instrument A relative to the second online monitoring instrument B and the third online monitoring instrument C is calculated. The updated posterior probability is corrected based on the deviation. The corrected posterior probability is used as the basis for the final tamper confidence score. The results of this cross-validation are stored in the knowledge graph associated with the online monitoring instruments.

[0033] In the above embodiments, specifically, when the tampering confidence score reaches the first confidence threshold and the cross-validation consistency of the multiple online monitoring instruments indicates that at least two online monitoring instruments are simultaneously abnormal, the data points marked as normal in the most recent maintenance cycle of the current online monitoring instrument are extracted from the knowledge graph associated with the online monitoring instruments, and the source-end calculation rules are self-optimized and updated. The specific process of the self-optimization update is as follows: using the incremental Bayesian linear regression method, the current slope coefficient and intercept coefficient of the source calculation rule are used as the mean of the prior distribution, the normal data points extracted from the knowledge graph associated with the online monitoring instrument are used as observation samples, the posterior distribution is calculated, and the updated slope coefficient and intercept coefficient are obtained based on the posterior distribution. After the update is completed, the updated slope coefficient and intercept coefficient are stored as new source-end calculation rules, and the baseline parameters of the corresponding online monitoring instruments in the online monitoring instrument association knowledge graph are updated synchronously.

[0034] It should be understood that the above embodiments are one or more embodiments of the present invention, and there are many other embodiments and variations based on the present invention; any variations and modifications made by those skilled in the art through the present invention without making pioneering innovations are all within the protection scope of the present invention.

Claims

1. A method for determining the accuracy of water pollution monitoring data based on multivariate statistics, characterized in that, The method includes: S1. Real-time acquisition of result data and source parameters output by the online water pollution source monitoring instrument in measurement mode; the result data is the final monitoring data output by the online monitoring instrument; the source parameters include the original measurement value, the slope of the calibration curve, and the intercept of the calibration curve; S2. After data cleaning and denoising, the monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept at the same measurement time are organized according to the time window to construct a multi-dimensional time series feature vector. S3. Based on multivariate statistical analysis, a regression model is established using the multidimensional time-series feature vector to connect monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept. The regression model matches the fitted parameters with the source parameters reported by the online monitoring instrument. Based on the matching results, the source-end calculation rules actually followed by the online monitoring instrument are determined. The source-end calculation rules characterize the true mathematical logic of calculating monitoring data from raw measurement values ​​within the online monitoring instrument. S4. During the real-time determination process, according to the determined source-end calculation rules, the real-time collected raw measurement values ​​are substituted into the source-end calculation rules to calculate the theoretical concentration value. The theoretical concentration value is then dynamically compared with the monitoring data output in real time by the online monitoring instrument to calculate the degree of deviation. When the degree of deviation exceeds a preset threshold, it is automatically determined that the current monitoring data is suspected of being tampered with.

2. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 1, characterized in that, S1 collects real-time data from the online water pollution source monitoring instrument in measurement mode, including the output data and source parameters, such as: Using IoT monitoring technology, at the end of each complete measurement cycle of the online water pollution source monitoring instrument, three sets of data streams are read simultaneously: The first data stream is the monitoring data finally output by the online monitoring instrument, which includes at least one of the concentration values ​​of chemical oxygen demand, ammonia nitrogen, total phosphorus, and total nitrogen; the second data stream is the raw measurement value directly collected by the online monitoring instrument during the measurement process, which is the raw absorbance value or luminescence intensity value reflecting the concentration of the analyte during spectrophotometric determination. The raw measurement value exists internally in the online monitoring instrument in the form of an electrical signal or a digital quantized value and is not directly displayed externally; the third data stream is the slope and intercept of the calibration curve used by the online monitoring instrument in the current measurement cycle, where the slope and intercept of the calibration curve are parameter values ​​stored in the internal register of the online monitoring instrument after the most recent calibration operation. All three data streams are accompanied by a precise timestamp and an identifier indicating the current working status of the online monitoring instrument. The working status identifier includes at least measurement status, calibration status, cleaning status, and standby status. Only when the working status identifier is measurement status is the data set marked as valid and stored in the local cache queue. When the working status identifier is not measurement status, the data set is marked as invalid and discarded.

3. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 2, characterized in that, S2 constructs a multi-dimensional temporal feature vector, including: The data preprocessing module reads raw records marked as valid data in batches from the local cache queue. The raw records at least include timestamps, monitoring data, raw measurement values, calibration curve slopes, and calibration curve intercepts. The data preprocessing module sequentially performs a first cleaning operation, a second cleaning operation, and a third cleaning operation on the raw records. The first cleaning operation involves traversing all raw records and removing records where the monitoring data is null or less than or equal to zero. The second cleaning operation involves removing records where the measured raw values ​​exceed the preset effective measurement range of the online monitoring instrument, where the preset effective measurement range is defined by the lower and upper limits of the online monitoring instrument's range. The third cleaning operation involves removing records where the slope of the calibration curve is less than or equal to zero, and records where the intercept of the calibration curve is null. After the first, second, and third cleaning operations, the data preprocessing module arranges the original records retained after cleaning in ascending order of timestamps to form a continuous time-series dataset. The feature vector construction module extracts the cleaned monitoring data, original measurement values, calibration curve slope, and calibration curve intercept corresponding to each timestamp from the time series dataset, and uses a dynamic adaptive time window algorithm to organize the sliding window. The dynamic adaptive time window algorithm automatically adjusts the window size according to the degree of data fluctuation, including: calculating the coefficient of variation (CV) of the monitoring data of a preset number of data points before the current time. When the CV is less than a first fluctuation threshold, it is determined to be a stable period, and a first window length is used; when the CV is greater than or equal to the first fluctuation threshold and less than a second fluctuation threshold, it is determined to be a normal fluctuation period, and a second window length is used; when the CV is greater than or equal to the second fluctuation threshold, it is determined to be a period of severe fluctuation, and a third window length is used; wherein the first window length is less than the second window length, and the second window length is less than the third window length; the feature vector construction module takes as input all data points within the current dynamic window and outputs as a set of multi-dimensional time-series feature vectors.

4. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 3, characterized in that, S3 utilizes the aforementioned multidimensional time-series feature vector to establish a regression model between monitoring data, raw measurement values, calibration curve slope, and calibration curve intercept, including: The monitoring data and original measurement values ​​of the multidimensional time series feature vectors under all valid measurement states within a preset time window are extracted to form an initial fitting dataset; a robust fitting operation based on outlier points based on confidence ellipse is performed on the initial fitting dataset. The robust fitting operation for outlier points based on confidence ellipses includes: A two-dimensional feature space is constructed, and the mean vector and covariance matrix of the two-dimensional feature space are calculated. The Mahalanobis distance of each data point is calculated based on the mean vector and the covariance matrix. According to a preset quantile threshold of 2 degrees of freedom for the chi-square distribution, data points whose Mahalanobis distance is greater than the preset quantile threshold are identified as outliers and removed from the initial fitted dataset to obtain a robust fitted dataset. The mean vector, covariance matrix, and preset quantile threshold used for outlier removal are recorded as reference benchmarks for subsequent determinations. Using the robust fitted dataset, a univariate linear regression is performed using the least squares method to obtain the fitting slope and fitting intercept. The fitting slope and fitting intercept are compared with the data sample... The module compares the average slope of the original calibration curve and the average intercept of the original calibration curve reported by the online monitoring instrument within the same preset time window. If the deviation between the fitted parameter and the reported parameter is less than a preset first deviation threshold, it is determined that the regression model accurately reflects the internal calculation logic of the online monitoring instrument, and the regression model is identified as the source calculation rule of the online monitoring instrument. At the same time, the mean vector, the covariance matrix, and the preset quantile critical value used for outlier removal are stored as auxiliary parameters of the source calculation rule. If the deviation between the fitted parameter and the reported parameter is greater than or equal to the preset first deviation threshold, the back-analysis process is triggered.

5. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 4, characterized in that, The reverse analysis process includes: Using the mean vector, the covariance matrix, and the preset quantile threshold, outlier removal is performed on the dataset within the current time window to obtain a robust dataset; for the robust dataset obtained after outlier removal, a theoretical inverse value sequence is constructed: ， Where L is the length of the time window, C u For monitoring data reported by online monitoring instruments, u represents the index corresponding to the time window length, and k u and b u The slope and intercept of the calibration curve reported by the online monitoring instrument; the theoretical back-calculated value sequence A calc The rounding process yields the rounded theoretical back-calculated value sequence, while the original measurement value sequence reported by the online monitoring instrument is also obtained. The empirical distribution consistency test is used to compare the empirical cumulative distribution functions of the rounded theoretical inverse value sequence and the measured original value sequence, and the maximum vertical distance D between the two is calculated. ks : ， Among them, F calc (x) is the empirical cumulative distribution function of the theoretical inverse value sequence, Fi raw (x) is the empirical cumulative distribution function of the values ​​reported by the online monitoring instrument, where x represents the original measurement value variable; Based on the critical value of the empirical distribution consistency test, if the maximum vertical distance is less than or equal to the critical value, it is determined that the theoretical back-calculated value sequence and the value sequence reported by the online monitoring instrument follow the same distribution, confirming that the online monitoring instrument follows the slope k of its built-in calibration curve. u and calibration curve intercept b u The determined linear relationship is used as the source calculation rule of the online monitoring instrument, and the recorded mean vector and covariance matrix are stored as verification parameters of the source calculation rule. If D KS If the value exceeds the threshold, a biased distribution estimation based on resampling is employed: multiple resampled samples are drawn with replacement from the original data, each sample being the same size as the original sample, and D is recalculated for each resampled sample. KS Statistic, to obtain D KS The resampling distribution; if the original D KS If the value is above the preset high percentile of the resampling distribution, a systematic bias is determined to exist, triggering the formula library matching process.

6. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 5, characterized in that, The real-time determination process in S4 includes: After the source-end calculation rules are determined, the regression model enters the daily real-time judgment process; The daily real-time judgment process adopts a cumulative deviation sliding detection algorithm based on cumulative summation; The cumulative deviation sliding detection algorithm includes: defining the residual e at time t. t =C real,t -C calc,t , where C real,t For the real-time monitoring data output by online monitoring instruments, C calc,t The theoretical concentration value calculated according to the source-end calculation rules; cumulative concentration is defined. and statistics : ， ， Where μ e This is the mean estimate of the residuals during the historical normal period, where K is the allowable drift amplitude parameter, and t-1 represents the time of the previous measurement period; The cumulative deviation sliding detection algorithm maintains a sliding window with a window length set to a preset number of measurement cycles; for each newly arrived data point within the window, the cumulative sum statistic is updated. Define upper control limit H upper =β·σ e and lower control limit H lower =-β·σ e Where β represents the preset control limit coefficient, σ e This represents the standard deviation estimate of the residuals during historical normal periods; when >H upper When a positive cumulative bias is detected, it indicates that the monitoring data is consistently higher than the theoretical value; when <H lower When a negative cumulative deviation is detected, it indicates that the monitoring data is consistently lower than the theoretical value. When the cumulative deviation exceeds the control limit, the deviation confirmation mechanism is triggered: the cumulative and statistical values ​​must remain outside the control limit for a subsequent preset number of measurement periods before the data is finally confirmed as abnormal; if the cumulative and statistical values ​​fall back to within the control limit within a preset number of periods, the current exceedance is marked as a random disturbance and no abnormal alarm is output. When the cumulative deviation exceeds the control limit, after confirmation for a preset number of consecutive cycles, the real-time comparison and judgment module outputs a multi-dimensional anomaly alarm signal. The multi-dimensional anomaly alarm signal includes cumulative and statistical quantities, average deviation, duration, and deviation direction.

7. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 6, characterized in that, The real-time comparison and judgment module outputs a multi-dimensional anomaly alarm signal. A dynamic scoring mechanism based on Bayesian inference for tampering confidence is introduced to score the multi-dimensional anomaly alarm signal. This Bayesian tampering confidence scoring mechanism includes: The cumulative sum statistic and upper control limit are extracted from the multidimensional anomaly alarm signal, and the cumulative deviation intensity is calculated using the cumulative sum statistic and upper control limit; the Mahalanobis distance of the current data point relative to the established benchmark is calculated from the mean vector and the covariance matrix stored in S3 to obtain the outlier deviation degree; the residual time series autocorrelation is calculated from the deviation direction sequence of the multidimensional anomaly alarm signal. Construct a Bayesian probability model, treating data anomaly detection as a posterior probability inference problem. Let H0 represent the hypothesis that the data is real and has not been tampered with, and H1 represent the hypothesis that the data has been tampered with. The Bayesian probability model integrates multiple sources of evidence to calculate the likelihood ratio and posterior probability, whereby the multiple sources of evidence include the cumulative bias intensity. m The outlier deviation, the residual time series autocorrelation, and the consistency of cross-validation of multiple online monitoring instruments; Through the cumulative deviation intensity E cusu Outlier deviation E mahal Residual time series autocorrelation E autocorr Cross-validation of consistency among multiple online monitoring instruments cross Calculate the likelihood ratio: ， Wherein, P(·|H1) and P(·|H0) are the probability density functions observed under the tampering hypothesis and the truth hypothesis, respectively, and LR represents the likelihood ratio. The probability density function uses the mean vector and the covariance matrix as historical benchmarks to estimate the kernel density. The posterior probability is calculated as follows: ， Wherein, P(H1) is the prior tampering probability, which is updated online using the Bayesian update formula based on the result of each judgment; The tampering confidence score is output based on the posterior probability P(H1|E). When the tampering confidence score exceeds the first confidence threshold, a high-confidence tampering alarm is output. When the tampering confidence score is between the second confidence threshold and the first confidence threshold, a medium-confidence tampering warning is output. When the tampering confidence score is lower than the second confidence threshold, the current judgment result is temporarily stored in the pending confirmation queue, and a comprehensive judgment is made after the accumulation of evidence for a preset number of data periods.

8. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 7, characterized in that, The cross-validation consistency of the multiple online monitoring instruments is calculated according to the following process: When the posterior probability of the first online monitoring instrument A reaches the preset trigger threshold, the multi-online monitoring instrument cross-validation module obtains the original monitoring data of the second online monitoring instrument B and the third online monitoring instrument C with the same monitoring index on the same cross section within the same time window; and independently calculates the posterior probability of the second online monitoring instrument B and the posterior probability of the third online monitoring instrument C according to the cumulative deviation sliding detection algorithm and the Bayesian tampering confidence scoring mechanism described in S4, respectively. Based on the posterior probability I of the first online monitoring instrument A The posterior probability I of the second online monitoring instrument B The posterior probability I of the third online monitoring instrument C Calculate the anomaly consistency index I of the same cross section. consensus =(I A +I B +I C ) / 3; The same cross-sectional anomaly consistency index I consensus As the consensus of the cross-validation of the multiple online monitoring instruments E cross The output value, i.e., E cross =I consensus ; The consistency of the cross-validation of the multiple online monitoring instruments is fed back into the Bayesian probability model to update the likelihood ratio of the first online monitoring instrument A. The deviation of the first online monitoring instrument A relative to the second online monitoring instrument B and the third online monitoring instrument C is calculated. The updated posterior probability is corrected based on the deviation. The corrected posterior probability is used as the basis for the final tamper confidence score. The results of this cross-validation are stored in the knowledge graph associated with the online monitoring instruments.

9. The method for determining the accuracy of water pollution monitoring data based on multivariate statistics according to claim 8, characterized in that, When the tampering confidence score reaches the first confidence threshold and the cross-validation consistency of the multiple online monitoring instruments indicates that at least two online monitoring instruments are simultaneously abnormal, the data points marked as normal in the most recent maintenance cycle of the current online monitoring instrument are extracted from the knowledge graph associated with the online monitoring instruments, and the source-end calculation rules are self-optimized and updated. The specific process of the self-optimization update is as follows: using the incremental Bayesian linear regression method, the current slope coefficient and intercept coefficient of the source calculation rule are used as the mean of the prior distribution, the normal data points extracted from the knowledge graph associated with the online monitoring instrument are used as observation samples, the posterior distribution is calculated, and the updated slope coefficient and intercept coefficient are obtained based on the posterior distribution. After the update is completed, the updated slope coefficient and intercept coefficient are stored as new source-end calculation rules, and the baseline parameters of the corresponding online monitoring instruments in the online monitoring instrument association knowledge graph are updated synchronously.

Images