Photovoltaic station regional surface deformation monitoring and prediction integration method

Abstract

The invention discloses a photovoltaic station regional surface deformation monitoring and prediction integration method. The method comprises the following steps: acquiring a photovoltaic station surface deformation monitoring result; dividing the surface deformation area of the photovoltaic station according to the annual surface deformation rate; performing slope unit division on the photovoltaic station by using a positive and negative DEM method to obtain a slope unit-based photovoltaic station earth surface deformation grade; extracting the displacement time sequence in each subdivision deformation unit by adopting a complete random sampling technology, and calculating an average displacement time sequence in each subdivision unit; performing signal decomposition on the average displacement time sequences to obtain trend, seasonal and periodic components of each average displacement time sequence; key rainfall characteristic parameters for various displacement components of the photovoltaic station are obtained by adopting an MIC correlation analysis technology; and constructing a photovoltaic station deformation unit displacement combined prediction model based on a DHKELM method as a core, and realizing accurate prediction of the deformation condition of each region of the photovoltaic station.

Description

Technical Field

[0001] This invention belongs to the technical field of geological disaster prediction and forecasting, and particularly relates to an integrated method for monitoring and predicting surface deformation in photovoltaic power station areas. Background Technology

[0002] With the advancement of the "dual carbon" goals, the photovoltaic power generation industry has developed rapidly, making the construction and safe operation of large-scale photovoltaic power plants crucial. However, photovoltaic power plants are typically located in areas with complex geological conditions, such as mountains, hills, and coal mining subsidence areas. These sites cover vast areas and face risks of geological disasters such as landslides and ground subsidence, seriously threatening their structural safety and stable power generation. Therefore, large-scale, high-precision, and continuous monitoring and prediction of surface deformation for photovoltaic power plants has significant economic and safety implications.

[0003] Traditional methods for monitoring surface deformation, such as GPS and leveling, while highly accurate, suffer from limitations including high cost, low efficiency, and difficulty in achieving large-area, dense monitoring, failing to meet the needs of continuous, full-site monitoring of photovoltaic power plants. In recent years, Synthetic Aperture Radar Interferometry (InSAR) technology, particularly Small Baseline Set Assemblies (SBAS-InSAR), has become an important tool for surface deformation monitoring due to its wide-area, high-precision, and all-weather monitoring capabilities. However, current technologies face several bottlenecks: First, the massive deformation point data acquired by InSAR lacks organic integration with the engineering geological units of photovoltaic power plants, making it difficult to directly use for zonal risk assessment; second, surface deformation time series contain complex trend, seasonal, and periodic signals, which traditional methods struggle to adaptively and accurately separate; and finally, existing prediction models often fail to fully consider the complex nonlinear relationships between photovoltaic power plant deformation and external driving factors such as rainfall, and their generalization ability and prediction accuracy need improvement.

[0004] Therefore, there is an urgent need for a comprehensive approach that integrates advanced monitoring technology, intelligent data processing, and high-precision prediction models to achieve full-chain automation and intelligence in monitoring, analyzing, and predicting surface deformation across the entire photovoltaic power station area. Summary of the Invention

[0005] The technical problem to be solved by this invention is to provide an integrated method for monitoring and predicting surface deformation in photovoltaic power station areas, in order to address the following issues in existing technologies: the massive amount of deformation point data acquired lacks organic integration with the engineering geological units of photovoltaic power stations, making it difficult to directly use for zonal risk assessment; the surface deformation time series contains complex trend, seasonal, and periodic signals, which traditional methods cannot adaptively and accurately separate; existing prediction models often fail to fully consider the complex nonlinear relationship between photovoltaic power station deformation and external driving factors such as rainfall, and the generalization ability and prediction accuracy of the models need to be improved.

[0006] Technical solution of the present invention: An integrated method for monitoring and predicting surface deformation in photovoltaic power station areas, the method comprising: Step 1: Collect Sentinel-1 C-band SAR images, DEM data, and precise orbit data for the photovoltaic power station during the specified monitoring period to obtain the surface deformation monitoring results of the photovoltaic power station; Step 2: Based on the obtained surface deformation monitoring results, the Som-Kmeans two-step clustering method is used to divide the surface deformation area of ​​the photovoltaic power station according to the annual surface deformation rate, and the area is divided into three categories: small deformation, medium deformation and large deformation. Step 3: Based on the clustering results of the deformation area of ​​the photovoltaic power station, the forward and reverse DEM methods are used to divide the photovoltaic power station into slope units, and the clustering results of the deformation area of ​​the photovoltaic power station are further subdivided to obtain the result of the classification of the surface deformation level of the photovoltaic power station based on the slope units. Step 4: Based on the results of the subdivision of the deformation area of ​​the photovoltaic power station, and combined with the monitoring results of the surface deformation of the photovoltaic power station, the displacement time series of 50 locations in each subdivided deformation unit is extracted using a completely random sampling technique, and the average displacement time series in the subdivided unit is calculated. Step 5: Use the OVMD adaptive decomposition method to decompose the average displacement time series in each subdivided deformation unit to obtain the trend, seasonal and periodic components of each average displacement time series. Step 6: Use MIC correlation analysis technology to analyze the correlation between various rainfall characteristic parameters in the region and the trend, seasonal and periodic components of each average displacement time series, and obtain the key rainfall characteristic parameters for various displacement components of photovoltaic power stations. Step 7: Based on key rainfall characteristic parameters, construct a photovoltaic power station deformation unit displacement combination prediction model with the DHKELM method as the core to achieve accurate prediction of deformation in various areas of the photovoltaic power station.

[0007] The beneficial effects of this invention are: This invention acquires the time series of surface displacement at photovoltaic power plants using SBAS-InSAR technology, and combines it with Som-kmeans two-step clustering. It then uses OVMD adaptive decomposition to decompose the nonlinear displacement time series into trend displacement and seasonal displacement, and analyzes the correlation between displacement components and influencing factors using average MIC correlation analysis. Based on DHKELM kernel-based displacement prediction algorithms, it effectively achieves accurate prediction of surface displacement at photovoltaic power plants, providing crucial data support and technical assurance for early warning and prevention of geological disasters at photovoltaic power plants. This results in more effective and accurate prediction of surface deformation in the study area, providing strong data support for early warning and prevention of geological disasters at photovoltaic power plants.

