GNSS (Global Navigation Satellite System) monitoring data abnormal value detection method and device based on CEEMD (Continuous Empirical Mode Decomposition) and local abnormal factors, electronic equipment and storage medium

Abstract

The invention discloses a GNSS (Global Navigation Satellite System) monitoring data abnormal value detection method based on CEEMD (Compensation Empirical Mode Decomposition) and a local abnormal factor. The method comprises the following steps: carrying out preprocessing, gross error removal, interpolation processing and mirror image expansion on a three-dimensional direction coordinate time sequence of a GNSS monitoring terminal; cEEMD decomposition is carried out on the extended time sequence; determining a boundary layer of the noiseless signal component and the noise signal component according to a cross validation model; reconstructing a time sequence according to the boundary layer; determining an abnormal point according to a K-nearest neighbor algorithm and a local abnormal factor algorithm; abnormal points in the residual term are removed, interpolation processing is carried out, and an accurate residual term is obtained; and reconstructing a GNSS coordinate time sequence according to the precise residual term and the trend term. The method is especially suitable for processing a data set with high complexity and great local density difference, the accuracy of a processing result is improved so as to solve the problem of accurate detection and elimination of the gross error of the GNSS monitoring data, and the invention further discloses a device for realizing the method, electronic equipment and a computer readable storage medium.

Description

Technical Field

[0001] The present invention relates to the technical field of GNSS data processing, and particularly to a method, device, electronic device and storage medium for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor. Background Art

[0002] In recent years, with the continuous development and popularization of the Global Navigation Satellite System (GNSS) technology, low-cost GNSS receivers have received extensive attention and use from users due to their low cost and low power consumption. However, low-cost receivers are vulnerable to the influence of observation environment errors. Due to the complexity and uncertainty of the deformation monitoring working environment, GNSS signals are easily strongly interfered, resulting in inevitable gross errors in the deformation monitoring time series, and the gross errors show the characteristic of a large variation range. These large-scale gross errors often lead to the inability to accurately reflect the actual deformation of the deformed body or the actual change of the water level height, and will have an inestimable impact on its stability analysis and evaluation and even disaster warning and emergency response.

[0003] The detection and elimination of gross errors in GNSS sequences are the key prerequisites for GNSS coordinate time series processing and a hot issue in the processing of measurement data at home and abroad. At present, this technology can be divided into two major types: statistical gross error detection methods and non-statistical gross error detection methods. The statistical method constructs a gross error discrimination statistic based on the error distribution characteristics of the observed data and detects gross errors through hypothesis testing. The three most commonly used statistics are three times the mean square error, the interquartile range (IQR), and the median absolute deviation (MAD). The method identifies outliers by judging whether the residual value between the GNSS monitoring point data and the sample mean is greater than three times the mean square error, and it has a significant effect when the prior model is accurate and the gross error data is small, but it strongly depends on the normal distribution of the monitoring data, is overly sensitive to outliers, and has insufficient robustness. The IQR method, as a robust statistical method, replaces the mean and standard deviation with the median and the standard interquartile range to construct a robust Z-ratio score statistic to detect gross errors, and has strong adaptability to data deviating from the model. However, when the dispersion is large, the sensitivity to small-scale gross errors decreases due to the expansion of the interquartile range. In addition, the analysis methods based on statistics are generally sensitive to the trend term of the time series. It is necessary to fit the trend by least squares (LSQ) and then analyze the residuals. However, LSQ itself has weak anti-gross error ability, its residuals are easily contaminated, and the fitting effect depends on the model matching degree, ultimately affecting the gross error detection effect. Therefore, existing methods need to optimize the model assumptions and anti-error strategies in combination with data characteristics to balance sensitivity and robustness.

[0004] In recent years, in response to the demand for error suppression of GNSS coordinate time series, spatial filtering techniques represented by empirical mode decomposition (EMD) and its improved algorithms have been widely used. As a non-parametric adaptive signal processing method, EMD extracts the inherent time-scale characteristics of the signal, decomposes the original time series into several intrinsic mode functions (IMFs) from high frequency to low frequency and a residual term, and realizes denoising by removing the first few IMFs containing high-frequency noise. However, its decomposition process is vulnerable to endpoint effects and mode mixing interference. To this end, ensemble empirical mode decomposition (EEMD) introduces Gaussian white noise perturbation and multiple averaging strategies to suppress mode mixing through the statistical characteristics of noise. The further optimized complementary ensemble empirical mode decomposition (CEEMD) uses the superposition method of positive and negative white noise pairs in EEMD, eliminates residual noise while retaining the advantages of EEMD, significantly improves the GNSS signal reconstruction accuracy and calculation efficiency, and alleviates the endpoint effect. However, in the environment of water level monitoring, the water level time series data is complex and variable, and the stability and accuracy of anomaly detection cannot meet the requirements. Summary of the Invention

