A method, device and equipment for diagnosing and measuring the overall variability of a hydrological sequence

Abstract

The present invention discloses a method, device and equipment for diagnosing and measuring the overall variability of a hydrological sequence. The original hydrological sequence is decomposed into multiple IMF sequences and a trend term sequence; it is judged whether there is a trend variation point in the trend term sequence. If so, the trend variability of the trend term sequence is calculated; the multiple IMF sequences are superimposed and recombined to obtain an IMFs sequence, and it is judged whether there is a mean variation point in the IMFs sequence. If so, the means of the subsequences on both sides of the mean variation point and the mean variability of the IMFs sequence are calculated; the IMFs sequence is subjected to mean normalization processing according to the means of the subsequences on both sides of the mean variation point, and it is judged whether there is a variance variation point in the mean-normalized IMFs sequence. If so, the variance variability of the mean-normalized IMFs sequence is calculated; the comprehensive variability of the original hydrological sequence is calculated according to the trend variability, mean variability and variance variability. The present invention can obtain complete variability information and fully reflect the specific change degree of the hydrological sequence.

Description

Technical Field

[0001] The present invention relates to the field of hydrological and hydraulic calculations, and more specifically, provides a method, device, and equipment for diagnosing and measuring the overall variability of hydrological sequences. Background Art

[0002] Under the influence of human production and life, the water in nature is constantly changing. In addition, changes in the natural environment can also cause mutations in hydrological elements. These changes are often reflected in the hydrological sequence through changes in sequence characteristic values, that is, hydrological variations, making the hydrological phenomenon usually show a non-stationary change trend in time series.

[0003] In order to diagnose variations in non-stationary time series, relevant scholars have proposed a series of research methods, such as the Pettitt mutation point test method, the heuristic segmentation algorithm (BG algorithm), the Bayesian estimation method, the STARS method, the iterative cumulative sum of squares algorithm (ICSS), and the SIC change point test theory. Although these methods can diagnose the variation type from a certain aspect and have certain advantages, the result information obtained by a single diagnostic method is often very limited and cannot obtain complete variation information. In addition, there is less research on the degree of hydrological variation. Although some scholars have proposed some variation degree calculation methods based on the correlation coefficient, most of them calculate the variation scale starting from a certain type of variation and cannot fully reflect the specific change degree of the hydrological sequence. Currently, there is no research method that can separately measure different variation types. Summary of the Invention

[0004] Aiming at the problems existing in the prior art, the present invention provides a method, device, and equipment for diagnosing and measuring the overall variability of hydrological sequences, which can obtain complete variation information and fully reflect the specific change degree of the hydrological sequence.

[0005] In order to solve the above technical problems, the present invention is realized through the following technical solutions:

[0006] A method for diagnosing and measuring the overall variability of hydrological sequences includes:

[0007] Obtain the original hydrological sequence, and decompose the original hydrological sequence into a plurality of IMF sequences and a trend item sequence;

[0008] Judge whether there is a trend variation point in the trend item sequence. If so, calculate the trend variability of the trend item sequence;

[0009] Stack and recombine a plurality of the IMF sequences to obtain an IMFs sequence, and judge whether there is a mean variation point in the IMFs sequence. If so, calculate the means of the subsequences on both sides of the mean variation point and the mean variability of the IMFs sequence;

[0010] Perform mean normalization on the IMFs sequence according to the means of the subsequences on both sides of the mean mutation point, and determine whether there is a variance mutation point in the mean-normalized IMFs sequence. If so, calculate the variance mutation degree of the mean-normalized IMFs sequence.

[0011] Calculate the comprehensive mutation degree of the original hydrological sequence according to the trend mutation degree, the mean mutation degree, and the variance mutation degree.

[0012] Further, the decomposition of the original hydrological sequence into multiple IMF sequences and a trend term sequence includes:

[0013] 1) Add a white noise signal ω(t) with the same time scale and conforming to the normal distribution to the original hydrological sequence x(t) to obtain a new hydrological sequence y(t) after noise reduction processing:

[0014] y(t) = x(t) + ω(t)

[0015] 2) Perform EMD decomposition on the new hydrological sequence y(t) to obtain multiple component sequences;

[0016] 3) Set the number of EMD decomposition times as j, repeat steps 1) - 2) j times, and average each of the component sequences obtained by cycling to obtain multiple IMF sequences and a trend term sequence.

[0017] 3. According to a method for diagnosing and measuring the overall mutation degree of a hydrological sequence as described in claim 1, wherein the determination of whether there is a trend mutation point in the trend term sequence includes:

[0018] Use the self-segmented linear fitting method to perform segmented trend line fitting on the trend term sequence;

[0019] According to the χ 2 test method, determine whether the segmented trend line satisfies significance. If it satisfies significance, there is a trend mutation point.

[0020] Further, the determination of whether there is a mean mutation point in the IMFs sequence includes:

[0021] 1) Set the IMFs sequence as M(t), the number of its sequence points as n, and i as a possible mutation point of M(t), where l ≤ i ≤ n - l, and l is the minimum truncation length of M(t);

[0022] 2) Calculate the combined standard deviation σ D (i), and the calculation formula is as follows:

[0023]

[0024] where: μ r (i) and μ l (i) are the means of the sequence segments on both sides of M(i), respectively, and σ r (i) and σ l (i) are the standard deviations of the sequence segments on both sides of M(i), respectively; n r and n l are the number of sequence points on the left and right sides of point i, respectively;

[0025] 3) Construct a test statistic T(i) to quantitatively represent the mean difference between the sequence segments on the left and right sides of point i. The formula is as follows:

[0026]

[0027] where: the value of T(i) represents the degree of mean variation. The larger its value, the greater the amplitude of the mean jump at point i. Calculate the corresponding statistic T(i) in turn within the value range of i to obtain the test statistic sequence T(t) corresponding to M(t);

[0028] 4) Set the maximum value of the test statistic sequence T(t) as T m , and its significance level P is expressed as:

[0029]

[0030] η = 4.19ln(n) - 11.54

[0031] where v = n - 2; δ = 0.4; S x (a, b) is the incomplete beta function;

