A method for estimating snow water equivalent change based on inversion error

Abstract

The application discloses a snow water equivalent change estimation method based on inversion error, comprising the following steps: S1, obtaining a precipitation data set, performing difference value processing on the precipitation data set, and obtaining a plurality of difference value data sets; S2, filtering all the difference value data sets through snow cover data with the same space-time resolution to obtain filtered difference value data sets; S3, fusing all the filtered difference value data sets through a triangular cap method based on Bayes to obtain a fusion data set; and S4, estimating the snow water equivalent change through the fusion data set. Compared with the prior art, the application simplifies the snow water equivalent change estimation process, estimates through the difference value between relatively mature and reliable precipitation data sets and soil water data sets without considering the snow depth and snow density. Compared with the original estimation path, the spatial resolution is higher, reaching 0.1*0.1. The application provides a completely independent snow water equivalent change numerical reference for the product produced through the original technical path.

Description

Technical Field

[0001] This invention belongs to the field of snow water equivalent prediction technology, specifically relating to a method for estimating snow water equivalent changes based on inversion error. Background Technology

[0002] Snow water equivalent (SWE) is an important parameter for assessing snowmelt resources, and SWE estimation is a crucial aspect of snow cover research in high-altitude and cold mountainous areas. Basic methods for SWE estimation include in-situ measurement, microwave remote sensing estimation, reanalysis, and surface modeling. Based on these methods, SWE reconstruction and empirical regression methods have been developed.

[0003] In-situ measurements are relatively direct and reliable, including both direct and indirect measurements. Direct measurements measure snow column mass to determine snow radiation emission intensity (SWE), while indirect measurements are based on the attenuation of gamma radiation by water molecules. By analyzing the difference in gamma radiation levels before and after snowfall, and using optimization coefficients, the difference can be converted into real-time SWE. However, because in-situ measured SWE is insufficient to characterize the average value within its region, and because snow cover exhibits spatial heterogeneity, especially in complex terrain, this method is not suitable for obtaining SWE values ​​efficiently over large areas. Furthermore, the complex terrain of high-altitude and cold mountainous areas makes it difficult for in-situ observations to cover the entire region, resulting in missing measurement areas and increasing data uncertainty.

[0004] Microwave remote sensing methods developed over the past three decades are used to estimate snowfall precipitation (SWE). This method utilizes the relationship between electromagnetic radiation and snow cover characteristics, estimating snowfall based on the difference in light temperature (TB) at different frequencies. It has advantages such as global coverage, high temporal resolution, and long estimation history. However, it also has drawbacks, including potentially large measurement errors in mountainous and complex terrain; estimates limited to shallow snow layers (<150 mm); and low spatial resolution, thus failing to accurately capture spatial variations caused by mixed land cover and topography. Some studies have proposed improving the quality of satellite-based SWE estimation over large areas by assimilating satellite data with ground-based observation data.

[0005] The SWE derived from reanalysis and surface modeling relies on meteorological data-driven snow cover models, which can simulate complex snowmelt processes in cold environments. However, these models require multiple atmospheric variables as input, while meteorological observation data for high-altitude mountainous areas is scarce. Therefore, due to inaccurate input parameters, the estimation results in high-altitude mountainous areas have not significantly improved compared to other methods. Furthermore, this method is data-intensive and time-consuming when applied on a large scale.

[0006] In addition to the methods mentioned above, other studies have developed SWE reconstruction algorithms and empirical regression algorithms, which estimate SWE based on remotely sensed snow area or snow depth.

[0007] The SWE reconstruction algorithm utilizes snow cover area and the energy required for pure snowmelt to backtrack and integrate the SWE from the snowmelt date to the date the SWE reaches its maximum value, calculating the SWE for the snow season both temporally and spatially. This reconstruction method can provide spatiotemporal patterns of SWE in large mountainous areas. The algorithm simulates the peak snow cover from the snowmelt period, therefore it can only be evaluated after snowfall and is only valid during the snowmelt period. Because the algorithm does not consider precipitation during the meltwater season, the estimation results are biased. Furthermore, the accuracy of this algorithm depends on the accuracy of many factors, including snow cover data, snowfall date, SWE peak value, and the snowmelt estimation method.