[0005] To solve the above problems existing in the prior art, the present invention provides a method, device, electronic device and storage medium for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor. The technical problems to be solved by the present invention are realized through the following technical solutions:

[0006] The first aspect of the embodiment of the present invention provides a method for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor, including the following steps:

[0007] Preprocess the three-dimensional direction coordinate time series of the GNSS monitoring terminal to obtain a preprocessed time series;

[0008] Remove gross errors, perform interpolation processing and mirror extension on the preprocessed time series to obtain an extended time series;

[0009] Perform CEEMD decomposition on the extended time series to obtain a decomposition representation of the extended time series;

[0010] Divide the extended time series into an odd sequence and an even sequence, use the even sequence as a verification sample, and the odd sequence as a filtering sample, and determine the boundary layer between the noise-free signal component and the noise signal component according to the cross-validation model;

[0011] Delete the data of the mirror extension part in the decomposed representation of the extended signal to obtain the decomposed representation of the original signal. In the decomposed representation of the original signal, reconstruct the trend term of the time series according to the demarcation layer, and reconstruct the residual term of the time series for the remaining IMF components;

[0012] Determine the outliers in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm;

[0013] Remove the outliers in the residual term and perform interpolation processing to obtain an accurate residual term;

[0014] Reconstruct the GNSS coordinate time series according to the accurate residual term and the trend term.

[0015] In one embodiment of the present invention, the preprocessing includes: coordinate transformation and mean removal processing.

[0016] In one embodiment of the present invention, the obtaining of the extended time series by performing gross error removal, interpolation processing, and mirror extension on the preprocessed time series includes:

[0017] Remove the data points within the range of μ±3σ in the preprocessed time series according to the 3σ criterion to obtain an intermediate time series;

[0018] Perform interpolation processing on the missing values of the time series of the intermediate time series using the KKF method, and then perform mirror extension to obtain an extended time series.

[0019] In one embodiment of the present invention, the obtaining of the decomposed representation of the extended time series by performing CEEMD decomposition on the extended time series includes:

[0020] Generate interference signals: Add multiple pairs of complementary Gaussian white noise pairs to the extended time series to obtain corresponding multiple groups of interference time series;

[0021] Perform EMD decomposition on multiple groups of the interference time series to obtain the components of the current layer;

[0022] Determine the residual of the current layer according to the components of the current layer, and return the residual of the current layer as the extended time series to be decomposed to the step of generating interference signals until the components of the first preset number of layers and the residual of the first preset number of layers are obtained;

[0023] Determine the decomposed representation of the extended signal according to the components of each layer and the residual of the first preset number of layers.

[0024] In one embodiment of the present invention, the dividing the extended time series into an odd sequence and an even sequence, using the even sequence as a verification sample, and using the odd sequence as a filtering sample, and determining the demarcation layer between the noise-free signal component and the noise signal component according to the cross-validation model includes:

[0025] Represent the extended time series as a composite function with time as the independent variable and displacement as the dependent variable, and divide the composite function into an odd sequence and an even sequence;

[0026] Randomly sample the even sequence to obtain verification samples;

[0027] Perform CEEMD decomposition on the odd sequence to obtain IMF components of a second preset layer and a residual vector;

[0028] Take the sum of the IMF components of the second preset layer accumulated in sequence as a filtering value, and calculate the variance of the verification samples with respect to the filtering value accumulated each time;

[0029] Take the layer corresponding to the minimum variance as the boundary layer between the noise-free signal component and the noise signal component.

[0030] In an embodiment of the present invention, deleting the data of the mirror extension part in the decomposition representation of the extended time series to obtain the decomposition representation of the original signal, and in the decomposition representation of the original signal, reconstructing the trend term of the time series according to the boundary layer, and reconstructing the residual term of the time series with the remaining IMF components, includes:

[0031] Delete the data of the mirror extension part in the decomposition representation of the extended time series to obtain the decomposition representation of the original signal;

[0032] In the decomposition representation of the original signal, reconstruct the trend term of the time series according to the low-frequency IMF components from the boundary layer to the last layer of the second preset layer, and reconstruct the residual term of the time series with the remaining IMF components.

[0033] In an embodiment of the present invention, determining the outlier points in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm includes:

[0034] Determine the k-nearest neighbors and k-distance of any sample point in the residual term according to the K-nearest neighbor algorithm;

[0035] Determine the local reachability distance according to the k-distance of the arbitrary sample point;

[0036] Determine the local reachability density according to the local reachability distance;

[0037] Determine the local outlier factor according to the local reachability density;

[0038] Determine the outlier points according to the local outlier factor and the factor threshold.

[0039] The second aspect of the embodiments of the present invention provides a GNSS monitoring data outlier detection device based on CEEMD and local outlier factor, including:

[0040] A first preprocessing module, configured to preprocess the three-dimensional direction coordinate time series of the GNSS monitoring terminal to obtain a preprocessed time series;

[0041] A second preprocessing module, configured to perform gross error removal, interpolation processing, and mirror extension on the preprocessed time series to obtain an extended time series;

[0042] A decomposition module, configured to perform CEEMD decomposition on the extended time series to obtain a decomposition representation of the extended time series;

[0043] A cross-validation module, configured to divide the extended time series into an odd sequence and an even sequence, use the even sequence as a validation sample, and the odd sequence as a filtering sample, and determine the boundary layer between the noise-free signal component and the noise signal component according to the cross-validation model;

[0044] A reconstruction module, configured to delete the data of the mirror extension part in the decomposition representation of the extended signal to obtain a decomposition representation of the original signal. In the decomposition representation of the original signal, reconstruct the trend term of the time series according to the boundary layer, and reconstruct the residual term of the time series for the remaining IMF components;

[0045] A calculation module, configured to determine the outliers in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm;

[0046] A denoising module, configured to remove the outliers in the residual term and perform interpolation processing to obtain an accurate residual term;

[0047] An accurate reconstruction module, configured to reconstruct the GNSS coordinate time series according to the accurate residual term and the trend term.

[0048] The third aspect of the embodiments of the present invention provides an electronic device, including a memory, a processor, and a computer program stored on the memory and executable on the processor. When the processor executes the program, it implements a GNSS monitoring data outlier detection method provided in the first aspect of the embodiments of the present invention based on CEEMD and local outlier factor.

[0049] The fourth aspect of the embodiments of the present invention provides a computer-readable storage medium, on which a computer program is stored. When the computer program is executed by a processor, it implements a GNSS monitoring data outlier detection method provided in the first aspect of the embodiments of the present invention based on CEEMD and local outlier factor.

[0050] Advantages of the present invention:

[0051] In terms of trend extraction and noise separation, the present invention abandons the traditional gross error detection algorithms based on statistical analysis and least squares. The obtained trend is more accurate and clear, and the separated noise is also purer and free of impurities, greatly improving the accuracy of data processing.

[0052] Regarding the boundary determination between the signal component and the noise component, the present invention uses the cross-validation method, which has better robustness and higher reliability, providing reliable data support for subsequent accurate analysis, and using mirror extension to have a certain inhibitory effect on the endpoint effect of modal decomposition filtering.

[0053] The present invention introduces the LOF method based on local density, which has unique advantages for skewed data with uneven distribution. Whether it is an isolated outlier or a density outlier, it can be accurately captured, especially suitable for processing datasets with high complexity and large differences in local density.

[0054] Moreover, estimating the noise by means of IMF components further strengthens the anti-trend performance of the LOF method, enabling it to still play a stable role in a complex and changing data environment, providing strong technical support for the data processing requirements in various fields.

[0055] Other features and advantages of the present invention will be described in the subsequent specification, and part of them will become obvious from the specification or be understood by implementing the present invention. The objectives and other advantages of the present invention can be achieved and obtained through the structures specifically pointed out in the written specification, claims, and drawings.

[0056] The technical solutions of the present invention will be further described in detail below through the drawings and embodiments. Description of the Drawings

[0057] The drawings are used to provide further understanding of the present invention and constitute a part of the specification. They are used together with the embodiments of the present invention to explain the present invention and do not constitute a limitation to the present invention. In the drawings:

[0058] Figure 1 It is a schematic flow chart of a method for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor provided by an embodiment of the present invention;

[0059] Figure 2 It is a schematic diagram of the reconstructed GNSS coordinate time series after the method for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor provided by an embodiment of the present invention is applied to the processing of GNSS buoy terminal water level monitoring data;

[0060] Figure 3A schematic diagram of a GNSS monitoring data outlier detection device based on CEEMD and local anomaly factors provided by an embodiment of the present invention. DETAILED DESCRIPTION

[0061] The present invention will be further described in detail below with reference to specific examples, but the embodiments of the present invention are not limited thereto.

[0062] like Figure 1 As shown, a first aspect of an embodiment of the present invention provides a method for detecting outliers in GNSS monitoring data based on CEEMD and a local anomaly factor, comprising the following steps:

[0063] Step 11: preprocess the three-dimensional direction coordinate time series of the GNSS monitoring terminal to obtain a preprocessed time series.

[0064] Specifically, the three-dimensional direction coordinate time series of the GNSS monitoring terminal is obtained, and the preprocessed time series of the three-dimensional direction or vertical component is obtained through coordinate conversion and de-meaning processing.

[0065] Here, the three-dimensional direction coordinate time series is the data generated when the GNSS monitoring terminal performs water level monitoring.

[0066] Step 12: performing gross error removal, interpolation processing, and mirror expansion on the preprocessed time series to obtain an expanded time series.

[0067] The specific steps of step 12 include steps 121 to 124:

[0068] In step 121 , preliminary gross error detection is performed on the preprocessed time series according to the 3σ criterion, and GNSS data points within the range of ±3σ of the mean μ are eliminated to obtain an intermediate time series.

[0069] Step 122: Use Kriging Kalman Filter (KKF) method to interpolate the missing values of the intermediate time series to obtain an interpolated time series.

[0070] Step 123: Extend the beginning and end of the interpolated time series by one-quarter of the total length of the series using a mirror extension method. Mirror extension can suppress the influence of endpoint effects on modal decomposition accuracy during the subsequent CEEMD filtering process.

[0071] Step 13: Perform CEEMD decomposition on the extended time series to obtain a decomposition representation of the extended signal.

[0072] Step 13 includes steps 131 to 133:

[0073] Step 131, generate interference signals: Add multiple pairs of complementary Gaussian white noise pairs to the extended time series x(t) to obtain corresponding multiple groups of interference time series.

[0074] In this step, the added complementary Gaussian white noise pairs are positive noise ω j (t) and negative noise -ω j (t), where j = 1, 2, …, q represents the number of groups of noise pairs, and t represents time. Any group of obtained interference time series is:

[0075]

[0076] Step 132, perform EMD decomposition on multiple groups of interference time series to obtain the components of the current layer.

[0077] Perform EMD decomposition on the interference time series after the first addition of interference. A group of IMF components obtained is denoted as and Take the mean of this group of components to effectively cancel the error introduced by the noise and obtain the final IMF component IMF1 of the first layer:

[0078]

[0079] The residual corresponding to the first layer decomposition is: Res1(t) = x(t) - IMF1.

[0080] Step 133, then regard the residual Res1(t) after the first layer decomposition as the signal to be decomposed, process it according to Step 131, add multiple pairs of complementary Gaussian white noise pairs to Res1(t), and then perform decomposition according to the method of Step 132 to calculate and obtain the second layer IMF2: Then the residual after the second layer decomposition is: Res2(t) = Res1(t) - IMF2. Continue to regard Res2(t) as the signal to be decomposed, process it according to the methods of Step 131 and Step 132, and continue to decompose the obtained results in a loop until the components IMF n and the residual Res n (t) of the nth layer required are obtained.

[0081] At this time, the decomposition representation of the signal is obtained: That is, the decomposition representation of the extended time series. i = 1, 2, …, n represents the serial number of the decomposition layer, and the decomposition layer is determined according to actual requirements.

[0082] Step 14: Divide the extended time series into an odd sequence and an even sequence. Take the even sequence as the validation sample and the odd sequence as the filtering sample, and determine the boundary layer between the noise-free signal component and the noise signal component according to the cross-validation model. In this step, the cross-validation model is used to determine the boundary IMF between the relatively pure signal component and the pure noise component.

[0083] The specific steps of Step 14 include Step 141 - Step 145:

[0084] Step 141: Represent the extended time series as a composite function with time as the independent variable and displacement as the dependent variable, and divide the composite function into an odd sequence and an even sequence.

[0085] The extended time series is represented as a function with time as the independent variable and displacement as the dependent variable, i.e., x(t(i)), where i = 1, 2, …, m, representing the serial number of the sample, and the total number of samples is m. Divide it into an odd sequence x1(t(i ′ )) and an even sequence x2(t(i ″ ))). Among them, i ′ = 2n - 1, n = 1, 2, …, m1, i ″ = 2n, n = 1, 2, …, m2.

[0086] When m is even, m1 = m2 = m / 2; when m is odd, m1 = (m + 1) / 2, m2 = (m - 1) / 2.

[0087] Step 142: Perform CEEMD decomposition on the odd sequence to obtain the IMF components of the second preset layer and the residual vector. Take the odd sequence samples as the filtering samples, and perform CEEMD decomposition on the odd sequence x1(t(i ′ )) to obtain the IMF components of the second preset layer with b layers and the residual vector R. a = 1, 2, …, b, representing the serial number of the decomposition layer of the odd sequence, and the total number of layers is b.

[0088] Step 143: Perform random sampling on the even sequence to obtain the validation samples. Considering the statistical characteristics and computational amount of random division, select M different random divisions r s (s = 1, 2, … M) each time of calculation. r represents a certain random division of the even sequence. Each random division is a randomly selected sample size m test = 0.2m, that is, 20% of the total number of samples m, to obtain the validation samples x test (t(i ″ ))).

[0089] Step 144: Take the sum of the IMF components of the second preset layer accumulated in sequence as the filtering value, and calculate the variance of the validation samples with respect to each accumulated filtering value.

[0090] In this step, perform b accumulations. The first accumulation is to use the sum of the IMF components from the first layer to the b-th layer as the first filtering value f. ′ , and at the same time, calculate the individual variance of each of the M verification samples selected this time with respect to the first filtering value f. ′ The individual variance of any one verification sample with respect to the filtering value is : Obtain M individual variances, calculate the average value var of the M individual variances. a , This average value var a is used as the variance of the first layer.