[0032] 5) If P ≥ P0, where P0 is the significance critical value and its value range is 0.5 - 0.95, it indicates that the mean difference is significant. Taking point i as the segmentation point, M(t) is divided into two sub-sequences, and point i is the mean variation point;

[0033] 6) Perform the above operations on the two sub-sequences obtained in 5) respectively to check whether there are still mean variation points in the sub-sequences, and finally obtain c IMFs sub-sequences and (c - 1) mean variation points.

[0034] Furthermore, the mean normalization process of the IMFs sequence according to the means of the sub-sequences on both sides of the mean variation point includes:

[0035] 1) Let an IMFs sub-sequence before processing be m(t), t = 1, 2, 3, …, is the mean of m(t);

[0036] 2) The processed IMFs sub-sequence is m*(t), and its calculation formula is as follows:

[0037]

[0038] 3) Combine the processed multi-segment IMFs subsequences, and the mean-zeroed IMFs sequence is obtained.

[0039] Further, determining whether there is a variance mutation point in the mean-zeroed IMFs sequence includes:

[0040] 1) Define the entire mean-zeroed IMFs sequence as Y t , where t = 1, 2, 3,..., n, and its mean is 0, and the overall variance is Assume that there is a variance mutation point at r in the sequence, and construct the test statistic S max :

[0041]

[0042] In the formula, L(θ) is the maximum likelihood function; λ is the degree of freedom during calculation. When calculating for the entire sequence, λ = 1, and when calculating for the variance mutation point, λ = 2; n is the sequence length; is the variance of the subsequence before the r point, is the variance of the subsequence after the r point;

[0043] 2) For C(α) > 1, where α is the significance level, if S max ≥ C(α), then the variance mutation at the r point is significant, and the r point is a variance mutation point;

[0044] 3) Apply the bisection method to divide Y t into two subsequences at the variance mutation point, and then perform variance mutation tests on each subsequence according to the above steps to obtain all the variance mutation points of Y t .

[0045] Further, the trend mutation degree of the trend term sequence, the mean mutation degree of the IMFs sequence, and the variance mutation degree of the mean-zeroed IMFs sequence are calculated by using the correlation integral method, specifically as follows:

[0046] 1) Define the sequence for which the mutation degree is to be calculated as Z t , where t = 1, 2, 3,..., N, then the vertical distances between any two points in Z t form a new sequence D s :

[0047] D s = ||Z i - Z j ||

[0048] where \(s = 1, 2, \ldots, N(N - 1) / 2\); \(i = 1, 2, \ldots, N - 1\); \(j = i + 1, i + 2, \ldots, N\);

[0049] 2) Define the statistic \(I(N, R)\) to measure the non - stationarity of \(Z\) t :

[0050]

[0051] where \(N\) is the total number of points; \(R\) is the reference distance, determined by the statistical distribution law of \(D\) s ; \(H\) is the unit step function, defined as follows:

[0052]

[0053] 3) If there are change points in the test sequence, calculate the non - stationarity degrees \(I_1(N_1, R_1)\) and \(I_2(N_2, R_2)\) of the sequences before and after the change points respectively, and obtain the change degree \(U\) corresponding to the trend term sequence:

[0054] \(U=\vert I_1(N_1, R_1)-I_2(N_2, R_2)\vert\).

[0055] Furthermore, calculating the comprehensive change degree of the original hydrological sequence according to the trend change degree, the mean change degree and the variance change degree includes:

[0056] Adding the trend change degree, the mean change degree and the variance change degree to obtain the comprehensive change degree of the original hydrological sequence.

[0057] A device for diagnosing and measuring the overall change degree of a hydrological sequence includes:

[0058] A decomposition module, configured to obtain an original hydrological sequence and decompose the original hydrological sequence into a plurality of IMF sequences and a trend term sequence;

[0059] A trend change degree calculation module, configured to determine whether there is a trend change point in the trend term sequence, and if so, calculate the trend change degree of the trend term sequence;

[0060] A mean change degree calculation module, configured to stack and recombine a plurality of the IMF sequences to obtain an IMFs sequence, determine whether there is a mean change point in the IMFs sequence, and if so, calculate the means of the subsequences on both sides of the mean change point and the mean change degree of the IMFs sequence;

[0061] A variance change degree calculation module, configured to perform mean zeroing on the IMFs sequence according to the means of the subsequences on both sides of the mean change point, determine whether there is a variance change point in the mean - zeroed IMFs sequence, and if so, calculate the variance change degree of the mean - zeroed IMFs sequence;

[0062] A comprehensive variability calculation module, configured to calculate the comprehensive variability of the original hydrological sequence according to the trend variability, the mean variability, and the variance variability.

[0063] A device includes a memory, a processor, and a computer program stored in the memory and executable on the processor. When the processor executes the computer program, the steps of the method for diagnosing and measuring the overall variability of a hydrological sequence are implemented.

[0064] Compared with the prior art, the present invention has at least the following beneficial effects:

[0065] The method for diagnosing and measuring the overall variability of a hydrological sequence provided by the present invention starts with the variability diagnosis of the hydrological sequence, constructs a hydrological variability diagnosis system combining multiple methods, determines the possible types of variability and variability nodes through the tests of three types of variability, namely trend, mean, and variance. Moreover, the tests of each type of variability are carried out separately, and there is no mutual influence between the tests of different types of variability, avoiding the detection errors caused by the coupling and intersection of variability information, making the variability diagnosis result more credible, and at the same time being able to obtain relatively complete variability information of the hydrological sequence. Then, the variability of different types of variability is further measured, and the influence of different types of variability on the overall variability is comprehensively considered to quantitatively represent the overall variability degree of the hydrological sequence, providing a reliable method for measuring the variability of the hydrological sequence.

[0066] To make the above objects, features, and advantages of the present invention more obvious and understandable, the following specific preferred embodiments are given, and in conjunction with the accompanying drawings, the detailed description is as follows. Description of the Drawings

[0067] To more clearly illustrate the technical solutions in the specific embodiments of the present invention, the following will briefly introduce the drawings required for the description of the specific embodiments. Obviously, the following drawings are some embodiments of the present invention. For those of ordinary skill in the art, without creative efforts, other drawings can be obtained based on these drawings.