[0008] This invention systematically solves the core problem in monitoring and predicting surface deformation of photovoltaic power plants through a highly integrated technical solution. To address the lack of organic integration between massive deformation point data and engineering geological units, this study first uses SBAS-InSAR technology to acquire surface deformation monitoring results from photovoltaic power plants. Then, a two-step Som-Kmeans clustering method is employed, using initial clustering with a Som neural network and secondary subdivision with Kmeans to divide the entire power plant into small, medium, and large deformation regions based on annual deformation rates. Simultaneously, using forward and reverse DEM methods, including DEM depression filling, slope and flow direction calculation, reverse DEM construction, and overlaying of watershed lines and slope toe lines, slope units are delineated. Finally, the deformation region clustering results are overlaid with the slope unit delineation results to achieve spatial correlation between deformation data and engineering geological units, providing a foundation for zoning risk assessment. To address the difficulty in accurately separating complex surface deformation time series signals, the average displacement time series of subdivided deformation units is used as input. An OVMD adaptive decomposition method is employed, using K-value optimization to determine the optimal modal decomposition number and α-value optimization to select the best penalty parameter, combined with FFT. Spectral analysis and power spectral density calculation are used to separate trend displacement (low frequency) and seasonal displacement (high frequency) according to frequency thresholds, and further identify the annual cycle and long cycle components in seasonal displacement to achieve adaptive and accurate decomposition of complex signals. To address the problem of insufficient consideration of the nonlinear relationship between deformation and external driving factors and low prediction accuracy, MIC correlation analysis technology is first used to screen key driving factors by taking monitoring time as the influencing factor for trend displacement and multi-dimensional rainfall characteristic parameters (cumulative rainfall, effective rainfall, etc.) as the influencing factors for seasonal and periodic displacement. Then, a combined prediction model with DHKELM as the core is constructed. Deep features are extracted by ELM-AE stacking, and the nonlinear fitting ability is enhanced by combining a hybrid kernel function of linear kernel and radial basis kernel. The trend, seasonal and periodic displacement components of each sub-unit are predicted separately. Finally, the three types of components are superimposed to obtain high-precision and high-generalization deformation prediction results for the entire station. Attached Figure Description

[0009] Figure 1 This is a schematic diagram of the process of the present invention. Detailed Implementation

[0010] This invention provides an integrated method for monitoring and predicting surface deformation in photovoltaic power station areas, comprising: Step 1: Collect Sentinel-1 C-band SAR images, DEM data, and precise orbit data for the photovoltaic power station during the specified monitoring period. Use SBAS-InSAR technology to obtain the surface deformation monitoring results of the photovoltaic power station. The implementation steps include:

[0011] Step 1.1, Data Acquisition: Sentinel-1A Single Look Complex (SLC) data, precise orbit data, and reference DEM data for the study area during the study period were obtained from publicly available platforms.

[0012] Step 1.2: Obtain surface deformation information: (1) Based on the data obtained above, the Small Baseline Subset (SBAS-InSAR) method was used to process the surface deformation information of the study area. The N+1 Sentinel-1A Single Look Complex (SLC) SAR images were sorted in chronological order: Select and register the super master image to generate M interferograms, where M satisfies the following conditions: ; (2) In The two SAR images acquired at time 1 generate the first j Interferogram, with flat and topographic phases removed respectively, azimuth coordinate x y The differential interference phase at the range coordinate r pixel can be expressed by the following formula. ; In the formula, The center wavelength of the radar. and They are respectively t b and t a Time relative to t 0 Cumulative deformation of radar line-of-sight direction at any given time. and They are respectively derived from deformation variables and The resulting change in phase deformation value.

[0013] (3) Any interferogram in equation (2.2) j The deformation phase can be represented by the average rate over the corresponding time period. v j This means, that is: ; Equation (2.2) can then be expressed as: ; in Denotes the end time corresponding to the j-th interferogram. This represents the starting time corresponding to the j-th interferogram. v j Represents any interferogram j The average rate of deformation phase within the corresponding time period, This represents the deformation phase of the j-th interferogram.

[0014] Equation (2.4) can be written in matrix form: ; Where A is an M×N coefficient matrix, v The deformation rate vector, Let A be the matrix composed of the phases of M differential interferograms. When M ≥ N, the rank of the coefficient matrix A is N, which can be solved using the least squares criterion according to equation (2.6).

[0015] ; in Let A be the transpose of the coefficient matrix A. For matrix The inverse matrix.

[0016] Find Then, the deformation phase, deformation amount, and corresponding deformation rate vector for each imaging time period can be calculated. ).

[0017] When M < N, the equation has infinitely many solutions. The minimum norm solution of the deformation rate vector is obtained by using the singular value decomposition method. The rate is integrated over each time period to obtain the deformation of each time period, thus obtaining the time series of surface displacement of the photovoltaic power station.

[0018] Solve equation (2.5) using the singular value decomposition method. ; In the formula, U is A T An M×M orthogonal matrix within A; R is an M×M matrix within A. T An M×N orthogonal matrix in A; A T The diagonal matrix of elements on the inner diagonal of A is C. Let A be the transpose of matrix R. Let the rank of A be H, then A... T The first H eigenvalues ​​of A are not 0, and the remaining eigenvalues ​​are all 0. Therefore, the above equation can be changed to: In the formula, It is the l-th element in C; , Let l be the l-th element in U and R.

[0019] Therefore, an estimated value can be obtained. Replace at average rate The formula is ; In the formula, G is a matrix where G(j′, B) is 0. v G is the deformation rate vector. G can be obtained by calculating it using singular value decomposition. The minimum norm solution is the deformation rate vector. Assume the duration of each time period is... Then the deformation of each time period It can be obtained through the deformation rate vector Duration of Time Period The product of is obtained, that is

[0020] ; Based on the deformation of the photovoltaic power station at different time periods, a time series of surface displacement of the photovoltaic power station can be constructed. Assume the initial displacement is... Then the displacement d at the k-th time point k It is the sum of the deformation variables of the first k time periods, i.e.

[0021] ; Step 2: Based on the obtained surface deformation monitoring results, the Som-Kmeans two-step clustering method is used to divide the surface deformation areas of the photovoltaic power station according to the annual surface deformation rate. These areas are divided into three categories: small deformation, medium deformation, and large deformation. The steps include: Step 2.1, First-stage clustering of displacement time series based on SOM: (1) Extract the deformation rate values ​​of all pixels in the entire photovoltaic power station area from the annual average deformation rate results generated in step 1 to form a sample set. V={V 1 ,V 2 ,V 3 ，...,V n }， n represents the total number of pixels, which serves as the input data for the first stage of SOM-based displacement time series clustering.

[0022] (2) Initialize the SOM network parameters and set the initial weights. W g Learning rate or Winning Areas s Number of times of learningS .