[0091] The second accumulation is to use the sum of the IMF components from the second layer to the b-th layer as the second filtering value. At the same time, calculate the individual variance of each of the M verification samples selected this time (each time pre-selecting verification samples again, randomly selecting M different verification samples each time, that is, randomly selecting different M verification samples each time) with respect to the second filtering value, obtain M individual variances, calculate the average value of the M individual variances, and this average value is used as the variance of the second layer.

[0092] Then perform the third accumulation and calculate the variance of the third layer, and so on for accumulative calculation until the b-th accumulation is completed and the variance of the b-th layer is obtained.

[0093] Step 145: Use the layer corresponding to the minimum variance as the boundary layer between the noise-free signal component and the noise signal component.

[0094] In this step, among the b variances calculated in step 145, the a corresponding to the minimum variance is the boundary layer. a = 1, 2,..., b, representing the serial number of the decomposition layer of the odd sequence, and the total number of layers is b. At this time, the sum of the low-frequency IMF components from the a-th layer to the b-th layer is the layer where the CEEMD-filtered signal is located.

[0095] Step 15: Delete the data of the mirror extension part in the decomposition representation of the extended time series to obtain the decomposition representation of the original signal. In the decomposition representation of the original signal, reconstruct the trend term of the time series according to the boundary layer, and reconstruct the residual term of the time series with the remaining IMF components.

[0096] In this step, specifically, delete the data of that mirror extension part in the decomposition representation of the extended time series to obtain the decomposition representation of the original signal. Then, in the decomposition representation of the original signal, reconstruct the trend term of the time series according to the low-frequency IMF components from the boundary layer to the second preset layer and the residual of the last layer, and reconstruct the residual term of the time series with the remaining IMF components.

[0097] The low-frequency IMF components are used to accumulate and reconstruct the GNSS sequence to obtain its signal layer, that is, the trend term:

[0098]

[0099] The remaining components are reconstructed into the noise layer, that is, the residual term:

[0100]

[0101] For the boundary determination of the signal component and the noise component, the present invention adopts the cross-validation method. Compared with the traditional indicators that rely on the correlation coefficient, energy density, etc., this method has better robustness and higher reliability, laying a solid foundation for subsequent accurate analysis.

[0102] Step 16, determine the outliers in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm. The specific steps of Step 16 include Step 161 - Step 165:

[0103] Step 161, determine the k-nearest neighbor points and the k-distance of any sample point in the residual term according to the K-nearest neighbor algorithm.

[0104] In this step, for any sample point p in the residual term of Step 15, select the k sample points closest to it, that is, the k-nearest neighbor points N k (p). Among them, for any two sample points p and v, the distance metric method uses the Euclidean distance: d(p,v) = ‖p - v‖.

[0105] Then, calculate the k-distance k-dist(p). For point p, its k-distance is defined as the maximum value of the distance from point p to its k-th nearest neighbor point o, that is:

[0106] k-dist(p) = max{d(p,o)|o ∈ N k (p)}

[0107] Step 162, determine the local reachability distance according to the k-distance of any sample point.

[0108] In this step, calculate the local reachability distance reach-dist(p,o). The reachability distance of point p and its neighbor o ∈ N k (p) is defined as the larger value of the k-distance of point o and the distance between point p and point o.

[0109] reach-dist(p,o) = max{k-dist(o), d(p,o)}

[0110] Step 163, determine the local reachability density according to the local reachability distance.

[0111] Calculate the Local Reachability Density (LRD). The local reachability density is an indicator used to measure the density degree of point p relative to its neighborhood. Specifically, it is characterized as the reciprocal of the average local reachability distance of the data points within the k-distance neighborhood of point p. Its physical meaning is to indirectly reflect whether the density of point p is lower than that of its neighborhood points, that is, to describe the relative density degree between this point and its neighborhood points. The calculation formula for the local reachability density is as follows:

[0112]

[0113] It should be noted that LRD k (p) reflects the degree of deviation of the density of point p from that of its neighbor points. The larger its value, the higher the density around point p.

[0114] Step 164: Determine the local outlier factor according to the local reachability density.

[0115] In this step, calculate the local outlier factor LOF k (p). LOF is defined as the average value of the ratio of the local reachability density of point p to the local reachability densities of the data points within its k-distance neighborhood.

[0116]

[0117] Step 165: Determine the outlier points (i.e., outliers) according to the local outlier factor and the factor threshold 1.

[0118] LOF k (p) is used to characterize the degree of deviation of the data point p. When LOF k (p) ≈ 1, the density of p is similar to that of its neighbors, indicating that p is a normal point; LOF k (p) > 1, the density of p is significantly lower than that of its neighbors, indicating that p is a potential outlier. The larger LOF k (p) is, the more likely p is an outlier.