[0068] Figure 1 It is a flowchart of the method for diagnosing and measuring the overall variability of a hydrological sequence according to the present invention.

[0069] Figure 2 It is a schematic diagram of the overall empirical mode decomposition method in the method for diagnosing and measuring the overall variability of a hydrological sequence according to the present invention.

[0070] Figure 3 It is a schematic diagram of the annual runoff sequence of Dingjiagou Hydrological Station in the Wuding River Basin from 1960 to 2010 in the embodiment of the present invention.

[0071] Figure 4 It is a schematic diagram of the decomposition result of the original hydrological sequence applying the ensemble empirical mode decomposition method in the embodiment of the present invention.

[0072] Figure 5 It is a schematic diagram of the trend variation diagnosis result in the embodiment of the present invention.

[0073] Figure 6 It is a schematic diagram of the mean variation diagnosis result in the embodiment of the present invention.

[0074] Figure 7 It is a schematic diagram of the variance variation diagnosis result in the embodiment of the present invention.

[0075] Figure 8 It is a schematic diagram of the calculation result of the overall variation degree of the hydrological sequence in the embodiment of the present invention. Specific embodiments

[0076] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings. Apparently, the described embodiments are some but not all of the embodiments of the present invention. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts shall fall within the scope of protection of the present invention.

[0077] As a specific embodiment of the present invention, as Figure 1 shown, the present invention provides a method for diagnosing and measuring the overall variation degree of a hydrological sequence, which specifically includes the following steps:

[0078] S1. Obtain the original hydrological sequence, and decompose the original hydrological sequence into a plurality of IMF sequences and a trend item sequence;

[0079] S2. Judge whether there is a trend variation point in the trend item sequence. If so, calculate the trend variation degree of the trend item sequence;

[0080] S3. Stack and recombine the plurality of IMF sequences to obtain an IMFs sequence, and judge whether there is a mean variation point in the IMFs sequence. If so, calculate the means of the subsequences on both sides of the mean variation point and the mean variation degree of the IMFs sequence;

[0081] S4. Perform mean zeroing processing on the IMFs sequence according to the means of the subsequences on both sides of the mean variation point, and judge whether there is a variance variation point in the IMFs sequence after the mean zeroing processing. If so, calculate the variance variation degree of the IMFs sequence after the mean zeroing processing;

[0082] To avoid the influence of the mean mutation information on the variance mutation diagnosis, it is necessary to perform zero-mean processing on the IMFs sequence of the fluctuation term.

[0083] S5. Calculate the comprehensive mutation degree of the original hydrological sequence according to the trend mutation degree, the mean mutation degree and the variance mutation degree.

[0084] On the basis of the above embodiments, as a more preferred embodiment, regarding decomposing the original hydrological sequence into multiple IMF sequences and a trend term sequence in step S1, the ensemble empirical mode decomposition method (EEMD) is specifically adopted, which specifically includes the following steps:

[0085] 1) Add a white noise signal ω(t) with the same time scale and conforming to the normal distribution to the original hydrological sequence x(t) to obtain a new hydrological sequence y(t) after noise reduction processing:

[0086] y(t) = x(t) + ω(t)

[0087] Among them, the intensity of the white noise signal ω(t) can be determined by setting its standard deviation;

[0088] 2) Perform EMD decomposition on the new hydrological sequence y(t) to obtain multiple component sequences;

[0089] In this embodiment, more preferably, regarding performing EMD decomposition on the new hydrological sequence y(t) in step 2) to obtain multiple component sequences, it specifically includes the following steps:

[0090] 2.1) Extract the extreme points in the new hydrological sequence y(t). For the obtained maximum points and minimum points, use the cubic spline interpolation method to obtain their interpolation functions respectively. Perform interpolation processing on the maximum points to obtain the upper envelope Y of the sequence y(t) max(t) , and perform interpolation processing on the minimum points to obtain the lower envelope Y of the sequence y(t) min(t) ;

[0091] 2.2) Calculate the mean p(t) of the upper envelope Y of the new hydrological sequence y(t) max(t) and the lower envelope Y of the new hydrological sequence y(t) min(t) :

[0092] p(t) = (y max(t) + y min(t) ) / 2

[0093] 2.3) Calculate the difference between the new hydrological sequence y(t) and the mean p(t) to obtain the difference time series d(t):

[0094] d(t) = y(t) - p(t)

[0095] 2.4) Let \(a\) be the number of extreme points of the difference time series \(d(t)\), \(b\) be the number of zero-crossing points of the difference time series \(d(t)\), and \(u(t)\) be the mean value of the outer envelope of the difference time series \(d(t)\). If:

[0096] \(\vert a - b\vert\leq1\)

[0097] \(u(t)=0\)

[0098] Then the difference time series \(d(t)\) is an intrinsic mode component sequence of the new hydrological sequence \(y(t)\), denoted as \(I n (t)\); otherwise, re- interpolate and fit the upper and lower envelopes of the new hydrological sequence \(y(t)\).

[0099] 2.5) Subtract \(I n (t)\) from the new hydrological sequence \(y(t)\) to obtain the remainder sequence \(R n (t)\):

[0100] \(R n (t)=y(t)-I n (t)\)

[0101] 2.6) Determine whether the remainder sequence \(R n (t)\) is monotonic. If it is not monotonic, then take \(R n (t)\) as the new original sequence, and successively obtain \(I_2(t)\), \(I_3(t)\), \(\cdots\), \(I n (t)\) according to the above steps until \(R n (t)\) is a monotonic sequence. In this way, \(y(t)\) is decomposed into multiple component sequences:

[0102]

[0103] 3) Set the number of EMD decompositions as \(j\), repeat steps 1) - 2) \(j\) times, and average each of the component sequences obtained in the loop to obtain multiple IMF sequences and a trend term sequence.

[0104] The empirical mode decomposition method (EMD) has significant advantages in processing complex non - linear signals. It can split a complex time series into a trend term and a fluctuation term, and for an unknown time series, it can directly perform decomposition without prior analysis or processing. Simply put, EMD is like a sophisticated data screening machine that can automatically divide the mixed data into several parts, each containing different data information, and in this process, no human intervention and processing are required. And EEMD is an improvement of the EMD method. On the premise of retaining its advantages, by adding a white noise signal to the original time series, it solves the mode mixing problem in EMD decomposition, making the separation of different types of signals more complete.