[0023] (3) Output neurons W g ={w1,w2,…,w g}, g=1,2,3,…, L , L To determine the number of output neurons, calculate the normalized pixel deformation rate value. With all output neurons W g European distance d g The neuron with the smallest distance is selected as the winning neuron.

[0024] ; In the formula, This represents the deformation rate value of the nth normalized pixel. .

[0025] (4) Taking the winning neuron as the center, the winning neighborhood is derived based on the neighborhood radius. Neurons in the winning neighborhood all have the opportunity to adjust their weights.

[0026] (5) The weights of all neurons in the winning neighborhood are adjusted as follows: ; In the formula, This indicates that after the adjustment, the first g The weight vector of each output neuron. With training time t The learning rate changes.

[0027] (6) Determine whether the training count has reached the set number of learning counts or the learning rate has decreased to the set value. or If the condition is met, training stops; otherwise, the next round of clustering begins.

[0028] Step 2.2: Second-stage clustering of surface deformation rate values ​​of photovoltaic power plants based on K-means.

[0029] (1) Input all normalized pixel deformation rate samples Select cluster centers C = {C1, C2, ..., C} Y The cluster number Y (Y=3) is fed into K-means to calculate the deformation rate value of each pixel in the dataset. With cluster center C g The distance is used to assign the deformation rate value of each pixel to the category S belonging to the nearest cluster center. kThree categories were obtained: {S1, S2, S3}, which correspond to three types of regions with small deformation, medium deformation, and large deformation, respectively, according to the magnitude of the surface deformation rate.

[0030] ; in, C is the normalized deformation rate vector of the nth pixel. g For the g-th cluster center, To find the cases where g=1, 2, 3 The minimum value.

[0031] (2) Recalculate the cluster center of each pixel's deformation rate value category. .

[0032] ; In the formula, After recalculation, the first g Cluster centers for each pixel's deformation rate value category. M represents the category of the deformation rate value of the i-th pixel. Sg For category S g Total number of pixel deformation rate values.

[0033] (3) Determine whether the distance between the new cluster center and the original cluster center is less than the set value, or whether the iteration has reached the maximum number of iterations. If the conditions are met, stop the iteration and output the clustering results and cluster centers; otherwise, repeat steps (1) and (2) for iterative calculation.

[0034] (4) Based on the clustering results, the surface deformation state of the photovoltaic power station is divided into small deformation region, medium deformation region and large deformation region.

[0035] Step 3: Based on the clustering results of the deformable areas of the photovoltaic power station, the forward and reverse DEM methods are used to divide the photovoltaic power station into slope units. Combined with the slope unit division, the clustering results of the deformable areas of the photovoltaic power station are further subdivided. The implementation steps include:

[0036] Step 3.1: Based on the original DEM data obtained in Step 1, perform depression filling processing. Iteratively determine each grid cell. If it is a depression point, fill its elevation to the elevation of its lowest outlet neighborhood, thereby generating a hydrologically continuous digital elevation model without false depressions.

[0037] Step 3.2: Extract the preliminary boundary of the slope unit based on the natural water flow characteristics of the terrain, and calculate the slope of the forward DEM according to the formula. ; In the formula z Represents the elevation values ​​of the DEM raster. , They are x Direction (east-west) and y The rate of change of elevation in the north-south direction is calculated using the following formula: ; ; In the formula, cellsize is the raster resolution of the DEM. The elevation of the center grid. It is the elevation of the neighboring rasters surrounding the central raster.

[0038] The D8 algorithm is used to determine the water flow direction in a forward DEM. It is assumed that water will only flow towards the cell with the largest elevation drop among its eight neighboring cells. By comparing the elevation differences between the central cell and its eight neighboring cells, the direction with the largest elevation difference is identified as the water flow direction for that cell. After determining the water flow direction, the catchment area is calculated using the formula.

[0039] ; Indicates the ()th in the DEM e,f The catchment area of ​​the ()th grid; upstream represents the set of all "upstream grids", that is, the area of ​​the terrain where water can naturally flow towards the ()th grid. e,f All adjacent and more distant graticles of the graticles; To count the number of all upstream grids that can be merged into the target grid (each upstream grid is counted as 1, and the sum is the total number of upstream grids).

[0040] Step 3.3: Construct the reverse DEM. Calculate the reverse DEM according to equation (4.5). ; in It is the elevation value of the corresponding raster in the reverse DEM. It is the maximum elevation value of the forward DEM within the photovoltaic power station area. This is the original elevation value of each grid cell in the forward DEM. Using this formula, we can obtain the reverse elevation value corresponding to each grid cell, and thus generate a complete reverse DEM.

[0041] Step 3.4: Based on the reverse DEM obtained in (3), calculate the water flow direction and catchment area of ​​the reverse DEM using the same analysis method as the forward DEM. Finally, divide and optimize the slope units by spatially superimposing the watershed line extracted from the forward DEM (as the top and lateral boundaries of the slope unit) with the slope toe line extracted from the reverse DEM (as the bottom boundary of the slope unit). After superposition, these boundaries will divide the terrain of the photovoltaic power station into independent closed regions, each of which is a preliminary slope unit, thus obtaining the slope unit division vector map of the photovoltaic power station.

[0042] Step 3.5: Based on the clustering results in Step 2, the surface deformation state of the photovoltaic power station is divided into small deformation area, medium deformation area, and large deformation area. Combined with the slope unit division results of the photovoltaic power station, the clustering results of the deformation area of ​​the photovoltaic power station are further subdivided to obtain the result of the surface deformation level classification of the photovoltaic power station based on the slope unit.

[0043] Step 4: Based on the subdivision results of the photovoltaic power station deformation area and combined with the surface deformation monitoring results of the photovoltaic power station, a completely random sampling technique is used to extract the displacement time series of 50 locations within each subdivided deformation unit, and the average displacement time series within the subdivided unit is calculated; the implementation steps include: Step 4.1: Based on the results of the classification of surface deformation levels of photovoltaic power stations obtained in Step 3, extract all SBAS pixels covered inside each slope unit polygon.

[0044] Step 4.2: Assign a unique number from 1 to F to all pixels within each subdivision unit. Using a computer random number generator, generate 50 unique random integers within the range [1, F]. The points corresponding to these 50 random numbers constitute a random sample. The sampling probability is...