[0119] Step 17: Remove the outlier points (i.e., outliers) in the residual term and perform interpolation processing to obtain the accurate residual term.

[0120] After determining the outlier points according to the above steps, remove the outlier points. To ensure the continuity and usability of the data, the KKF method is used again for missing value interpolation.

[0121] Step 18: Reconstruct the GNSS coordinate time series according to the accurate residual term and the trend term.

[0122] In this embodiment, as Figure 2As shown in the figure, it is the GNSS coordinate time series reconstructed from the precise residual term and the trend term of the water level monitoring data of the GNSS buoy terminal processed according to the above steps. The plus signs in the figure indicate the eliminated gross errors, and the dots indicate the time series after removing the gross errors. It should be noted that for water level deformation monitoring, the water level fluctuates greatly and the time series is extremely unstable. Through the method of the present invention, the trend can be more accurately identified and the noise can be reduced, the outliers can be more accurately detected, and finally more accurate time series data can be obtained, providing more accurate data support for water level deformation analysis.

[0123] As Figure 3 shown, the second aspect of the embodiment of the present invention provides an outlier detection device for GNSS monitoring data based on CEEMD and local outlier factor, including:

[0124] The first preprocessing module 21 is used to preprocess the three-dimensional direction coordinate time series of the GNSS monitoring terminal to obtain a preprocessed time series;

[0125] The second preprocessing module 22 is used to remove gross errors, perform interpolation processing and mirror extension on the preprocessed time series to obtain an extended time series;

[0126] The decomposition module 23 is used to perform CEEMD decomposition on the extended time series to obtain a decomposition representation of the extended time series;

[0127] The cross-validation module 24 is used to divide the extended time series into an odd sequence and an even sequence, use the even sequence as the verification sample, use the odd sequence as the filtering sample, and determine the boundary layer between the noise-free signal component and the noise signal component according to the cross-validation model;

[0128] The reconstruction module 25 is used to delete the data of the mirror extension part in the decomposition representation of the extended signal to obtain the decomposition representation of the original signal. In the decomposition representation of the original signal, the trend term of the time series is reconstructed according to the boundary layer, and the residual terms of the time series are reconstructed from the remaining IMF components;

[0129] The calculation module 26 is used to determine the outliers in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm;

[0130] The denoising module 27 is used to remove the outliers in the residual term and perform interpolation processing to obtain a precise residual term;

[0131] The precise reconstruction module 28 is used to reconstruct the GNSS coordinate time series according to the precise residual term and the trend term.

[0132] In an embodiment of the present invention, the preprocessing includes: coordinate transformation and mean removal processing.

[0133] In one embodiment of the present invention, the preprocessed time series is subjected to gross error removal, interpolation processing, and mirror extension to obtain an extended time series, including:

[0134] According to the 3σ criterion, data points within the range of μ±3σ in the preprocessed time series are removed to obtain an intermediate time series;

[0135] The missing values in the intermediate time series are interpolated using the KKF method and then mirror extended to obtain an extended time series.

[0136] In one embodiment of the present invention, the extended time series is subjected to CEEMD decomposition to obtain the decomposition representation of the extended time series, including:

[0137] Generate interference signals: Add multiple pairs of complementary Gaussian white noise pairs to the extended time series to obtain corresponding multiple groups of interference time series;

[0138] Perform EMD decomposition on the multiple groups of interference time series to obtain the current layer components;

[0139] Determine the current layer residual based on the current layer components, and return the current layer residual as the extended time series to be decomposed to the step of generating interference signals until the components of the first preset number of layers and the residuals of the first preset number of layers are obtained;

[0140] Determine the decomposition representation of the extended signal according to the components of each layer and the residuals of the first preset number of layers.

[0141] In one embodiment of the present invention, the extended time series is divided into an odd sequence and an even sequence. The even sequence is used as the verification sample, and the odd sequence is used as the filtering sample. The boundary layer between the noise-free signal component and the noise signal component is determined according to the cross-validation model, including:

[0142] Represent the extended time series as a composite function with time as the independent variable and displacement as the dependent variable, and divide the composite function into an odd sequence and an even sequence;

[0143] Randomly sample the even sequence to obtain the verification sample;

[0144] Perform CEEMD decomposition on the odd sequence to obtain the IMF components of the second preset number of layers and the residual vector;

[0145] Take the sum of the IMF components of the second preset number of layers accumulated in sequence as the filtering value, and calculate the variance of the verification sample with respect to each accumulated filtering value;

[0146] Take the layer corresponding to the minimum variance as the boundary layer between the noise-free signal component and the noise signal component.