[0105] On the basis of the above embodiments, as a more preferred embodiment, regarding the determination of whether there is a trend change point in the trend item sequence in step S2, it specifically includes the following steps:

[0106] 1) Use the self-segmented linear fitting method to perform segmented trend line fitting on the trend item sequence;

[0107] In this embodiment, more preferably, regarding the use of the self-segmented linear fitting method to perform segmented trend line fitting on the trend item sequence in step 1), it specifically includes the following steps:

[0108] 1.1) Set the trend item sequence as F(t i ), i = 1, 2,... n, let the number of segments of F(t) be k and the minimum truncation length be s;

[0109] 1.2) Perform segmented regression line fitting on the trend item sequence F(t i ) according to the least squares method, and obtain a set of fitting regression line equations as:

[0110] y i (j) = a j t i + b j , (j - 1)s < i ≤ n

[0111] In the formula, j = 1, 2,... k, representing the serial number of the segmented regression line; a j is the regression coefficient of the j-th segment of the regression line; b j is the vertical axis intercept of the j-th segment of the regression line;

[0112] 1.3) Calculate the sum of squared deviations S between the regression fitting line of the trend item sequence F(t i ) and the sequence values. The set of fitting regression line equations with the minimum S value is the best fitting result. The calculation formula for the sum of squared deviations S is as follows:

[0113]

[0114] 1.4) Let the correlation coefficients of the subsequences on both sides of a trend change point be r1 and r2 respectively. The difference between the two can reflect the degree of trend difference of the sequences on both sides of the inflection point. Construct a statistic χ 2 for hypothesis testing:

[0115]

[0116] In the formula, n1 and n2 are the lengths of the subsequences on both sides of the trend change point respectively, and n1 > 3, n2 > 3; z i is calculated through the correlation coefficient r i , and the formula is as follows:

[0117]

[0118] 2) According to χ 2 The test method is used to determine whether the piecewise trend line satisfies significance. If it satisfies significance, there is a trend change point. More specifically, select the significance level α, and compare the χ 2 test statistic with the upper quantile value If it is considered that there is a significant difference between r1 and r2, that is, the piecewise trend is significant, and the intersection point of the piecewise trend line is the trend change point.

[0119] On the basis of the above implementation manner, as a more preferred implementation manner, regarding determining whether there is a mean change point in the IMFs sequence in step S3, it specifically includes the following steps:

[0120] 1) Set the IMFs sequence as M(t), the number of its sequence points is n, and i is the possible change point of M(t), where 1 ≤ i ≤ n - 1, and 1 is the minimum truncation length of M(t), and its function is to ensure the effectiveness of the statistical test;

[0121] 2) Calculate the combined standard deviation σ D (i) at the i-th point, and the calculation formula is as follows:

[0122]

[0123] In the formula: μ r (i) and μ l (i) are the means of the sequence segments on both sides of M(i) respectively, and σ r (i) and σ l (i) are the standard deviations of the sequence segments on both sides of M(i) respectively; n r and n l are the numbers of sequence points on the left and right sides of the i-th point respectively;

[0124] 3) Construct a test statistic T(i) to quantitatively represent the mean difference between the sequence segments on the left and right sides of the i-th point, and the formula is as follows:

[0125]

[0126] In the formula: the numerical value of T(i) represents the degree of mean variation. The larger its value, the greater the amplitude of the mean jump at the i-th point. Calculate the corresponding statistic T(i) in turn within the value range of i, that is, obtain the test statistic sequence T(t) corresponding to M(t);

[0127] 4) Set the maximum value of the test statistic sequence T(t) as T m , and its significance level P is expressed as:

[0128]

[0129] η = 4.19ln(n) - 11.54

[0130] where v = n - 2; δ = 0.4; S x (a, b) is the incomplete beta function;

[0131] 5) If P ≥ P0, where P0 is the significance critical value and its value range is 0.5 to 0.95, it indicates that the mean difference is significant. Taking point i as the segmentation point, M(t) is divided into two sub - sequences, and point i is the mean mutation point;

[0132] 6) Perform the above operations on the two sub - sequences obtained in 5) respectively to check whether there are still mean mutation points in the sub - sequences. Eventually, c IMFs sub - sequences and (c - 1) mean mutation points will be obtained.

[0133] On the basis of the above - mentioned implementation manner, as a more preferred implementation manner, regarding the mean zero - normalization processing of the IMFs sequence according to the means of the sub - sequences on both sides of the mean mutation point in step S4, it specifically includes the following steps:

[0134] 1) Let an IMFs sub - sequence before processing be m(t), t = 1, 2, 3, …, is the mean of m(t);

[0135] 2) The processed IMFs sub - sequence is m*(t), and its calculation formula is as follows:

[0136]

[0137] 3) Combine the processed multiple IMFs sub - sequences to obtain the mean zero - normalized IMFs sequence.

[0138] On the basis of the above - mentioned implementation manner, as a more preferred implementation manner, regarding the judgment of whether there is a variance mutation point in the IMFs sequence after the mean zero - normalization processing in step S4, it specifically includes the following steps:

[0139] 1) Define the entire IMFs sequence after mean zero - normalization processing as Y t , t = 1, 2, 3, …, n, and its mean is 0, and the overall variance is Assume that there is a variance mutation point in the sequence at r, and construct the test statistic S max :

[0140]

[0141]

[0142]

[0143] In the formula, L(θ) is the maximum likelihood function; λ is the degree of freedom during calculation. When calculating for the entire sequence, λ = 1, and when calculating for the variance mutation point, λ = 2; n is the sequence length; is the variance of the subsequence before the r point, is the variance of the subsequence after the r point;

[0144] 2) For C(α) > 1, where α is the significance level, if S max ≥ C(α), then the variance mutation at the r point is significant, and the r point is a variance mutation point;