[0008] In recent years, many empirical regression models for snow depth-SWE have been developed. These regression-based methods require paired snow depth-SWE ground measurements to calibrate the parameters, and then use them to estimate SWE with only snow depth measurements available. This greatly simplifies SWE calculations and has led to their widespread application in regional studies. However, the SWE value at a specific time is estimated independently for each snow depth value, independent of previous snow depth values. Therefore, if the regression method directly converts snow depth time series to SWE without considering the settlement of new snow or the compaction of snowdrifts, the temporal evolution of SWE will be incorrect. This problem has been addressed by considering the temporal evolution of snow depth in the conversion model, but to date, no method for estimating snow density has been generalized to large areas or globally. Summary of the Invention

[0009] To address the aforementioned shortcomings in existing technologies, this invention provides a method for estimating snow water equivalent changes based on inversion error, which solves the problem of large errors in the conversion of snow water equivalent to snow density in existing technologies.

[0010] To achieve the aforementioned objectives, the present invention employs the following technical solution: a method for estimating snow water equivalent variation based on inversion error, comprising the following steps:

[0011] S1. Obtain the precipitation dataset, perform interpolation on the precipitation dataset, and obtain several interpolation datasets;

[0012] S2. Filter all difference datasets using snow cover data with the same spatiotemporal resolution to obtain the filtered difference dataset;

[0013] S3. The filtered difference datasets are merged using the Bayesian triangular hat method to obtain the merged dataset.

[0014] S4. Estimate the change in snow water equivalent by fusing datasets.

[0015] Further: In S1, the method for processing the difference in the precipitation dataset is specifically as follows:

[0016] All precipitation datasets are interpolated to the same spatiotemporal resolution and equivalent units. The precipitation datasets include a first precipitation dataset based on direct precipitation measurements and a second precipitation dataset based on remote sensing soil moisture observations. The difference dataset is obtained by subtracting the values ​​of the second precipitation dataset from the values ​​of the first precipitation dataset.

[0017] Further: In S2, the method for filtering all difference datasets using snow cover data with the same spatiotemporal resolution is specifically as follows:

[0018] The spatiotemporal difference data is filtered based on the snow cover of the grid, and the specific expression is as follows:

[0019]

[0020] In the formula, S(x,y,t) represents the snow cover at time t and spatial location (x,y), with a value range of [0,100], D(x,y,t) represents the difference data at time t and spatial location (x,y), and D′(x,y,t) represents the filtered difference data at time t and spatial location (x,y).

[0021] Furthermore: In S3, the method for fusing all filtered difference datasets using the Bayesian triangular hat method is specifically as follows:

[0022] S31. Calculate the probability density function for all difference datasets, where the i-th difference dataset NS i The probability density function p(NS) i NS t Specifically:

[0023]

[0024] In the formula, NS t ε represents the true value of the difference dataset. i Let σ be zero-mean white noise in the i-th difference dataset. i Let L(-) be the error covariance of the product of the i-th difference dataset, and L(-) be the likelihood function.

[0025] ε i =NS i -NS t

[0026] S32. Calculate the maximum probability maxL(NS) of the true value in the difference dataset. t NS i NS j (, its expression is:

[0027]

[0028] In the formula, NS j For the j-th difference dataset;

[0029] S33, In order to obtain the true value NS of the difference dataset t The maximum likelihood value is given by defining the cost function J, which is expressed as follows:

[0030]

[0031] S34. Set the first change in the cost function to zero to obtain the true value NS of the difference dataset. t Its specific expression is:

[0032]

[0033] S35, Define NS t =w i NS i +w j NS j Then define the weight w i and w j The specific expression is:

[0034]

[0035] For multiple difference datasets, the true value NS of the difference datasets t =w1NS1+…+w N NS N The weight w is defined for each difference dataset. k This is obtained by minimizing a similar cost function;

[0036]

[0037] For each difference dataset, a weight w is defined. k Error covariance σ i It is obtained by solving the relative uncertainty of each difference dataset through the triangular hat method.

[0038] Further: In S35, the method for solving the relative uncertainty of each difference dataset using the triangular hat method is specifically as follows:

[0039] S351. Randomly select one dataset from N difference datasets and calculate the difference Y between it and the remaining difference datasets. i,M Its specific expression is:

[0040] Y i,M =NS i -NS R =ε i -ε R

[0041] In the formula, NS R Y is a difference dataset randomly selected from N difference datasets. i,M Let S be a difference matrix with (N-1) rows and M columns, i = 1, 2, ..., N-1, where M is the total number of time samples in each NS product. The expression for its covariance matrix S is as follows:

[0042] S = cov(Y) i,M )

