Satellite data-based sea ice thickness reanalysis data correction method
By constructing a spatiotemporally weighted linear regression model and combining satellite remote sensing data with reanalysis data, the problem of spatiotemporal continuity and physical consistency of sea ice thickness data in the Antarctic region was solved, achieving high-precision data correction and improving the accuracy and reliability of the data.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- TIANJIN UNIV
- Filing Date
- 2026-02-10
- Publication Date
- 2026-06-19
AI Technical Summary
In existing technologies, the fusion or correction methods for sea ice thickness data fail to fully consider the spatiotemporal correlation and physical differences between satellite remote sensing products and reanalysis products. This results in insufficient spatiotemporal continuity and physical consistency of correction results in complex ice-affected areas such as Antarctica, making it difficult to meet the requirements for high-resolution, large-scale data correction.
A spatiotemporal weighted linear regression model based on satellite data was constructed. By spatial matching and temporal correlation analysis between satellite remote sensing data and reanalysis data, the regression coefficient matrix was calculated to correct the reanalysis data. The data was then corrected by combining the climatological field of satellite remote sensing data.
It significantly improves the accuracy and reliability of sea ice thickness data, reduces systematic errors, provides higher-quality basic data products, and supports climate diagnostics and scientific research.
Smart Images
Figure CN122241052A_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of remote sensing geosciences technology, specifically relating to a method for correcting sea ice thickness reanalysis data based on satellite data. Background Technology
[0002] As a key component of the polar climate system, sea ice thickness variations directly affect energy and mass exchange at the sea-atmosphere interface, thus playing a significant role in regulating global and regional climate. Accurate sea ice thickness information is crucial for scientific research and operational applications such as climate diagnostics, numerical model validation, and polar environmental monitoring.
[0003] Currently, operational methods for obtaining sea ice thickness mainly fall into two categories: First, direct inversion products based on satellite remote sensing, such as the sea ice thickness data provided by the European Space Agency's Soil Moisture and Ocean Salinity Mission (SMOS). These products directly reflect observed signals and have high spatiotemporal resolution; however, in high-latitude regions, the data exhibits significant uncertainty and systematic bias due to factors such as sensor performance, atmospheric conditions, and inversion algorithms. Second, reanalysis products generated based on numerical models and data assimilation systems, such as the sea ice thickness reanalysis data provided by the Copernicus Marine Environment Monitoring Service (CMEMS). These products have good physical consistency and spatiotemporal continuity, but their accuracy is constrained by model physics, parameterization schemes, and the assimilated observational data (usually only including sea ice concentration, not directly including thickness observations), leading to systematic errors in reflecting the true sea ice thickness distribution.
[0004] It is worth noting that, due to differences in data sources and generation mechanisms, satellite remote sensing products and reanalysis products often exhibit significant differences in their characterization of sea ice thickness, especially in complex ice-affected regions such as Antarctica. Although verification and correction can be achieved through on-site observation data from ships and ice, the spatiotemporal coverage of measured data is extremely limited, making it difficult to directly support the need for large-scale, high-resolution data correction.
[0005] Existing technologies for sea ice thickness data fusion or correction primarily focus on the Arctic region and typically employ simple bias subtraction or spatial interpolation methods, failing to adequately consider the spatiotemporal correlations and physical differences between heterogeneous data sources. Traditional methods often neglect temporal correlations or rely solely on spatial proximity for interpolation, resulting in insufficient spatiotemporal continuity and physical consistency in the correction results. Therefore, there is an urgent need to develop an efficient correction method that can effectively integrate the advantages of multi-source observations, consider spatiotemporal correlations, and is applicable to large-scale gridded data. This will improve the accuracy and reliability of sea ice thickness reanalysis products and meet the demand for high-quality ice thickness data in climate research and applications. Summary of the Invention
[0006] The purpose of this invention is to provide a method for correcting sea ice thickness reanalysis data based on satellite data, aiming to solve the technical problems existing in the background art.
[0007] To achieve the above objectives, the present invention provides the following technical solution: This invention provides a method for correcting sea ice thickness reanalysis data based on satellite data, comprising the following steps: S1: Acquire and prepare daily average satellite remote sensing sea ice thickness products, daily average sea ice thickness reanalysis data products, and measured sea ice thickness data; S2: Preprocess the satellite remote sensing sea ice thickness product, spatially resample it to the same latitude and longitude grid as the reanalysis data product to achieve spatial matching; calculate the average value of the long-term series satellite remote sensing data and reanalysis data along the time dimension after spatial matching to obtain the corresponding remote sensing data climate field and reanalysis data climate field; spatially match the measured data with the reanalysis data grid. S3: For each grid point in the target area, based on the reanalysis data and satellite remote sensing data matched within the historical time period, a spatiotemporal weighted linear regression model is constructed to calculate the regression coefficient from the sea ice thickness anomaly value in the reanalysis data to the sea ice thickness anomaly value in the satellite remote sensing data, thereby obtaining a spatiotemporal regression coefficient matrix covering the target area. S4: For the time period to be corrected, the reanalysis data is corrected using the spatiotemporal regression coefficient matrix, and the corrected sea ice thickness field is obtained by combining the climatological field of satellite remote sensing data. S5: Using the measured data, evaluate the accuracy of the sea ice thickness data before and after correction.
[0008] Preferably, in step S2, the specific steps for spatially resampling the satellite remote sensing data to the latitude and longitude grid of the reanalysis data include: S21: Based on the polar projection parameters attached to the satellite remote sensing data, the projection coordinates of its grid points are inverted into latitude and longitude coordinates; S22: For each target point T(λ_T, φ_T) on the latitude and longitude grid of the reanalysis data, find the K nearest satellite remote sensing data points that have completed coordinate transformation within the preset search radius R; S23: Calculate the spherical distance d_k between the k-th remote sensing data point and the target point T, and perform a weighted average with the reciprocal of the p-th power of the distance as the weight to obtain the satellite remote sensing sea ice thickness interpolation H_T at the target point T, thus completing spatial matching; the interpolation formula is:
[0009]
[0010] Where H_k is the sea ice thickness value of the kth remote sensing data point.
[0011] Preferably, in step S2, the specific method for spatially matching the measured data with the reanalysis data grid is as follows: for each observation point in the measured data, calculate its spatial distance with all reanalysis data grid points, and determine the grid point with the smallest distance as the grid point that matches the measured point.
[0012] Preferably, in step S3, the specific steps for constructing the spatiotemporal weighted linear regression model include: S31: Combine the longitude, latitude, and time information of all reanalysis data grid points within the historical time period to construct a global spatiotemporal dataset; S32: For each target grid point P(λ_P, φ_P) within the target area, construct a query point using its spatial coordinates and the center time of the historical time period, and query its spatiotemporal neighbor points in the global spatiotemporal dataset; S33: A hierarchical random sampling strategy is adopted to select sample points for regression calculation from the spatiotemporal neighbor points found in the query; S34: For each sample point i, calculate the spatial distance difference Δx, Δy and the temporal distance difference Δt between it and the target point P, and substitute them into the Gaussian kernel function to calculate the spatiotemporal weight w of the sample point; S35: Extract the reanalysis sea ice thickness value and satellite remote sensing sea ice thickness value corresponding to each sample point i, and subtract the reanalysis climate state and satellite remote sensing climate state at the target point P respectively to obtain the reanalysis data anomaly sequence and the satellite remote sensing data anomaly sequence. S36: Using the outliers in the reanalysis data as the independent variable X, the outliers in the satellite remote sensing data as the dependent variable Y, and the spatiotemporal weight w as the weight, the coefficient a of the linear equation is solved by weighted least squares. This coefficient a is the regression coefficient a_ij of the target grid point P. The Gaussian kernel function mentioned in step S34 is:
[0013] Where Lx and Ly are preset spatial scale parameters, and Lt is a preset time scale parameter.
[0014] Preferably, the hierarchical random sampling strategy is as follows: the spatiotemporal neighboring points found are divided into N intervals according to their distance from the target point P from smallest to largest, and no more than M points are randomly selected from each distance interval to form the final sample point set for regression calculation.
[0015] Preferably, step S4 includes the following steps: S41: Obtain the reanalysis daily average data within the time period to be corrected, and calculate the reanalysis data climatology and remote sensing data climatology corresponding to that time period; read the pre-generated spatiotemporal regression coefficient matrix A; S42: For each day within the period to be corrected, calculate the anomalous field of the sea ice thickness data and its climate state for that day; S43: Multiply the anomaly field with the spatiotemporal regression coefficient matrix A grid-by-grid point to obtain the intermediate result field; S44: Add the intermediate result field to the remote sensing data climate field to obtain the corrected sea ice thickness field.
[0016] Preferably, step S44, "adding the intermediate result field to the remote sensing data climate field", specifically means: directly adding the value of each grid point in the intermediate result field to the climate reference value of the corresponding grid point in the remote sensing data climate field to obtain the corrected sea ice thickness value of that grid point.
[0017] Preferably, in step S5, the accuracy assessment is achieved by calculating the root mean square error, specifically by comparing the corrected sea ice thickness value matched to the same reanalysis data grid point with the measured value, calculating its root mean square error, and comparing it with the error of the original reanalysis data at that point.
[0018] Preferably, the satellite remote sensing sea ice thickness product is sourced from the European Space Agency's Soil Moisture and Ocean Salinity Mission (SMOS), and the sea ice thickness reanalysis data product is sourced from the Copernicus Marine Environment Monitoring Service (CMEMS).
[0019] Preferably, the method is specifically applied to the correction of sea ice thickness data in the Antarctic region.
[0020] The beneficial effects of this invention are as follows: By constructing a spatiotemporally weighted regression model, this invention effectively integrates the advantages of both satellite remote sensing observations and reanalysis data. On the one hand, this method overcomes the limitations of a single data source, utilizing the direct observation characteristics of satellite remote sensing data to correct model-based reanalysis products, achieving complementary advantages from heterogeneous data sources, significantly improving data accuracy, and reducing systematic errors from a single source. On the other hand, the model not only considers spatial proximity but also innovatively introduces temporal correlation, smoothly characterizing the decay of data bias with spatiotemporal distance through a Gaussian kernel function, which is more scientifically sound and physically grounded than traditional simple spatial interpolation or time averaging methods.
[0021] This invention fully utilizes the statistical relationships between heterogeneous data sources implicit in historical data, making the corrected product more statistically similar to actual observations. Validation with independent measured data shows that this method can effectively and systematically reduce the bias in sea ice thickness in reanalysis products, significantly improving the accuracy and reliability of the data. This provides higher-quality, more spatiotemporally consistent basic data products for scientific research such as climate diagnostics and forecast verification. Attached Figure Description
[0022] Figure 1 This is an overall flowchart of the method of the present invention; Figure 2 This is a spatial distribution map of sea ice thickness products obtained from satellite remote sensing and reanalysis in September 2012, as described in this embodiment of the invention. Figure 3 This is a spatial distribution map of sea ice thickness after the 2012 revision of this invention; Figure 4 This is a spatial distribution map of the measured sea ice thickness data in an embodiment of the present invention. Detailed Implementation
[0023] The technical solution of the present invention will now be described in detail with reference to the accompanying drawings and embodiments. It should be noted that the following embodiments are only used to explain the present invention and do not constitute a limitation on the scope of protection of the present invention.
[0024] This invention provides a method for correcting sea ice thickness reanalysis data based on satellite data, the method being as follows: Figure 1 As shown, it includes the following steps: S1: This embodiment takes the Antarctic region as an example. Download the daily sea ice thickness products released by the European Space Agency's Soil Moisture and Ocean Salinity Mission (SMOS) south of 60°S in September 2012 from publicly available data sources, the daily sea ice thickness reanalysis products provided by the Copernicus Marine Environment Monitoring Service (CMEMS) in the same region during the same period, and the Antarctic Sea Ice Thickness Field Observation Data Set released by the Antarctic Data and Archive Centre (AADC).
[0025] S2: Read an example SMOS file and obtain its projection information string. Use the Python toolkit pyproj.Proj to create a projection converter, perform backprojection calculations on the x, y coordinates of all SMOS grid points, and obtain the corresponding latitude and longitude arrays. Obtain the latitude and longitude grid of the CMEMS data.
[0026] For daily average SMOS data, the converted SMOS data (discrete point set) for the day is used to iterate through each point in the CMEMS grid, querying the four nearest SMOS points within a 25km radius. The distances between these four points and the CMEMS point are calculated, and an inverse distance-weighted (p=2) average is performed to obtain the SMOS interpolation thickness at that CMEMS point. This process is repeated for all CMEMS points for the day to form SMOS data aligned with the CMEMS grid. For example... Figure 2 The image shows a comparison of SMOS and CMEMS data, revealing morphological differences between the two products. The SMOS and CMEMS climatologies are obtained by averaging the generated SMOS latitude and longitude grid and the original CMEMS time series along the time dimension.
[0027] For AADC measured data, calculate the distance between each measured point and all CMEMS grid points, and find the nearest neighbor grid point as the verification point.
[0028] S3: This step aims to establish a spatiotemporal statistical correlation between EMS reanalysis data and SMOS remote sensing data. It calculates the specific regression coefficients for each spatial grid point in the target area, generating a coefficient matrix A that perfectly matches the target grid. Key parameters of the regression model are set, including spatial scale parameters Lx and Ly (both set to 0.25 degrees to control the attenuation range of spatial influence) and temporal scale parameter Lt (set to 7 days to control the attenuation range of temporal influence). Simultaneously, calculation strategy parameters are set: the initial number of spatiotemporal neighbor points M for each query (e.g., 100). These parameters collectively determine the size and structure of the sample set used for regression calculation at each grid point, aiming to balance computational efficiency and statistical robustness. For each grid point (i, j) within the target area, regression coefficient calculation is performed. The core process of this calculation includes: S31: Querying Spatiotemporally Neighboring Points. A three-dimensional query vector is constructed using the geographic coordinates (latitude and longitude) of the grid point and the selected representative historical time point (the average value of the time series). This vector is input into the global spatiotemporal search tree to quickly retrieve the M spatiotemporally nearest historical data points and obtain their distances. To optimize sample representativeness and avoid local overfitting, a hierarchical random sampling strategy is adopted. The M retrieved neighboring points are evenly divided according to their distance from the target point, from smallest to largest. From each distance interval, neighboring points are randomly selected to form the final sample points used for regression calculation. This strategy ensures that the samples are distributed across different spatiotemporal distances.
[0029] S32: Data Extraction and Weight Calculation. Extract the original sea ice thickness values from historical CMEMS and SMOS data based on the original time and spatial points corresponding to each sample point. Calculate the longitude difference Δx, latitude difference Δy, and time difference Δt between each sample point and the target grid point (i, j). Calculate the spatiotemporal weight w for each sample point using the Gaussian kernel function formula.
[0030] S33: Calculate outliers and perform weighted regression. Subtract the CMEMS and SMOS climatologies corresponding to the target grid point (i, j) from the original CMEMS and SMOS thickness values for each sample point to obtain two sets of outliers. Using the CMEMS outliers as independent variables and the SMOS outliers as dependent variables, weights are calculated using the previously calculated weights w, and the weighted least squares method is used to solve for the linear regression coefficients a_ij. This represents the average response of the SMOS anomaly at that specific grid point when the CMERS anomaly changes by one unit.
[0031] S4: This step utilizes the established correction model to correct the CMEMS reanalysis data for the target period. The specific implementation process is as follows: S41: Target Data and Climate State Preparation. Read the daily CMEMS reanalysis data and SMOS data for the target period to be revised (e.g., September 2012) and calculate the corresponding short-term climate state for that target period.
[0032] S42: Daily Correction Calculation. For each day within the target period, the daily anomaly field is calculated based on the daily CMEMS sea ice thickness field. This anomaly field reflects the deviation of the daily state from the current climatic mean. The corrected daily average sea ice thickness value is obtained by directly multiplying each coefficient in the regression coefficient matrix A by the anomaly value at the corresponding position in the anomaly field, and then adding the climatic baseline value at the corresponding position in the SMOS climatological field. Taking September 25, 2012 as an example, the result is shown in the figure below. Figure 3 As shown.
[0033] S5: Effect Evaluation. September 25, 2012, was selected as the date with the most AADC measured data points. The thickness values at these points were extracted from the original CMEMS data, SMOS data, and the data corrected according to this invention, and compared with the measured values. The results are shown in Table 1.
[0034] Table 1. Original data and corrected data results The results show that the method of the present invention corrects CMEMS data, and the corrected data significantly reduces errors compared to CMEMS, improving the consistency with the measured data.
[0035] The above embodiments only illustrate the basic principles and operating procedures of the present invention. The present invention is not limited to the above embodiments. Several improvements and modifications can be made to the present invention without departing from its principles, and all such improvements and modifications fall within the scope of the present invention. The scope of protection of the present invention is defined by the appended claims and their equivalents.
Claims
1. A method for correcting sea ice thickness reanalysis data based on satellite data, characterized in that, Includes the following steps: S1: Acquire and prepare daily average satellite remote sensing sea ice thickness products, daily average sea ice thickness reanalysis data products, and measured sea ice thickness data; S2: Preprocess the satellite remote sensing sea ice thickness product, spatially resample it to the same latitude and longitude grid as the reanalysis data product to achieve spatial matching; calculate the average value of the long-term series satellite remote sensing data and reanalysis data along the time dimension after spatial matching to obtain the corresponding remote sensing data climate field and reanalysis data climate field; spatially match the measured data with the reanalysis data grid. S3: For each grid point in the target area, based on the reanalysis data and satellite remote sensing data matched within the historical time period, a spatiotemporal weighted linear regression model is constructed to calculate the regression coefficient from the sea ice thickness anomaly value in the reanalysis data to the sea ice thickness anomaly value in the satellite remote sensing data, thereby obtaining a spatiotemporal regression coefficient matrix covering the target area. S4: For the time period to be corrected, the reanalysis data is corrected using the spatiotemporal regression coefficient matrix, and the corrected sea ice thickness field is obtained by combining the climatological field of satellite remote sensing data. S5: Using the measured data, evaluate the accuracy of the sea ice thickness data before and after correction.
2. The method according to claim 1, characterized in that, In step S2, the specific steps for spatially resampling satellite remote sensing data to the latitude and longitude grid of the reanalysis data include: S21: Based on the polar projection parameters attached to the satellite remote sensing data, the projection coordinates of its grid points are inverted into latitude and longitude coordinates; S22: For each target point T(λ_T, φ_T) on the latitude and longitude grid of the reanalysis data, find the K nearest satellite remote sensing data points that have completed coordinate transformation within the preset search radius R; S23: Calculate the spherical distance d_k between the k-th remote sensing data point and the target point T, and perform a weighted average with the reciprocal of the p-th power of the distance as the weight to obtain the satellite remote sensing sea ice thickness interpolation H_T at the target point T, thus completing spatial matching; the interpolation formula is: ； ； Where H_k is the sea ice thickness value of the kth remote sensing data point.
3. The method according to claim 1, characterized in that, In step S2, the specific method for spatially matching the measured data with the reanalysis data grid is as follows: for each observation point in the measured data, calculate its spatial distance with all reanalysis data grid points, and determine the grid point with the smallest distance as the grid point that matches the measured point.
4. The method according to claim 1, characterized in that, In step S3, the specific steps for constructing the spatiotemporal weighted linear regression model include: S31: Combine the longitude, latitude, and time information of all reanalysis data grid points within the historical time period to construct a global spatiotemporal dataset; S32: For each target grid point P(λ_P, φ_P) within the target area, construct a query point using its spatial coordinates and the center time of the historical time period, and query its spatiotemporal neighbor points in the global spatiotemporal dataset; S33: A hierarchical random sampling strategy is adopted to select sample points for regression calculation from the spatiotemporal neighbor points found in the query; S34: For each sample point i, calculate the spatial distance difference Δx, Δy and the temporal distance difference Δt between it and the target point P, and substitute them into the Gaussian kernel function to calculate the spatiotemporal weight w of the sample point; S35: Extract the reanalysis sea ice thickness value and satellite remote sensing sea ice thickness value corresponding to each sample point i, and subtract the reanalysis climate state and satellite remote sensing climate state at the target point P respectively to obtain the reanalysis data anomaly sequence and the satellite remote sensing data anomaly sequence. S36: Using the reanalysis data outliers as the independent variable X, the satellite remote sensing data outliers as the dependent variable Y, and the spatiotemporal weight w as the weight, the weighted least squares method is used to solve the linear equation. The coefficient a is the regression coefficient a_ij of the target grid point P; The Gaussian kernel function mentioned in step S34 is: ； Where Lx and Ly are preset spatial scale parameters, and Lt is a preset time scale parameter.
5. The method according to claim 4, characterized in that, The hierarchical random sampling strategy is as follows: the spatiotemporal neighboring points are divided into N intervals according to their distance from the target point P from smallest to largest, and no more than M points are randomly selected from each distance interval to form the final sample point set for regression calculation.
6. The method according to claim 1, characterized in that, The specific steps of step S4 include: S41: Obtain the reanalysis daily average data within the time period to be corrected, and calculate the reanalysis data climatology and remote sensing data climatology corresponding to that time period; read the pre-generated spatiotemporal regression coefficient matrix A; S42: For each day within the period to be corrected, calculate the anomalous field of the sea ice thickness data and its climate state for that day; S43: Multiply the anomaly field with the spatiotemporal regression coefficient matrix A grid-by-grid point to obtain the intermediate result field; S44: Add the intermediate result field to the remote sensing data climate field to obtain the corrected sea ice thickness field.
7. The method according to claim 6, characterized in that, Step S44, "adding the intermediate result field to the remote sensing data climate field", specifically means: directly adding the value of each grid point in the intermediate result field to the climate reference value of the corresponding grid point in the remote sensing data climate field to obtain the corrected sea ice thickness value of that grid point.
8. The method according to claim 1, characterized in that, In step S5, the accuracy assessment is achieved by calculating the root mean square error. Specifically, the corrected sea ice thickness value matched to the same reanalysis data grid point is compared with the measured value, the root mean square error is calculated, and the error of the original reanalysis data at that point is compared with the error of the original reanalysis data.
9. The method according to claim 1, characterized in that, The satellite remote sensing sea ice thickness product is from the European Space Agency's Soil Moisture and Ocean Salinity Mission (SMOS), and the sea ice thickness reanalysis data product is from the Copernicus Marine Environment Monitoring Service (CMEMS).
10. The method according to claim 1, characterized in that, The method is specifically applied to the correction of sea ice thickness data in the Antarctic region.