[0145] 3) Apply the bisection method to divide Y t into two subsequences at the variance mutation point, and then perform variance mutation tests on each subsequence according to the above steps to obtain all the variance mutation points of Y t .

[0146] On the basis of the above embodiments, as a more preferred embodiment, for the calculation methods of the trend mutation degree of the trend term sequence, the mean mutation degree of the IMFs sequence, and the variance mutation degree of the IMFs sequence after mean normalization, the correlation integral method is respectively used for calculation, specifically as follows:

[0147] 1) Define the sequence for which the mutation degree is to be calculated as Z t , t = 1, 2, 3,..., N, then the vertical distance between any two points in Z t constitutes a new sequence D s :

[0148] D s = ||Z i - Z j ||

[0149] In the formula, s = 1, 2,..., N(N - 1) / 2; i = 1, 2,..., N - 1; j = i + 1, i + 2,..., N;

[0150] 2) Define the statistic I(N, R) to measure the non - stationarity degree of Z t :

[0151]

[0152] In the formula, N is the total number of points; R is the reference distance, generally determined by the statistical distribution law of D s ; H is the unit step function, defined as follows:

[0153]

[0154] 3) If there are mutation points in the test sequence, calculate the non-stationarity degrees I1(N1, R1) and I2(N2, R2) of the sequences before and after the mutation points respectively, and obtain the mutation degree U corresponding to the trend term sequence:

[0155] U = |I1(N1, R1) - I2(N2, R2)|.

[0156] The correlation integral method is a method for calculating mutation scales through probability statistics. When the mutation points are known, calculate the non-stationarity degrees of the subsequences on both sides of the mutation points respectively, and then measure the variability through the change in non-stationarity before and after the mutation, so as to classify and grade the mutations through quantitative values.

[0157] On the basis of the above embodiments, as a more preferred embodiment, regarding calculating the comprehensive mutation degree of the original hydrological sequence according to the trend mutation degree, the mean mutation degree, and the variance mutation degree in step S5, specifically:

[0158] Add the trend mutation degree, the mean mutation degree, and the variance mutation degree to obtain the comprehensive mutation degree of the original hydrological sequence.

[0159] The present invention decomposes the original hydrological sequence to obtain subsequences containing different mutation information, then conducts mutation tests on different subsequences, calculates the scales of different types of mutations respectively, and finally obtains the hydrological mutation degree considering multiple mutations comprehensively. Specifically, the low-frequency component (trend term sequence) of the original hydrological sequence is separated, and the trend mutation detection is carried out on the trend term sequence alone. When the trend term component is separated, the high-frequency component (IMFs sequence) of the original sequence is subjected to subsequent jump mutation identification to obtain the mean mutation information of the fluctuation term (IMFs) sequence. The mean of the IMFs sequence is normalized so that the overall mean of the processed time series is 0, and the variance mutation test is carried out, thus eliminating the influence of the trend term and the mean jump component on the variance mutation test. That is to say, through the trend mutation test, the trend analysis of the trend term sequence is carried out to obtain its detailed trend change law. For the high-frequency component of the sequence, the jump mutation points are identified through mean mutation and variance mutation, avoiding the detection error caused by the simultaneous action of several different mutations, so as to obtain a more accurate mutation diagnosis result. In this way, the mutation degrees of different mutation types can be calculated respectively, and the influence degree of different mutations on the original hydrological sequence can be considered comprehensively to obtain the overall mutation degree of the hydrological sequence.

[0160] Embodiment

[0161] Step 1. Select the Wuding River Basin as the experimental area. The rainfall during the flood season in this basin features heavy rain with short rainfall duration and high rainfall intensity. The Wuding River Basin has a dry climate, low vegetation coverage rate, intense soil erosion activities, and an increasingly serious soil and water loss situation. It is one of the main sources of coarse sediment deposition in the middle reaches of the Yellow River. The data used in this invention is sourced from the Dingjiagou Hydrological Station on the main stream of the Wuding River. The time span of the data is from 1960 to 2010, and this data segment contains a 51-year annual runoff sequence. In this embodiment, a method for diagnosing the overall variability of hydrological sequences is applied to detail the steps for diagnosing the variability type of the annual runoff sequence at this station and the calculation process of hydrological variability. The annual runoff sequence of the Dingjiagou Hydrological Station is as Figure 3 shown.

[0162] Step 2. Apply the Ensemble Empirical Mode Decomposition (EEMD) method to decompose the original annual runoff sequence into a trend term sequence (Residual) and four Intrinsic Mode Function subsequences (IMF1, IMF2, IMF3, and IMF4). The decomposed subsequences are as Figure 4 shown. It can be seen from the figure that the trend term sequence (Residual) generally shows a concave decreasing trend, and the attenuation rate gradually slows down.

[0163] Step 3. Use the piecewise linear fitting method to fit the trend line of the trend term sequence. Set the minimum truncation length to 15. Finally, it is found that when the number of segments is 2, there is a trend mutation point that satisfies the significance test. The parameter values of the piecewise trend line obtained are shown in Table 1 below. It can be seen from the table that 1985 is the trend mutation point, and the linear regression coefficients before and after the mutation point are -0.206 and -0.0677 respectively. Then, a significance test for the trend mutation point is carried out at a significance level of α = 0.01. By calculation, the χ 2 value is 16.58, which is much larger than the upper quantile value satisfying the significance test. The trend mutation diagnosis result is as Figure 5 shown. The trend term sequence generally shows a decreasing trend, but the decreasing trend of the sequence before 1985 is significantly greater than that after 1985.

[0164] Table 1. Parameter values of the piecewise trend line using the piecewise linear fitting method

[0165]

[0166] Step 4. Combine the four IMF subsequences to obtain the IMFs sequence, and then identify the mean mutation information. First, set the significance critical value P0 = 0.9, and calculate according to the mean mutation test method described above. Finally, two mean mutation points that satisfy the significance test are obtained, namely i = 1984 and i = 1971, and the calculation results of the mutation information at each mutation point are shown in Table 2 below. AsFigure 6 As shown, there are two mean jump variations in the IMFs sequence. These two mean variations divide the IMFs sequence into three subsequences, and their means are 0.355, -1.154, and 0.552 from left to right in turn.