[0147] In one embodiment of the present invention, data of the mirror extension part is deleted in the decomposed representation of the extended time series to obtain the decomposed representation of the original signal. In the decomposed representation of the original signal, the trend term of the time series is reconstructed according to the demarcation layer, and the residual term of the time series is reconstructed by the remaining IMF components, including:

[0148] Data of the mirror extension part is deleted in the decomposed representation of the extended time series to obtain the decomposed representation of the original signal;

[0149] In the decomposed representation of the original signal, the trend term of the time series is reconstructed according to the low-frequency IMF components from the demarcation layer to the last layer of the second preset number of layers, and the residual term of the time series is reconstructed by the remaining IMF components.

[0150] In one embodiment of the present invention, the outlier points in the residual term are determined according to the K-nearest neighbor algorithm and the local outlier factor algorithm, including:

[0151] The k-nearest neighbors and k-distance of any sample point in the residual term are determined according to the K-nearest neighbor algorithm;

[0152] The local reachability distance is determined according to the k-distance of any sample point;

[0153] The local reachability density is determined according to the local reachability distance;

[0154] The local outlier factor is determined according to the local reachability density;

[0155] The outlier points are determined according to the local outlier factor and the factor threshold.

[0156] The third aspect of the embodiments of the present invention provides an electronic device, including a memory, a processor, and a computer program stored on the memory and executable on the processor. When the processor executes the program, it implements the method for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor provided by the embodiments of the present invention.

[0157] The fourth aspect of the embodiments of the present invention further provides a computer-readable storage medium, on which a computer program is stored. When the computer program is executed by a processor, it implements the steps of the method for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor provided by the embodiments of the present invention.

[0158] Among them, the memory may include a random access memory (RAM), and may also include a non-volatile memory (NVM), such as at least one disk memory. Optionally, the memory may also be at least one storage device located far from the aforementioned processor.

[0159] The above-mentioned processor may be a general-purpose processor, including a Central Processing Unit (CPU), a Network Processor (NP), etc.; it may also be a Digital Signal Processor (DSP), an Application Specific Integrated Circuit (ASIC), a Field-Programmable Gate Array (FPGA), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware devices.

[0160] The method provided by the embodiments of the present invention can be applied to an electronic device. Specifically, the electronic device may be: a desktop computer, a portable computer, a smart mobile terminal, a server, etc. There is no limitation here. Any electronic device that can implement the present invention belongs to the protection scope of the present invention.

[0161] For the device / electronic device embodiments, since they are basically similar to the method embodiments, the description is relatively simple. For the relevant parts, please refer to the partial description of the method embodiments.

[0162] The present invention is described with reference to the flowcharts and / or block diagrams of methods, devices (apparatus), and computer program products according to the embodiments of the present invention. It should be understood that each flow and / or block in the flowchart and / or block diagram, as well as the combination of flows and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to the processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing devices to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing devices generate a device for implementing the functions specified in Figure 1 one flow or multiple flows and / or blocks Figure 1 one block or multiple blocks.

[0163] These computer program instructions can also be stored in a computer-readable memory that can direct a computer or other programmable data processing device to work in a specific manner, so that the instructions stored in the computer-readable memory generate a manufactured article including an instruction device, and the instruction device implements the functions specified in Figure 1 one flow or multiple flows and / or blocks Figure 1 one block or multiple blocks.

[0164] These computer program instructions can also be loaded onto a computer or other programmable data processing apparatus, so that a series of operation steps are executed on the computer or other programmable apparatus to produce a computer-implemented process, thereby providing instructions for implementing the functions specified in one process or a plurality of processes and / or blocks Figure 1 one process or a plurality of processes and / or blocks Figure 1 steps for the functions specified in one block or a plurality of blocks.

[0165] Obviously, those skilled in the art can make various changes and modifications to the present invention without departing from the spirit and scope of the present invention. Thus, if these modifications and variations of the present invention fall within the scope of the claims of the present invention and their equivalent technologies, the present invention is also intended to include these changes and modifications.

Claims

1. A method for detecting outliers in GNSS monitoring data based on CEEMD and local outlier factor, characterized in that, It includes the following steps: Preprocess the three-dimensional direction coordinate time series of the GNSS monitoring terminal to obtain a preprocessed time series; Perform gross error removal, interpolation processing, and mirror extension on the preprocessed time series to obtain an extended time series; Perform CEEMD decomposition on the extended time series to obtain a decomposition representation of the extended time series; Divide the extended time series into an odd sequence and an even sequence, use the even sequence as a verification sample, and the odd sequence as a filtering sample, and determine the boundary layer between the noise-free signal component and the noise signal component according to the cross-validation model; Delete the data of the mirror extension part in the decomposition representation of the extended signal to obtain a decomposition representation of the original signal. In the decomposition representation of the original signal, reconstruct the trend term of the time series according to the boundary layer, and the remaining IMF components reconstruct the residual term of the time series; Determine the outliers in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm; Remove the outliers in the residual term and perform interpolation processing to obtain an accurate residual term; Reconstruct the GNSS coordinate time series according to the accurate residual term and the trend term.