[0043] In the formula, cov(·) is the covariance operator, introducing an unknown N*N noise covariance matrix R, whose relationship with S is:

[0044]

[0045] In the formula, J N-1,N For the same matrix, its specific expression is:

[0046]

[0047] The definition of matrix R is:

[0048]

[0049] In the formula, σ ij =cov(ε i , ε j ), i=1,…,N, j=1,…,N;

[0050] S352. Establish the objective function of the Kuhn-Tucker theorem, minimize the objective function under constraints to obtain a set of free parameter solutions, and obtain the relative uncertainty of each difference dataset.

[0051] The objective function of the Kuhn-Tucker theorem is as follows:

[0052]

[0053] In the formula, The constraint functions are as follows:

[0054]

[0055] In the formula, Q is a diagonal matrix.

[0056] The beneficial effects of this invention are as follows: This invention provides a method for estimating the change in snow water equivalent based on inversion error. Addressing the characteristic that errors between products continuously evolve over time, this invention introduces snow cover data based on optical remote sensing as a basis for determining whether there is snow accumulation in a region and excluding the difference portion not caused by snow accumulation. Furthermore, it uses a Bayesian-based generalized triangular hat fusion method to fuse all datasets to extract the difference portion belonging to the change in snow water equivalent. Compared with existing technologies, this invention has the following advantages:

[0057] (1) Simplify the process of estimating snow water equivalent changes by estimating the difference between relatively mature and reliable precipitation datasets and soil water datasets without having to consider snow depth and snow density.

[0058] (2) Compared with the original estimation path, the spatial resolution is higher, reaching 0.1°*0.1°.

[0059] (3) Provide a completely independent numerical reference for the snow water equivalent variation of products produced through the original technical route. Attached Figure Description

[0060] Figure 1 This is a flowchart of a method for estimating the change in snow water equivalent based on inversion error according to the present invention. Detailed Implementation

[0061] The specific embodiments of the present invention are described below to enable those skilled in the art to understand the present invention. However, it should be understood that the present invention is not limited to the scope of the specific embodiments. For those skilled in the art, various changes are obvious as long as they are within the spirit and scope of the present invention as defined and determined by the appended claims. All inventions utilizing the concept of the present invention are protected.

[0062] like Figure 1 As shown, in one embodiment of the present invention, a method for estimating the change in snowmelt equivalent based on inversion error includes the following steps:

[0063] S1. Obtain the precipitation dataset, perform interpolation on the precipitation dataset, and obtain several interpolation datasets;

[0064] S2. Filter all difference datasets using snow cover data with the same spatiotemporal resolution to obtain the filtered difference dataset;

[0065] S3. The filtered difference datasets are merged using the Bayesian triangular hat method to obtain the merged dataset.

[0066] S4. Estimate the change in snow water equivalent by fusing datasets.

[0067] In this embodiment, the present invention is based on the coupling relationship between snow cover and soil hydrology. During winter, snow covers the surface, and soil moisture remains moist with low flux, resulting in minimal evaporation or water replenishment. At this time, the bottom-up inversion method for precipitation based on soil moisture fails to detect precipitation signals due to insignificant changes in soil moisture, while the top-down inversion method can monitor precipitation processes. In spring, as snow melts, the soil transitions from a moist, low-flux state to a moist, high-flux state, which continues until the soil is dry or spring precipitation continues to moisten it. At this point, due to increased soil moisture, the bottom-up inversion method identifies this soil wetting as precipitation, while the top-down inversion method fails to detect precipitation. Therefore, during the snow cover and snowmelt periods, the inversion errors of precipitation products from these two different inversion methods can be used to estimate changes in snowmelt equivalent in high-altitude mountainous areas.

[0068] In S1, the method for processing the difference in the precipitation dataset is as follows:

[0069] All precipitation datasets are interpolated to the same spatiotemporal resolution and equivalent units. The precipitation datasets include a first precipitation dataset based on direct precipitation measurements and a second precipitation dataset based on remote sensing soil moisture observations. The difference dataset is obtained by subtracting the values ​​of the second precipitation dataset from the values ​​of the first precipitation dataset.

[0070] In this embodiment, the same spatiotemporal resolution ensures the stability and rationality of operations between different datasets. This embodiment uses mm as the unit to directly represent snow water equivalent without requiring parameter transformations based on physical or meteorological assumptions. Utilizing the first precipitation dataset, which directly measures precipitation and has the ability to observe snowfall, and the second precipitation dataset, based on remote sensing soil moisture observations, which lacks timely and direct observation capabilities for snowfall, their difference can be considered as the amount of water resources covering the surface in a spatiotemporal environment with snow cover, i.e., snow water equivalent. This embodiment collects three sets of difference datasets to ensure the smooth fusion of the difference datasets.

[0071] In step S2, the method for filtering all difference datasets using snow cover data with the same spatiotemporal resolution is as follows:

[0072] The spatiotemporal difference data is filtered based on the snow cover of the grid, and the specific expression is as follows:

[0073]

[0074] In the formula, S(x,y,t) represents the snow cover at time t and spatial location (x,y), with a value range of [0,100], D(x,y,t) represents the difference data at time t and spatial location (x,y), and D′(x,y,t) represents the filtered difference data at time t and spatial location (x,y).

[0075] In this embodiment, a grid is established for the difference dataset. When the snow cover of the grid is 100, the difference data in the corresponding space-time is completely preserved. When there is snow cover in the grid but the snow cover rate is less than 100 (i.e., the snow cover rate is 0-100), the difference data is multiplied by the proportion of snow cover. When there is no snow cover in the grid (i.e., the snow cover rate is 0), the difference data in the corresponding space-time is deleted and set to 0. The specific principle of filtering the difference dataset is as follows:

[0076] When the entire grid cell is covered by snow, the estimated difference between the datasets is considered to be due to the snow cover, and this value can be regarded as the snow water equivalent. When a portion of the grid cell contains snow, the difference is multiplied by the proportion of the grid cell covered by snow. The result can be regarded as the snow water equivalent of the portion of the snow cover within the grid cell, while the difference in other parts of the grid cell is regarded as caused by the estimation error between different datasets. When there is no snow in the grid cell, all differences are regarded as caused by the estimation error between different datasets. Through this filtering rule, the error caused by snow can be initially separated from the overall difference between datasets.

[0077] In step S3, all filtered difference datasets are fused using a Bayesian-based triangular hat (TCH) method. The TCH method is used to evaluate the relative uncertainty of gridded datasets without requiring any prior knowledge. This method considers the cross-correlation between datasets without requiring the data to be independent, but the number of datasets must be greater than or equal to three. The Bayesian-based triangular hat (BTCH) method is used to integrate different datasets without using any reference truth value, thereby improving its accuracy. The algorithm for the Bayesian-based triangular hat (BTCH) method is as follows:

[0078] S31. Calculate the probability density function for all difference datasets, where the i-th difference dataset NS i The probability density function p(NS) i NS t Specifically:

[0079]

[0080] In the formula, NS t ε represents the true value of the difference dataset. i Let σ be zero-mean white noise in the i-th difference dataset.i Let L(-) be the error covariance of the product of the i-th difference dataset, and L(-) be the likelihood function.

[0081] ε i =NS i -NS t

[0082] S32. Calculate the maximum probability maxL(NS) of the true value in the difference dataset. t NS i NS j The specific expression is:

[0083]

[0084] In the formula, NS j For the j-th difference dataset;

[0085] S33, In order to obtain the true value NS of the difference dataset t The maximum likelihood value is given by defining the cost function J, which is expressed as follows:

[0086]

[0087] S34. Set the first change in the cost function to zero to obtain the true value NS of the difference dataset. t Its specific expression is:

[0088]

[0089] S35, Define NS t =w i NS i +w j NS j So, let's define the weight w. i and w j The specific expression is:

[0090]

[0091] For multiple difference datasets, the true value NS of the difference datasets t =w1NS1+…+w N NS N The weight w is defined for each difference dataset. k This is obtained by minimizing a similar cost function;

[0092]