[0167] Table 2. Calculation Results of Mean Variation Test

[0168]

[0169] Note: i represents the position of the mean variation point, that is, the year; μ l and μ r respectively represent the means of the sequence segments on the left and right sides of the variation point i; σ D is the combined standard deviation at the variation point i; T m is the quantified mean difference between the sequence segments on both sides of the variation point i; P is the significance level of the variation point i; P0 is the significance critical value.

[0170] Step 5. In order to conduct the next variance variation diagnosis, first, according to the results of the mean variation test, perform mean zeroing on the IMFs sequence to obtain the transformed time series, and then conduct a variance variation test on the processed IMFs sequence. It is calculated that there is a variance jump variation in the mean-zeroed IMFs sequence around 1981, and the standard deviations of the sequence segments on both sides of the variation point are 1.809 and 0.984 respectively. At the significant level of α = 0.05, the statistical test statistic S max = 9.28, which is greater than its corresponding test statistic critical value C(α) = 9.171, indicating that the variance variation is significant. Taking 1981 as the boundary, draw the standard deviations of the sequence segments on both sides of the variation point on the graph, and obtain the variance variation diagnosis graph as Figure 7 shown. It can be seen from the graph that the standard deviation of the mean-zeroed IMFs sequence has a significant difference before and after 1981, that is, a variance jump variation occurred around 1981.

[0171] Step 6. According to the variation identification results obtained in the above steps, use the relevant integral method to calculate the trend, mean, and variance variation degrees of the hydrological sequence respectively. The sum of the three is the overall variation degree of the hydrological sequence. When calculating, the reference distance is taken as R = μ D + σ D , that is, the sum of the mean and standard deviation of the sequence D s . Since D s follows the positive half-axis part of the standard normal distribution, according to the probability statistics theory, D sThe probability of <R> is 0.683. Therefore, 0.317 and 0.1 are taken as the variability demarcation points, and the variability classification table is shown in Table 3 below. The calculation results of the annual runoff series variability of Dingjiagou Station are shown in Table 4 below. As can be seen from Table 4, for the annual runoff series of Dingjiagou Station, its overall variability is 0.1667, belonging to the medium variability intensity, and the mean jump variability is the main variability component, having the greatest impact on the overall variability. The schematic diagram of the calculation results of the overall variability of the annual runoff series of Dingjiagou Station is as Figure 8 shown.

[0172] Table 3. Variability Classification Table

[0173]

[0174]

[0175] Table 4. Calculation Results of Hydrological Variability

[0176]

[0177] To sum up, the annual runoff series of Dingjiagou Hydrological Station from 1960 to 2010 has a triple composite variability type. First, a trend variability occurred in 1985, making the downward trend of the annual runoff series significantly slower. Second, after removing the influence of the trend component, there are three jump mutation points in the high-frequency component of the annual runoff series of Dingjiagou Station, namely two mean variability points and one variance variability point. The two mean variabilities occurred in 1971 and 1984 respectively. Among them, an attenuated jump mutation occurred in 1971, while a climbing jump mutation occurred in 1984. The variance variability occurred in 1981, and through the sequence standard deviations on both sides of the variance variability point, it can be seen that the sequence fluctuation from 1960 to 1981 is significantly more complex than that from 1981 to 2010. The interweaving of the three variability information makes the overall annual runoff series of this station show a downward trend. Based on the identification of the variability points, further variability calculations show that the mean jump variability is the main variability component of the annual runoff series of Dingjiagou Station. Calculating different types of variabilities separately makes the obtained overall variability more credible and can reflect the variability degree of the hydrological series.

[0178] The comprehensive variation diagnosis of hydrological sequences is of great significance for the attribution analysis of changes in hydrological elements. Errors in the variation diagnosis results will, to a certain extent, lead to the distortion of attribution analysis. The present invention first decomposes the annual runoff sequence by using the method of ensemble empirical mode decomposition, then identifies the variations of subsequences containing different variation information, constructs a hydrological variation diagnosis system combining multiple methods, and calculates the variation degrees of the test sequences by classification to quantitatively represent the variability of hydrological sequences. Through case tests, this method can accurately and effectively obtain the variation information of hydrological sequences by testing and identifying different variation types. Moreover, there is no mutual influence between the diagnoses of different variation types, avoiding the detection errors caused by the coupling and intersection of variation information, and the variation diagnosis results are accurate and reasonable. In addition, the scale calculation of hydrological variation can quantitatively represent the influence degree of different types of variation on the overall variation, reflect the overall variation degree of hydrological sequences, and provide a reliable method for measuring the variability of hydrological sequences.

[0179] The present invention provides a device for diagnosing and measuring the overall variation degree of hydrological sequences, including:

[0180] A decomposition module, configured to obtain an original hydrological sequence and decompose the original hydrological sequence into a plurality of IMF sequences and a trend item sequence.

[0181] A trend variation degree calculation module, configured to determine whether there is a trend variation point in the trend item sequence. If so, calculate the trend variation degree of the trend item sequence.

[0182] A mean variation degree calculation module, configured to stack and recombine a plurality of the IMF sequences to obtain an IMFs sequence, determine whether there is a mean variation point in the IMFs sequence. If so, calculate the means of the subsequences on both sides of the mean variation point and the mean variation degree of the IMFs sequence.

[0183] A variance variation degree calculation module, configured to perform mean normalization processing on the IMFs sequence according to the means of the subsequences on both sides of the mean variation point, determine whether there is a variance variation point in the mean-normalized IMFs sequence. If so, calculate the variance variation degree of the mean-normalized IMFs sequence.

[0184] A comprehensive variation degree calculation module, configured to calculate the comprehensive variation degree of the original hydrological sequence according to the trend variation degree, the mean variation degree, and the variance variation degree.

[0185] The device for diagnosing and measuring the overall variation degree of hydrological sequences provided by the present invention is used to implement the foregoing method for diagnosing and measuring the overall variation degree of hydrological sequences. Therefore, the specific implementation manners in the device for diagnosing and measuring the overall variation degree of hydrological sequences can be seen in the embodiment part of the method for diagnosing and measuring the overall variation degree of hydrological sequences in the foregoing text.