[0045] ; Step 4.3: For the 50 sampled points, extract their complete displacement time series. For each SAR image time... ( c =1, 2, ..., N+1 (where N+1 is the number of SAR images), calculate the average displacement value of these 50 points.

[0046] ; In the formula For this subdivided deformation element in The average cumulative displacement at time t. For the first w Each sample point in The cumulative displacement at each moment, arranged in chronological order. This constitutes the average displacement time series of the subdivided deformation unit. Then, for each subdivided deformation unit, the effective deformation points are randomly sampled and the average displacement time series is calculated to obtain the average displacement time series results of all subdivided deformation units.

[0047] Step 5: Using the OVMD adaptive decomposition method, the average displacement time series within each subdivided deformation unit is decomposed to obtain the trend, seasonality, and periodic components of each average displacement time series; the implementation steps include: Step 5.1: Calculate the average displacement time series of all subdivided deformation elements obtained in Step 4. As input, the OVMD (Optimized Variational Mode Decomposition) method is used for adaptive signal decomposition, specifically including: K-value optimization: Determine the optimal range for K and set the balance parameter. a Set to 1 temporarily.

[0048] Step 5.2, Initial VMD Decomposition (1) Constructing variational problems.

[0049] Based on the initial parameter settings (i.e.) and The VMD method is used to perform preliminary decomposition of the input subdivided deformation element time series to obtain... K There are several modal components (IMFs). Assume the displacement time series of the original landslide deformation characteristic points can be decomposed into... K The modal component, then the... k Each modal component is

[0050] ; in: It represents the instantaneous amplitude and is also called the envelope. The phase of the decreasing function; h =1,2,..., K . yes and Pure and harmonic signals have relatively slow changes in amplitude and frequency; t It's time.

[0051] To calculate the bandwidth of each mode, the analytic signals of each intrinsic mode component are obtained using the Hilbert transform. The spectrum of each obtained mode signal is then tuned to the corresponding "baseband," and the bandwidth of the mode is evaluated using the following formula: ; in: j represents the unit impulse function; a Represents the imaginary unit; * represents the convolution operator, ω h This represents the center frequency of the h-th mode; It is an exponential term representing the center frequency on a complex surface.

[0052] By calculating the L2 norm of the aforementioned modal width, for the original input signal... Ultimately, a constrained variational problem can be constructed, expressed as: ; ; in: Indicates the h-th mode; It is the partial derivative of the function; denoted by Dirac distribution, and st indicates that it is constrained by .

[0053] (9) Solve for variational modes.

[0054] 1) By introducing the Lagrange multiplier λ(t) and the penalty parameter α, the unconstrained problem is obtained, and then the optimal solution of the unconstrained problem is found by using the alternating multiplication operator algorithm. and The update formula is as follows:

[0055] ; 2) Minimize

[0056] ; in: , , They are , , Fourier transform.

[0057] (3) Find the optimal value based on the following two criteria. If the criteria are not met, perform iterative calculations until the criteria are met: 1) Judgment Criterion 1: Determine whether there is a moderate or high correlation between different IMF components after decomposition, that is, whether the distance correlation coefficient between different components is greater than 0.4.

[0058] 2) Judgment criterion 2: Determine whether there is mode aliasing between different components after decomposition, that is, whether the center frequency difference (DCF) between any two components is higher than half of the sum of their 3dB bandwidths (HSB).

[0059] Step 5.3: Optimizing the α value. After obtaining the optimal mode decomposition number... Then, the parameters were set using the logarithmic interval method. α The optimization range. Based on the determined optimal parameters. Find all values ​​within the optimization interval. α and The input time series is combined and decomposed according to various combinations. Based on the decomposition results, the approximate entropy, center frequency, 3dB bandwidth, and signal-to-noise ratio (SNR) between the reconstructed and original signals are calculated for each IMF component. The Composite Judgment Index (CDI) corresponding to different values ​​is then calculated. The specific definition of CDI is as follows:

[0060] ; In the formula: This represents the approximate entropy average of all components after all decompositions; It is an empty variable; t i Here are the logical variables for the Mode Aliasing Decision Matrix (DMMA), which is calculated based on the DCF and HSB values ​​between the decomposed IMF components. If mode aliasing exists between the decomposed components, then... ,otherwise When the CDI of the decomposed components reaches its minimum value, the corresponding α value is considered to be the optimal value. It is important to note that when multiple identical minimum CDI values ​​exist simultaneously, the value should be determined based on the maximum SNR between the reconstructed timing sequence and the original timing sequence. This minimizes the information loss of the reconstructed signal.