[0093] For each difference dataset, a weight w is defined. k Error covariance σ iIt is obtained by solving the relative uncertainty of each difference dataset through the Triangle Hat (TCH) method, which is specifically as follows:

[0094] S351, Due to the inability to obtain the true value NS of the difference dataset. t In order to solve for the true value NS of the difference dataset t Randomly select one dataset from N difference datasets and calculate the difference Y between it and the remaining difference datasets. i,M Its expression is as follows:

[0095] Y i,M =NS i -NS R =ε i -ε R

[0096] In the formula, NS R Y is a difference dataset randomly selected from N difference datasets. i,M Let S be a difference matrix with (N-1) rows and M columns, i = 1, 2, ..., N-1, where M is the total number of time samples in each NS product. The expression for its covariance matrix S is as follows:

[0097] S = cov(Y) i,M )

[0098] In the formula, cov(·) is the covariance operator, introducing an unknown N*N noise covariance matrix R, whose relationship with S is:

[0099]

[0100] In the formula, J N-1,N For the same matrix, its specific expression is:

[0101]

[0102] The definition of matrix R is:

[0103]

[0104] In the formula, σ ij =cov(ε i , ε j Let R be a matrix S, i = 1, ..., N, j = 1, ..., N. The covariance matrix S has N × (N+1) / 2 unknowns (the number of distinct elements in R), but only N × (N-1) / 2 equations (the number of distinct elements in S). Therefore, N "free" parameters must be reasonably determined to obtain a unique solution. To determine the N free parameters, the following objective function needs to be minimized according to the Kuhn-Tucker theorem.

[0105] S352, Galindo, and Palacio proposed the Kuhn-Tucker theorem, the objective function of which is:

[0106]

[0107] In the formula, The constraint functions are as follows:

[0108]

[0109] In the formula, Q is a diagonal matrix with elements σ on the diagonal. 11 ,…,σ NN These arrays can be calculated by iteratively minimizing the objective function of the Kuhn-Tucker theorem under initial conditions, thereby minimizing the objective function under constraints to obtain a set of free parameter solutions σ. 1N ,…,σ NN This refers to the relative uncertainty of each difference dataset.

[0110] While the difference caused by snow cover and the estimation errors between different datasets can be initially separated by using the snow cover data, estimation errors between different datasets still exist in the spatiotemporal context of snow cover. The triangular hat method is used to assess the uncertainty of each dataset, and posterior probabilities are used to weight the datasets according to the magnitude of the uncertainty, thereby achieving data fusion in the absence of actual measured data and eliminating estimation errors between different datasets in the spatiotemporal context of snow cover.

[0111] The beneficial effects of this invention are as follows: This invention provides a method for estimating the change in snow water equivalent based on inversion error. Addressing the characteristic that errors between products continuously evolve over time, this invention introduces snow cover data based on optical remote sensing as a basis for determining whether there is snow accumulation in a region and excluding the difference portion not caused by snow accumulation. Furthermore, it uses a Bayesian-based generalized triangular hat fusion method to fuse all datasets to extract the difference portion belonging to the change in snow water equivalent. Compared with existing technologies, this invention has the following advantages:

[0112] (1) Simplify the process of estimating snow water equivalent changes by estimating the difference between relatively mature and reliable precipitation datasets and soil water datasets without having to consider snow depth and snow density.

[0113] (2) Compared with the original estimation path, the spatial resolution is higher, reaching 0.1°*0.1°.

[0114] (3) Provide a completely independent numerical reference for the snow water equivalent variation of products produced through the original technical route.

[0115] In the description of this invention, it should be understood that the terms "center," "thickness," "upper," "lower," "horizontal," "top," "bottom," "inner," "outer," and "radial," etc., indicating orientation or positional relationships based on the orientation or positional relationships shown in the accompanying drawings, are only for the convenience of describing the invention and simplifying the description, and do not indicate or imply that the device or element referred to must have a specific orientation, or be constructed and operated in a specific orientation, and therefore should not be construed as a limitation of the invention. Furthermore, the terms "first," "second," and "third" are used for descriptive purposes only and should not be construed as indicating or implying the relative importance or the number of technical features implicitly specified. Therefore, a feature defined by "first," "second," and "third" may explicitly or implicitly include one or more of that feature.

Claims

1. A method for estimating snow water equivalent variation based on inversion error, characterized in that, Includes the following steps: S1. Obtain the precipitation dataset, perform interpolation on the precipitation dataset, and obtain several interpolation datasets; S2. Filter all difference datasets using snow cover data with the same spatiotemporal resolution to obtain the filtered difference dataset; S3. The filtered difference datasets are merged using the Bayesian triangular hat method to obtain the merged dataset. S4. Estimate the change in snow water equivalent by fusing datasets.

2. The method for estimating snow water equivalent variation based on inversion error according to claim 1, characterized in that, In S1, the method for processing the difference in the precipitation dataset is as follows: All precipitation datasets are interpolated to the same spatiotemporal resolution and equivalent units. The precipitation datasets include a first precipitation dataset based on direct precipitation measurements and a second precipitation dataset based on remote sensing soil moisture observations. The difference dataset is obtained by subtracting the values ​​of the second precipitation dataset from the values ​​of the first precipitation dataset.

3. The method for estimating snow water equivalent variation based on inversion error according to claim 2, characterized in that, In step S2, the method for filtering all difference datasets using snow cover data with the same spatiotemporal resolution is as follows: The spatiotemporal difference data is filtered based on the snow cover of the grid, and the specific expression is as follows: In the formula, S(x,y,t) represents the snow cover at time t and spatial location (x,y), with a value range of [0,100], D(x,y,t) represents the difference data at time t and spatial location (x,y), and D′(x,y,t) represents the filtered difference data at time t and spatial location (x,y).

4. The method for estimating snow water equivalent variation based on inversion error according to claim 3, characterized in that, In S3, the method for fusing all filtered difference datasets based on the Bayesian triangular hat method is as follows: S31. Calculate the probability density function for all difference datasets, where the i-th difference dataset NS i The probability density function p(NS) i |NS t Specifically: In the formula, NS t ε represents the true value of the difference dataset. i Let σ be zero-mean white noise in the i-th difference dataset. i Let L(-) be the error covariance of the product of the i-th difference dataset, and L(-) be the likelihood function. ε i =NS i -NS t S32. Calculate the maximum probability maxL(NS) of the true value in the difference dataset. t |NS i NS j The specific expression is: In the formula, NS j For the j-th difference dataset; S33, In order to obtain the true value NS of the difference dataset t The maximum likelihood value is given by defining the cost function J, which is expressed as follows: S34. Set the first change in the cost function to zero to obtain the true value NS of the difference dataset. t Its expression is as follows: S35, Define NS t =w i NS i +w j NS j Then define the weight w i and w j The specific expression is: For multiple difference datasets, the true value NS of the difference datasets t =w1NS1+…+w N NS N The weight w is defined for each difference dataset. k This is obtained by minimizing a similar cost function; For each difference dataset, a weight w is defined. k Error covariance σ i It is obtained by solving the relative uncertainty of each difference dataset through the triangular hat method.

5. The method for estimating snow water equivalent variation based on inversion error according to claim 4, characterized in that, In step S35, the method for solving the relative uncertainty of each difference dataset using the triangular hat method is as follows: S351. Randomly select one dataset from N difference datasets and calculate the difference Y between it and the remaining difference datasets. i,M Its specific expression is: Y i,M =NS i -NS R =ε i -ε R In the formula, NS R Y is a difference dataset randomly selected from N difference datasets. i,M Let S be a difference matrix with (N-1) rows and M columns, i = 1, 2, ..., N-1, where M is the total number of time samples in each NS product. The expression for its covariance matrix S is as follows: S=cov(Y i,M ) In the formula, cov(·) is the covariance operator, introducing an unknown N*N noise covariance matrix R, whose relationship with S is: In the formula, J N-1,N For the same matrix, its specific expression is: The definition of matrix R is: wherein, σ ij =cov(ε i , ε j ), i=1,…,N, j=1,…,N; S352. Establish the objective function of the Kuhn-Tucker theorem, minimize the objective function under constraints to obtain a set of free parameter solutions, and obtain the relative uncertainty of each difference dataset. The objective function of the Kuhn-Tucker theorem is as follows: In the formula, The constraint functions are as follows: In the formula, Q is a diagonal matrix.

Images