[0186] In one embodiment of the present invention, a computer device is provided, which includes a processor and a memory. The memory is used to store a computer program, and the computer program includes program instructions. The processor is used to execute the program instructions stored in the computer storage medium. The processor may be a central processing unit (CPU), or may also be other general-purpose processors, digital signal processors (DSPs), application specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. It is the computing core and control core of the terminal, and is suitable for implementing one or more instructions. Specifically, it is suitable for loading and executing one or more instructions to implement the corresponding method flow or corresponding function. The processor described in the embodiments of the present invention can be used to implement the operations of a method for diagnosing and measuring the overall variability of a hydrological sequence.

[0187] In one embodiment of the present invention, when a method for diagnosing and measuring the overall variability of a hydrological sequence is implemented in the form of a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on such an understanding, to implement all or part of the processes in the above-described embodiment methods of the present invention, it can also be completed by instructing relevant hardware through a computer program. The computer program can be stored in a computer-readable storage medium. When the computer program is executed by a processor, the steps of the above-described various method embodiments can be implemented. Among them, the computer program includes computer program code, and the computer program code can be in the form of source code, object code, executable file, or some intermediate form, etc. The computer-readable storage medium includes permanent and non-permanent, removable and non-removable media, and information storage can be achieved by any method or technology. The information can be computer-readable instructions, data structures, program modules, or other data.

[0188] The computer storage medium can be any available medium or data storage device that can be accessed by a computer, including but not limited to magnetic memories (such as floppy disks, hard disks, magnetic tapes, magneto-optical disks (MO), etc.), optical memories (such as CDs, DVDs, BDs, HVDs, etc.), and semiconductor memories (such as ROMs, EPROMs, EEPROMs, non-volatile memories (NAND FLASH), solid state drives (SSD)).

[0189] Those skilled in the art should understand that the embodiments of the present application can be provided as a method, a system, or a computer program product. Therefore, the present application can take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware aspects. Moreover, the present application can take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) that contain computer-usable program code.

[0190] The present application is described with reference to the flowcharts and / or block diagrams of the methods, apparatuses (systems), and computer program products of the embodiments of the present application. It should be understood that each flow and / or block in the flowchart and / or block diagram can be implemented by computer program instructions, and the combination of the flows and / or blocks in the flowchart and / or block diagram can also 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, such that the instructions executed by the processor of the computer or other programmable data processing devices generate a device for implementing the specified functions in Figure 1 one flow or multiple flows and / or blocks Figure 1 one block or multiple blocks.

[0191] 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, such that the instructions stored in the computer-readable memory generate a manufactured article including an instruction device that implements the specified functions in Figure 1 one flow or multiple flows and / or blocks Figure 1 one block or multiple blocks.

[0192] These computer program instructions can also be loaded onto a computer or other programmable data processing device, such that a series of operation steps are executed on the computer or other programmable device to generate a computer-implemented process, and thus the instructions executed on the computer or other programmable device provide steps for implementing the specified functions in Figure 1 one flow or multiple flows and / or blocks Figure 1 one block or multiple blocks.

[0193] Finally, it should be noted that the above embodiments are only examples for clearly explaining the present invention. The application scope of the method of the present invention is obviously not limited to a specific basin or region, and the selected calculation data is not a limitation on the type of hydrological sequence. For those of ordinary skill in the art, there may be some different forms of changes or modifications based on the above description. It is not necessary and impossible to list all implementation manners here. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present invention shall be included within the protection scope of the claims of the present invention.

Claims

1. A diagnostic measurement method for the overall variability of a hydrological sequence, characterized in that, Including: Obtain the original hydrological sequence, and decompose the original hydrological sequence into a plurality of IMF sequences and a trend item sequence; Judge whether there is a trend mutation point in the trend item sequence. If so, calculate the trend mutation degree of the trend item sequence; The judgment of whether there is a trend mutation point in the trend item sequence is specifically: Use the self-segmented linear fitting method to perform segmented trend line fitting on the trend item sequence; According to χ 2 Use the test method to determine whether the segmented trend line meets the significance. If it meets the significance, there is a trend change point; Stack and recombine the multiple IMF sequences to obtain an IMFs sequence, and judge whether there is a mean mutation point in the IMFs sequence. If so, calculate the mean values of the subsequences on both sides of the mean mutation point and the mean mutation degree of the IMFs sequence; Perform mean zeroing on the IMFs sequence according to the mean values of the subsequences on both sides of the mean mutation point, and judge whether there is a variance mutation point in the IMFs sequence after the mean zeroing. If so, calculate the variance mutation degree of the IMFs sequence after the mean zeroing; The judgment of whether there is a variance mutation point in the IMFs sequence after the mean zeroing is specifically: 1) Define the entire IMF sequence after mean normalization as Y t , where \(t = 1, 2, 3, \ldots, n\), and its mean is 0, and the overall variance is Suppose there is a variance change point at \(r\) in the sequence, and construct the test statistic \(S\) max : where n is the sequence length; is the variance of the subsequence before the r-th point, is the variance of the subsequence after the r-th point; 2) For C(α) > 1, where α is the significance level, if S max ≥ C(α), then the variance variation at point r is significant, and point r is a variance variation point; 3) Apply the bisection method to divide Y at the variance change point t into two sub - sequences, and then perform variance change tests on each sub - sequence according to the above steps to obtain all the variance change points of Y t ; Calculate the comprehensive mutation degree of the original hydrological sequence according to the trend mutation degree, the mean mutation degree, and the variance mutation degree.

2. The method for diagnosing and measuring the overall variability of a hydrological sequence according to claim 1, wherein The decomposition of the original hydrological sequence into a plurality of IMF sequences and a trend item sequence includes: 1) Add a white noise signal ω(t) with the same time scale and conforming to the normal distribution to the original hydrological sequence x(t) to obtain a new hydrological sequence y(t) after noise reduction processing: y(t) = x(t) + ω(t) 2) Perform EMD decomposition on the new hydrological sequence y(t) to obtain a plurality of component sequences; 3) Set the number of EMD decomposition times as j, repeat steps 1) to 2) j times, and average each component sequence obtained by the loop to obtain a plurality of IMF sequences and a trend item sequence.