[0061] Step 5.4 Finally, based on the optimized parameter settings (i.e. and The input subdivided deformation element average displacement time series is decomposed again using VMD. The decomposed data is then processed. K Spectral analysis was performed on each modal component. The modal components were transformed from the time domain to the frequency domain using Fast Fourier Transform (FFT) to obtain their spectral information. Based on the spectral information, the modal components were classified into trend displacements and seasonal displacements.

[0062] (1) Suppose that the OVMD decomposition yields K The modal components are respectively For modes Perform Discrete Fourier Transform (DFT): ; Where T represents the length of the time series. The sampling interval (e.g., SBAS time resolution).

[0063] (2) Calculate the power spectral density (PSD) ; in, Peak frequency This is the dominant frequency of the mode.

[0064] (3) Set a frequency threshold (generally We used 0.1 Hz (i.e., period > 10 sampling points) to distinguish between trend displacement and seasonal displacement. Frequency lower than... The modal components are classified as trend displacements with frequencies higher than [missing information]. The modal components are attributed to seasonal displacements.

[0065] (4) Synthesize the trend component and seasonal component. Sum the classified modal components to obtain the final trend component and initial seasonal component. Trend displacement It is the sum of all modal components belonging to the trend displacement, and the initial seasonal displacement. It is the sum of all modal components belonging to seasonal displacement.

[0066] ; ; In the formula, Indicates a trend of displacement. Indicates the initial seasonal displacement. This represents the i-th modal component. This indicates that the modal component is classified as a trend displacement. This indicates that the modal component is classified as a seasonal displacement.

[0067] (5) Based on the initial seasonal displacement obtained above, repeat the discrete Fourier transform (DFT) and power spectral density (PSD) calculations in (1) and (2) above, identify significant spectral peaks in the PSD, and record their corresponding periods. T iA periodic threshold is set, defining components with a period close to 365 days as strictly seasonal components, and components with a period between 1 and 5 years and no fixed length as periodic components. Finally, based on these thresholds, the various IMF modal components that make up the initial seasonal component are reclassified and synthesized to obtain the final seasonal component. and periodic components .

[0068] ; ; In the formula, Indicates the final seasonal component, Indicates periodic components, This represents the i-th modal component. This indicates that the modal component is classified as a seasonal displacement.

[0069] This indicates that the modal component is classified as a periodic displacement.

[0070] Step 6: Using MIC correlation analysis, analyze the correlation between various regional rainfall characteristic parameters and the trend, seasonal, and periodic components of each mean displacement time series, and propose key rainfall characteristic parameters applicable to various displacement components of photovoltaic power plants; the implementation steps include: Step 6.1, Selection of Displacement Influencing Factors: For trend-based displacement, monitoring time is selected as an influencing factor; for seasonal and periodic displacement, cumulative rainfall over 1, 3, 5, 7, 10, and 14 days, maximum single rainfall, and cumulative effective rainfall over 3, 5, 7, 10, and 14 days are selected as influencing factors. Cumulative and effective rainfall are calculated using equations 7.1 and 7.2. ; In the formula: H The cumulative rainfall over a specified period; n Represents the number of days included in the specified time period. H s For the first s Rainfall for the day.

[0071] ; In the formula: E Effective rainfall; R 0 represents the rainfall for that day; R s For the first sDaily rainfall for the day; q The rainfall infiltration coefficient for the study area is 0.8. This represents the coefficient related to the infiltration of rainfall on the previous day s, used to reflect the weight of the influence of rainfall from different days ago on the infiltration of the current effective rainfall.

[0072] Step 6.2: Perform correlation analysis Based on the selected influencing factors, the mean MIC analysis was used to analyze the correlation between displacement components and influencing factors, assuming that trend displacement, seasonal displacement, or periodic displacement are... V X The corresponding influencing factors are V Y .

[0073] (1) First, V X and V Y Sort in ascending order, then define a x 1 y 1 The grid G1 is used as a pair of partitions, one of which will V X Each sample point is divided into x 1 Another division will V Y Each sample point is divided into y 1 It contains a set of cells, and allows some of those cells to be empty.

[0074] (2) Find the mesh G 1 The probability mass distribution function of all cells in the array And obtain its maximum mutual information value. The eigenma values ​​at this time As shown in the following formula: ; In the formula: To fall within the grid G 1 The number of sample points in the J-th row and I-th column cell. Z The total number of sample points, Indicates that in variable V X The number of columns in a single grid division within the range of values. Indicates that in variable V Y The number of rows in a single grid division within the range of values.

[0075] Due to different grids G 1 This will lead to different Therefore, the global optimum is determined by exhaustive search of the feature matrix. Grid G 0 Then the variable V X and V Y MIC is ; In the formula: The maximum information coefficient is used to measure the information content of two variables. V X and V Y The degree of correlation between them This indicates that the variable is divided into x OK y Under the grid G ​​of the column, the feature matrix value is used to measure the correlation between displacement components and influence factors, and B(N) is the maximum grid area for the search. Indicates a specific grid G The maximum mutual information value under, Indicates the number of grid divisions x and y The smaller value in the equation is taken as the natural logarithm. This indicates that when the globally optimal value is locked... Grid G 0 The corresponding feature matrix value.

[0076] (4) Calculate the average MIC Sum the MIC values ​​of all variable pairs and then divide by the total number of variable pairs to obtain the average MIC value.

[0077] (5) Based on the calculated average MIC value, select the MIC correlation coefficient threshold as a reference value, delete the influence factors below the threshold, and retain the influence factors above the threshold as the key rainfall characteristic parameters of various displacement components of photovoltaic power stations.

[0078] Step 7: Based on the proposed key rainfall characteristic parameters, construct a photovoltaic power station deformation unit displacement combination prediction model based on the DHKELM method to achieve accurate prediction of deformation in various areas of the photovoltaic power station; the implementation steps include: Step 7.1: Constructing the model input dataset: (1) For trend displacement, the input feature is the monitoring time series T; the target variable is the trend displacement sequence extracted in step five; for seasonal / periodic displacement, the input feature is the key rainfall characteristic parameters screened in step six, and the target variable is the seasonal displacement sequence or periodic displacement sequence extracted in step five. Align the above features with the target variable and construct the original dataset according to time points. D ={ U u , V v},in U u For the input feature vector, V v This represents the corresponding displacement value.

[0079] (2) In order to eliminate the influence of data dimensions, various influencing factors and periodic displacements are normalized to the interval [0, l] according to different dimensions.

[0080] Input feature normalization: ; in, u u Represents the original values ​​of the characteristic variables. and These are the minimum and maximum values ​​of the feature, respectively.

[0081] Target variable normalization: ; in, v v Represents the original value of the target variable. and Find the minimum and maximum values ​​of the target variable, respectively.

[0082] Normalized dataset D The training set was divided chronologically in an 8:2 ratio. D train and prediction set D test .

[0083] Step 7.2: Train the DHKELM model: For each training set Perform the following training steps respectively: 1) Randomly initialize the input weights w in the E-th (E=1,2,...,u,u is the number of ELM-AEs) hidden layer. E and bias b E .

[0084] 2) By using the output feature matrix of the previous ELM-AE hidden layer as the input of the next ELM-AE hidden layer, the output feature matrix H of the E-th ELM-AE hidden layer is calculated sequentially. E and output weights .

[0085] ; ; In the formula: f(·) is the activation function; I is the identity matrix; Here, X represents the regularization coefficient, and X represents the original data matrix containing the influence factor features input to the E-th ELM-AE. H represents the output feature matrix of the E-th ELM-AE hidden layer. E The transpose of the matrix, This represents the output feature matrix of the (E-1)th hidden layer.

[0086] 3) The output feature matrix H of the last ELM-AE hidden layer u Used as the decision input for HKELM, and the output weights of HKELM are calculated. : ; In the formula G represents the output feature matrix of the last ELM-AE hidden layer, and G is the mixing kernel function matrix: ; Where μ is the mixing coefficient (0≤μ≤1), and These are the linear kernel and the radial basis kernel functions, respectively.

[0087] Step 7.3: For the prediction set D test Predicting from the samples in the dataset: For prediction set D test Perform forward propagation to obtain its deep features H test .make For E = 1 to u:

[0088] ; in, Represents the prediction set After the initial stage of forward propagation, the feature matrix is ​​given by E, where E is the index of the ELM-AE and u is the number of ELM-AEs. Indicating on the prediction set During the forward propagation process to acquire deep features, the deep feature matrix of the i-th layer is given by f(·), where f(·) is the activation function. w represents the depth feature matrix of the (i-1)th layer. E For the input weights, b E This represents the bias. Through the aforementioned forward propagation process, the prediction set is obtained. Deep features after multi-layer processing ;

[0089] Calculate the test kernel matrix G test ； and These are the linear kernel and the radial basis kernel functions, respectively.