2. The method according to claim 1, characterized in that, The preprocessing includes: coordinate transformation and mean removal processing.

3. The method according to claim 1, wherein The performing gross error removal, interpolation processing, and mirror extension on the preprocessed time series to obtain an extended time series includes: Remove the data points within the range of μ±3σ in the preprocessed time series according to the 3σ criterion to obtain an intermediate time series; Use the KKF method to perform interpolation processing on the missing values of the time series for the intermediate time series, and then perform mirror extension to obtain an extended time series.

4. The method according to claim 1, characterized in that The performing CEEMD decomposition on the extended time series to obtain a decomposition representation of the extended time series includes: Generate interference signals: Add multiple pairs of complementary Gaussian white noise pairs to the extended time series to obtain corresponding multiple groups of interference time series; Perform EMD decomposition on multiple groups of the interference time series to obtain the current layer components; Determine the current layer residual according to the current layer components, and return the current layer residual as the extended time series to be decomposed to the step of generating interference signals until the first preset number of layer components and the first preset number of layer residuals are obtained; Determine the decomposition representation of the extended signal according to each layer component and the first preset number of layer residuals.

5. The method according to claim 1, wherein The dividing the extended time series into an odd sequence and an even sequence, using the even sequence as a verification sample, and the odd sequence as a filtering sample, and determining the boundary layer between the noise-free signal component and the noise signal component according to the cross-validation model includes: Represent the extended time series as a composite function with time as the independent variable and displacement as the dependent variable, and divide the composite function into an odd sequence and an even sequence; Perform random sampling on the even sequence to obtain a verification sample; Perform CEEMD decomposition on the odd sequence to obtain the IMF components of the second preset number of layers and the residual vector; Use the sum of the IMF components of the second preset number of layers accumulated in sequence as the filtering value, and calculate the variance of the verification sample for each accumulated filtering value; Use the layer corresponding to the minimum variance as the boundary layer between the noise-free signal component and the noise signal component.

6. The method according to claim 5, wherein Deleting the data of the mirror extension part in the decomposition representation of the extended time series to obtain the decomposition representation of the original signal. In the decomposition representation of the original signal, reconstructing the trend term of the time series according to the boundary layer, and reconstructing the residual term of the time series with the remaining IMF components, including: Deleting the data of the mirror extension part in the decomposition representation of the extended time series to obtain the decomposition representation of the original signal; In the decomposition representation of the original signal, reconstructing the trend term of the time series according to the low-frequency IMF components from the boundary layer to the last layer of the second preset number of layers, and reconstructing the residual term of the time series with the remaining IMF components.

7. The method according to claim 1, wherein Determining the outliers in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm, including: Determining the k-nearest neighbors and the k-distance of any sample point in the residual term according to the K-nearest neighbor algorithm; Determining the local reachability distance according to the k-distance of the arbitrary sample point; Determining the local reachability density according to the local reachability distance; Determining the local outlier factor according to the local reachability density; Determining the outliers according to the local outlier factor and the factor threshold.

8. An outlier detection device for GNSS monitoring data based on CEEMD and local outlier factor, characterized in that, Including: A first preprocessing module for preprocessing the three-dimensional direction coordinate time series of the GNSS monitoring terminal to obtain a preprocessed time series; A second preprocessing module for removing gross errors, interpolation processing and mirror extension on the preprocessed time series to obtain an extended time series; A decomposition module for performing CEEMD decomposition on the extended time series to obtain the decomposition representation of the extended time series; A cross-validation module for dividing the extended time series into an odd sequence and an even sequence, using the even sequence as the validation sample and the odd sequence as the filtering sample, and determining the boundary layer between the noise-free signal component and the noise signal component according to the cross-validation model; A reconstruction module for deleting the data of the mirror extension part in the decomposition representation of the extended signal to obtain the decomposition representation of the original signal. In the decomposition representation of the original signal, reconstructing the trend term of the time series according to the boundary layer, and reconstructing the residual term of the time series with the remaining IMF components; A calculation module for determining the outliers in the residual term according to the K-nearest neighbor algorithm and the local outlier factor algorithm; A denoising module for removing the outliers in the residual term and performing interpolation processing to obtain an accurate residual term; An accurate reconstruction module for reconstructing the GNSS coordinate time series according to the accurate residual term and the trend term.

9. An electronic device, comprising a memory, a processor, and a computer program stored on the memory and executable on the processor, characterized in that, When the processor executes the program, it implements the GNSS monitoring data outlier detection method based on CEEMD and local outlier factor according to any one of claims 1 to 7.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed by the processor, it implements the GNSS monitoring data outlier detection method based on CEEMD and local outlier factor according to any one of claims 1 to 7.

Images