3. A method for diagnosing and measuring the overall variability of a hydrological sequence according to claim 1, characterized in that The judgment of whether there is a mean mutation point in the IMFs sequence includes: 1) Set the IMFs sequence as M(t), the number of its sequence points as n, and i as the possible mutation point of M(t), where 1 ≤ i ≤ n - 1, and 1 is the minimum truncation length of M(t); 2) Calculate the combined standard deviation σ D (i) at point i, and the calculation formula is as follows: Where: μ r (i) and μ l (i) are the means of the sequence segments on both sides of M(i), and σ r (i) and σ l (i) are the standard deviations of the sequence segments on both sides of M(i); n r and n l are the numbers of sequence points on the left and right sides of point i, respectively; 3) Construct a test statistic T(i) to quantitatively represent the mean difference between the sequence segments on the left and right sides of point i. The formula is as follows: In the formula: the numerical value of T(i) represents the degree of mean mutation. The larger its value, the greater the amplitude of the mean jump at point i. Calculate the corresponding statistic T(i) in the range of i values in turn, and obtain the test statistic sequence T(t) corresponding to M(t); 4) Set the maximum value of the sequence of the test statistic T(t) to T m , and its significance level P is expressed as: η = 4.19ln(n) - 11.54 where \(v = n - 2\); \(\delta=0.4\); \(S\) x (a, b) is the incomplete beta function; 5) If P ≥ P0, where P0 is the significance critical value, and its value range is 0.5 to 0.95, it means that the mean difference is significant. Taking point i as the segmentation point, divide M(t) into two subsequences, and point i is the mean mutation point; 6) Perform the above operations on the two subsequences obtained in step 5) respectively to check whether there is still a mean mutation point in the subsequences, and finally obtain c IMFs subsequences and (c - 1) mean mutation points.

4. A method for diagnosing and measuring the overall variability of a hydrological sequence according to claim 1, characterized in that, The mean zeroing of the IMFs sequence according to the mean values of the subsequences on both sides of the mean mutation point includes: 1) Let a subsequence of IMFs before processing be \(m(t)\), where \(t = 1, 2, 3, \ldots, n\). is the mean value of \(m(t)\); 2) The processed IMF subsequences are \(m^*(t)\), and its calculation formula is as follows: 3) Merge the processed multi-segment IMF subsequences to obtain the mean-zeroed IMF sequence.

5. A method for diagnosing and measuring the overall variability of a hydrological sequence according to claim 1, characterized in that The trend variability of the trend term sequence, the mean variability of the IMF sequence, and the variance variability of the mean-zeroed IMF sequence are calculated using the correlation integral method, specifically as follows: 1) Define the sequence for which the variability is to be calculated as Z t , where t = 1, 2, 3, …, N, then for Z t , the vertical distances between any two points form a new sequence D s : D s = ||Z i -Z j || where \(s = 1, 2, \ldots, N(N - 1) / 2\); \(i = 1, 2, \ldots, N - 1\); \(j = i + 1, i + 2, \ldots, N\); 2) Define a statistic I(N, R) to measure the non-stationarity of Z t : Where N is the total number of points; R is the reference distance, determined by the statistical distribution law of D s ; H is the unit step function, defined as follows: 3) If there is a breakpoint in the test sequence, calculate the non-stationarity degrees \(I_1(N_1, R_1)\) and \(I_2(N_2, R_2)\) of the sequences before and after the breakpoint respectively to obtain the variability \(U\) corresponding to the trend term sequence: \(U = |I_1(N_1, R_1) - I_2(N_2, R_2)|\).

6. A method for diagnosing and measuring the overall variability of a hydrological sequence according to claim 1, characterized in that, The calculation of the comprehensive variability of the original hydrological sequence according to the trend variability, the mean variability, and the variance variability includes: Add the trend variability, the mean variability, and the variance variability to obtain the comprehensive variability of the original hydrological sequence.

7. A diagnostic measurement device for the overall variability of a hydrological sequence, characterized in that, Including: A decomposition module for obtaining the original hydrological sequence and decomposing the original hydrological sequence into multiple IMF sequences and a trend term sequence; A trend variability calculation module for determining whether there is a trend breakpoint in the trend term sequence. If so, calculate the trend variability of the trend term sequence; The determination of whether there is a trend breakpoint in the trend term sequence is specifically: Use the self-segmented linear fitting method to perform segmented trend line fitting on the trend term sequence; According to χ 2 Use the test method to determine whether the piecewise trend line satisfies significance. If it satisfies significance, there is a trend change point; A mean variability calculation module for stacking and recombining multiple IMF sequences to obtain an IMF sequence, determining whether there is a mean breakpoint in the IMF sequence. If so, calculate the means of the subsequences on both sides of the mean breakpoint and the mean variability of the IMF sequence; A variance variability calculation module for performing mean-zeroing processing on the IMF sequence according to the means of the subsequences on both sides of the mean breakpoint, determining whether there is a variance breakpoint in the mean-zeroed IMF sequence. If so, calculate the variance variability of the mean-zeroed IMF sequence; The determination of whether there is a variance breakpoint in the mean-zeroed IMF sequence is specifically: 1) Define the entire sequence of IMFs after mean normalization as Y t , where \(t = 1, 2, 3, \ldots, n\), and its mean is 0, and the overall variance is Assume that there is a variance change point at \(r\) in the sequence, and construct the test statistic \(S\ max : Where n is the sequence length; is the variance of the subsequence before the r-th point, is the variance of the subsequence after the r-th point; 2) For C(α) > 1, where α is the significance level, if S max ≥ C(α), then the variance variation at point r is significant, and point r is a variance variation point; 3) Apply the bisection method to divide Y at the variance change point t into two sub-sequences, and then perform variance change tests on each sub-sequence according to the above steps to obtain all variance change points of Y t ; A comprehensive variability calculation module for calculating the comprehensive variability of the original hydrological sequence according to the trend variability, the mean variability, and the variance variability.

8. An apparatus, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the steps of a method for diagnosing and measuring the overall variability of a hydrological sequence as described in any one of claims 1 to 6.

Images