[0090] Calculate the final predicted value : ; Predicted values Inverse normalization yields predictions at the actual scale: ; in v max and v min These are the minimum and maximum values ​​of the original target variable before normalization.

[0091] Step 7.4: Based on the above steps, predict the trend, seasonal and periodic components of each subdivided deformation unit, and obtain the prediction results of the deformation of each area of ​​the photovoltaic power station by adding the trend term, seasonal term and periodic displacement.

Claims

1. An integrated method for monitoring and predicting surface deformation in photovoltaic power station areas, characterized in that: The method includes: Step 1: Collect Sentinel-1 C-band SAR images, DEM data, and precise orbit data for the photovoltaic power station during the specified monitoring period to obtain the surface deformation monitoring results of the photovoltaic power station; Step 2: Based on the obtained surface deformation monitoring results, the Som-Kmeans two-step clustering method is used to divide the surface deformation area of ​​the photovoltaic power station according to the annual surface deformation rate, and the area is divided into three categories: small deformation, medium deformation and large deformation. Step 3: Based on the clustering results of the deformation area of ​​the photovoltaic power station, the forward and reverse DEM methods are used to divide the photovoltaic power station into slope units, and the clustering results of the deformation area of ​​the photovoltaic power station are further subdivided to obtain the result of the classification of the surface deformation level of the photovoltaic power station based on the slope units. Step 4: Based on the results of the subdivision of the deformation area of ​​the photovoltaic power station, and combined with the monitoring results of the surface deformation of the photovoltaic power station, the displacement time series of 50 locations in each subdivided deformation unit is extracted using a completely random sampling technique, and the average displacement time series in the subdivided unit is calculated. Step 5: Use the OVMD adaptive decomposition method to decompose the average displacement time series in each subdivided deformation unit to obtain the trend, seasonal and periodic components of each average displacement time series. Step 6: Use MIC correlation analysis technology to analyze the correlation between various rainfall characteristic parameters in the region and the trend, seasonal and periodic components of each average displacement time series, and obtain the key rainfall characteristic parameters for various displacement components of photovoltaic power stations. Step 7: Based on key rainfall characteristic parameters, construct a photovoltaic power station deformation unit displacement combination prediction model with the DHKELM method as the core to achieve accurate prediction of deformation in various areas of the photovoltaic power station.

2. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 1, characterized in that: Methods for obtaining surface deformation monitoring results for photovoltaic power plants include: Step 1.1: Obtain Sentinel-1A single-view complex SLC data, precise orbit data, and reference DEM data for the study area during the study period from a public platform; Step 1.2: Obtain surface deformation information and obtain the surface deformation monitoring results of the photovoltaic power station; the steps include: Step 1.2.1: The surface deformation information of the study area is obtained by processing the small baseline subset method, and the N+1 Sentinel-1A single-look complex data SAR images are sorted in chronological order: Select and register the super master image to generate M interferograms, where M satisfies the following conditions: ； Step 1.2.2, in The two SAR images acquired at time 1 generate the first j Interferogram, with flat and topographic phases removed respectively, azimuth coordinate x y The differential interference phase at the range coordinate r pixel can be expressed by the following formula. ； In the formula, The center wavelength of the radar. and They are respectively t b and t a Time relative to t 0 Cumulative deformation of radar line-of-sight direction at any given time. and They are respectively derived from deformation variables and The resulting change in phase deformation value; Steps 1, 2, and 3, and any interference pattern in equation (2.2) j The deformation phase is represented by the average rate over the corresponding time period. v j This means, that is: ； Equation (2.2) can be expressed as: ； This represents the end time corresponding to the j-th interferogram. This represents the starting time corresponding to the j-th interferogram. v j Represents any interferogram j The average rate of deformation phase within the corresponding time period, This represents the deformation phase of the j-th interferogram; Equation (2.4) can be written in matrix form: ； A is an M×N coefficient matrix. v The deformation rate vector, Let A be the matrix composed of the phases of M differential interferograms; when M≥N, the rank of the coefficient matrix A is N, which can be solved using the least squares criterion according to equation (2.6); ； in Let A be the transpose of the coefficient matrix A. For matrix The inverse matrix; Find Then, the deformation phase, deformation amount, and corresponding deformation rate vector for each imaging time period are calculated. ); When M < N, the equation has infinitely many solutions. The minimum norm solution of the deformation rate vector is obtained by using the singular value decomposition method. The rate is integrated over each time period to obtain the deformation of each time period, thus obtaining the time series of surface displacement of the photovoltaic power station. Solve equation (2.5) using the singular value decomposition method. ； In the formula, U is A T An M×M orthogonal matrix within A; R is an M×M matrix within A. T An M×N orthogonal matrix in A; A T The diagonal matrix of elements on the inner diagonal of A is C. Let A be the transpose of matrix R; let the rank of A be H, then A T The first H eigenvalues ​​of A are not 0, and the remaining eigenvalues ​​are all 0. Therefore, the above formula is changed to: ； In the formula, It is the l-th element in C; , Let l be the l-th element in U and R; Therefore, the estimated value is obtained. Replace at average rate The formula is ； In the formula, G is a matrix where G(j′, B) is 0. v G is the deformation rate vector; G is calculated using singular value decomposition to obtain... The minimum norm solution is the deformation rate vector. Assume the duration of each time period is... Then the deformation of each time period Through deformation rate vector Duration of Time Period The product of these two is: ； Based on the deformation of the photovoltaic power station at different time periods, a time series of surface displacement of the photovoltaic power station is constructed, assuming the initial displacement is... Then the displacement d at the k-th time point k It is the sum of the deformation variables of the first k time periods, i.e. 。 3. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 1, characterized in that: The method for dividing the region into three categories—small deformation, medium deformation, and large deformation—as described in step 2 includes: Step 2.1, First-stage clustering of displacement time series based on SOM, including: Step 2.1.1: Extract the deformation rate values ​​of all pixels within the entire photovoltaic power station area from the annual average deformation rate results of the entire power station, and form a sample set. V={V 1 ,V 2 ,V 3 ，...,V n }， n is the total number of pixels, which serves as the input data for the first stage of SOM-based displacement time series clustering. Step 2.1.2: Initialize SOM network parameters and set initial weights. W g Learning rate η Winning Areas σ and number of times of learning S ; Step 2.1.3, Output Neuron W g ={w1,w2,…,w g }, g=1,2,3,…, L , L To determine the number of output neurons, calculate the normalized pixel deformation rate value. With all output neurons W g European distance d g The neuron with the smallest distance is selected as the winning neuron; ； In the formula, This represents the deformation rate value of the nth normalized pixel. ; Step 2.1.4: Using the winning neuron as the center, determine the winning neighborhood based on the neighborhood radius. Neurons within the winning neighborhood all have the opportunity to adjust their weights. Step 2.1.5: Adjust the weights of all neurons in the winning neighborhood: ； In the formula, This represents the adjusted weight vector of the g-th output neuron. The learning rate as a function of training time t; Step 2.1.6: Determine whether the training iterations have reached the set number of learning iterations or the learning rate has decreased to the set threshold. η If the condition is met, training stops; otherwise, the next round of clustering is performed. Step 2.2: Second-stage clustering of surface deformation rate values ​​for photovoltaic power plants based on K-means, including: Step 2.2.1: Input all normalized pixel deformation rate samples. Select cluster centers C = {C1, C2, ..., C} Y The cluster number Y (Y=3) is fed into K-means to calculate the deformation rate value of each pixel in the dataset. With cluster center C g The distance is used to assign the deformation rate value of each pixel to the category S belonging to the nearest cluster center. k Three categories were obtained: {S1, S2, S3}, which correspond to three types of regions with small deformation, medium deformation, and large deformation, respectively, according to the magnitude of the surface deformation rate. Step 2.2.2: Recalculate the cluster center for each pixel's deformation rate value category. ; Step 2.2.3: Determine whether the distance between the new cluster center and the original cluster center is less than the set value, or whether the iteration has reached the maximum number of iterations. If the conditions are met, stop the iteration and output the clustering result and cluster center; otherwise, repeat the iterative calculation. Step 2.2.4: Based on the clustering results, the surface deformation state of the photovoltaic power station is divided into small deformation region, medium deformation region and large deformation region.

4. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 1, characterized in that: Methods for obtaining the classification of surface deformation levels for photovoltaic power plants based on slope elements include: Step 3.1: Based on the original DEM data, perform depression filling processing, iteratively determine each grid cell, and if it is a depression point, fill the elevation to the elevation of the lowest outlet neighborhood, thereby generating a digital elevation model without false depressions and continuous in hydrology. Step 3.2: Extract the preliminary boundary of the slope unit based on the natural water flow characteristics of the terrain, and calculate the slope of the forward DEM; Step 3.3: Construct and calculate the reverse DEM; Step 3.4: Following the same analysis method as the forward DEM, calculate the water flow direction and catchment area of ​​the reverse DEM. Finally, perform slope unit division and optimization, including spatially superimposing the watershed line extracted from the forward DEM with the slope toe line extracted from the reverse DEM. After superposition, these boundaries divide the terrain of the photovoltaic power station into independent closed areas. Each closed area is a preliminary slope unit, thus obtaining the slope unit division vector map of the photovoltaic power station. Step 3.5: Based on the slope unit division results of the photovoltaic power station, further subdivide the clustering results of the deformation area of ​​the photovoltaic power station to obtain the result of the surface deformation level classification of the photovoltaic power station based on the slope unit.

5. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 1, characterized in that: Step 4 includes the following steps: Step 4.1: Classification of surface deformation levels of photovoltaic power stations based on slope units. For each slope unit polygon, extract all SBAS pixels covered inside. Step 4.2: Assign a unique number from 1 to F to all pixels in each subdivision unit. Use a computer random number generator to generate 50 non-repeating random integers in the range [1, F]. The points corresponding to these 50 random numbers constitute a random sample. Step 4.3: For the 50 sampled points, extract their complete displacement time series. For each SAR image time... , c =1, 2, ..., N+1, where N+1 is the number of SAR images, calculate the average displacement value of these 50 points; ； To subdivide deformation elements in The average cumulative displacement at time t. For the first w Each sample point in The cumulative displacement at each moment, arranged in chronological order. The average displacement time series constituting this subdivided deformation element { Then, for each subdivided deformation unit, the effective deformation points are randomly sampled and the average displacement time series is calculated to obtain the average displacement time series results of all subdivided deformation units.

6. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 1, characterized in that: Methods for obtaining the trend, seasonality, and periodic components of each average displacement time series include: Step 5.1: Calculate the average displacement time series of all subdivided deformation elements obtained in Step 4. As input, the OVMD optimized variational mode decomposition method is used for adaptive signal decomposition; Step 5.2: Perform initial VMD decomposition; Step 5.3: Optimize the α value to obtain the optimal mode decomposition number. Then, the parameters were set using the logarithmic interval method. α The optimization range is determined based on the optimal parameters. Find all within the optimization interval α and The input time series is combined and decomposed according to various combinations. Based on the decomposition results, the approximate entropy, center frequency, 3dB bandwidth, and signal-to-noise ratio (SNR) between the reconstructed signal and the original signal are calculated for each IMF component. Then, the comprehensive judgment index (CDI) corresponding to different values ​​is calculated. ； In the formula: This represents the approximate entropy average of all components after all decompositions; It is an empty variable; t i The logical variable for the Mode Aliasing Decision Matrix (DMMA) is calculated based on the DCF and HSB values ​​between the decomposed IMF components. If mode aliasing exists between the decomposed components, then... ,otherwise The optimal α value is the one where the CDI of the decomposed components reaches its minimum. When there are more than one identical minimum CDI value, the maximum SNR between the reconstructed timing sequence and the original timing sequence shall be used to determine the minimum CDI value. Step 5.4: Set the parameters according to the optimized parameters. and The input subdivided deformation element average displacement time series is decomposed again using VMD. The decomposed data is then processed. K Spectral analysis was performed on each modal component. The modal components were transformed from the time domain to the frequency domain using Fast Fourier Transform (FFT) to obtain spectral information. Based on the spectral information, the modal components were divided into trend displacements and seasonal displacements.

7. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 6, characterized in that: The implementation methods for step 5.4 include: Step 5.4.1: Assume that the decomposition obtained by OVMD... K The modal components are respectively For modes Perform Discrete Fourier Transform (DFT): ； T represents the length of the time series. The sampling interval; Step 5.4.2: Calculate the power spectral density (PSD): ； Peak frequency The dominant frequency of the mode; Step 5.4.3: Set a frequency threshold. To classify trend displacements and seasonal displacements; frequency lower than The modal components are classified as trend displacements with frequencies higher than [missing information]. The modal components are attributed to seasonal displacements; Step 5.4.4: Synthesize the trend component and seasonal component. Summate the classified modal components to obtain the final trend component and initial seasonal component, and the trend displacement. It is the sum of all modal components belonging to the trend displacement, and the initial seasonal displacement. It is the sum of all modal components belonging to seasonal displacement; ； ； Indicates a trend of displacement. Indicates the initial seasonal displacement. This represents the i-th modal component. This indicates that the modal component is classified as a trend displacement. This indicates that the modal component is classified as a seasonal displacement category; Step 5.4.5: Repeat the Discrete Fourier Transform (DFT) and Power Spectral Density (PSD) calculations, identify significant spectral peaks in the PSD, and record the corresponding periods. T i A periodic threshold is set, and the IMF modal components that make up the initial seasonal component are reclassified and synthesized according to the periodic threshold to obtain the final seasonal component. and periodic components ; ； ； In the formula, Indicates the final seasonal component, Indicates periodic components, This represents the i-th modal component. This indicates that the i-th modal component is classified as a seasonal displacement category. This indicates that the i-th modal component is classified as a periodic displacement.

8. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 1, characterized in that: Step 6 describes the method for obtaining key rainfall characteristic parameters for various displacement components of photovoltaic power plants, which includes: Step 6.1: For trend-based displacement, the monitoring time is selected as an influencing factor; for seasonal and periodic displacement, the cumulative rainfall over 1, 3, 5, 7, 10, and 14 days, the maximum single rainfall, and the cumulative effective rainfall over 3, 5, 7, 10, and 14 days are selected as influencing factors. The cumulative rainfall and effective rainfall are then calculated. ； In the formula: H The cumulative rainfall over a specified period; n Represents the number of days included in the specified time period. H s For the first s Rainfall for the day; ； In the formula: E Effective rainfall; R 0 represents the rainfall for that day; R s For the first s Daily rainfall for the day; q The rainfall infiltration coefficient for the study area is 0.

8. This represents the coefficient related to the infiltration of rainfall on the previous s-th day, used to reflect the weight of the infiltration impact of rainfall from different days ago on the current effective rainfall. Step 6.2, Correlation Analysis: The correlation between displacement components and influencing factors is analyzed using the mean MIC (Minimum Interval) method. Let the displacement be trend-based, seasonal, or periodic. V X The corresponding influencing factors are V Y ; First of all V X and V Y Sort in ascending order, then define a x 1 y 1 The grid G1 is used as a pair of partitions, one of which will V X Each sample point is divided into x 1 Another division will V Y Each sample point is divided into y 1 A set of cells, and allows some of those cells to be empty; Find the mesh G 1 The probability mass distribution function of all cells in the array And obtain the maximum mutual information value. The eigenma values ​​at this time As shown in the following formula: ； In the formula: To fall within the grid G 1 The number of sample points in the J-th row and I-th column cell. Z The total number of sample points, Indicates that in variable V X The number of columns in a single grid division within the range of values. This represents the number of rows in a single mesh division within the range of values ​​for variable VY; Due to different grids G 1 This will lead to different Therefore, the global optimum is determined by exhaustive search of the feature matrix. Grid G 0 Then the variable V X and V Y MIC is ； In the formula: For the maximum information coefficient, This represents dividing the variables into a grid G ​​with x rows and y columns, where B(N) is the maximum grid area for the search. This represents the maximum mutual information value under a specific grid G. This indicates taking the natural logarithm of the smaller of the grid division numbers x and y. This indicates that when the globally optimal value is locked... The feature matrix value corresponding to grid G0; Step 6.3: Sum the MIC values ​​of all variable pairs, then divide by the total number of variable pairs to obtain the average MIC value; Step 6.4: Based on the average MIC value, select the MIC correlation coefficient threshold as a reference value, delete the influencing factors below the threshold, and retain the influencing factors above the threshold as the key rainfall characteristic parameters of various displacement components of the photovoltaic power station.

9. The integrated method for monitoring and predicting surface deformation in photovoltaic power station areas according to claim 1, characterized in that: Step 7 includes the following steps: Step 7.1: Constructing the Model Input Dataset: For trend displacement, the input feature is the monitoring time series T; the target variable is the trend displacement sequence. For seasonal / periodic displacement, the input feature is key rainfall characteristic parameters, and the target variable is the seasonal displacement sequence or the periodic displacement sequence. Align the above features with the target variable and construct the original dataset according to time points. D ={ U u , V v },in U u For the input feature vector, V v The corresponding displacement values ​​are defined; various influencing factors and periodic displacements are normalized to the [0, l] interval according to different dimensions; the normalized dataset is then processed. D The training set was divided chronologically in an 8:2 ratio. D train and prediction set D test Step 7.2, Train the DHKELM model: For each training set Perform the following training steps respectively: Randomly initialize the input weights w in the E-th (E=1,2,...,u,u is the number of ELM-AEs) hidden layer. E and bias b E ; By using the output feature matrix of the previous ELM-AE hidden layer as the input of the next ELM-AE hidden layer, the output feature matrix H of the E-th ELM-AE hidden layer is calculated sequentially. E and output weights ; ； ； In the formula: f(·) is the activation function; I is the identity matrix; Here, X represents the regularization coefficient, and X represents the original data matrix containing the influence factor features input to the E-th ELM-AE. H represents the output feature matrix of the E-th ELM-AE hidden layer. E The transpose of the matrix, This represents the output feature matrix of the (E-1)th hidden layer; The output feature matrix Hu of the last ELM-AE hidden layer is used as the decision input of HKELM, and the output weights of HKELM are calculated. ; ； In the formula G represents the output feature matrix of the last ELM-AE hidden layer, and G is the mixing kernel function matrix: ； Where μ is the mixing coefficient (0≤μ≤1), and These are the linear kernel and the radial basis kernel functions, respectively. Step 7.3: Make predictions for the samples in the prediction set Dtest: For prediction set D test Perform forward propagation to obtain deep features H test ,make For E = 1 to u: ； in, Represents the prediction set After the initial stage of forward propagation, the feature matrix is ​​given by E, where E is the index of the ELM-AE and u is the number of ELM-AEs. Indicating on the prediction set During the forward propagation process to acquire deep features, the deep feature matrix of the i-th layer is given by f(·), where f(·) is the activation function. w represents the depth feature matrix of the (i-1)th layer. E For the input weights, and b E Indicates bias; Calculate the test kernel matrix G test : ； and These are the linear kernel and the radial basis kernel functions, respectively. Calculate the final predicted value : ； Predicted values Inverse normalization yields predictions at the actual scale: ； in v max and v min These are the minimum and maximum values ​​of the original target variable before normalization; Step 7.4: Predict the trend, seasonal and periodic components of each subdivided deformation unit, and obtain the prediction results of the deformation of each area of ​​the photovoltaic power station by adding the trend term, seasonal term and periodic displacement.

Images