Large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data

Abstract

The invention relates to the field of forest resource monitoring, and discloses a large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data, and the method comprises the following steps: S1, collecting multi-source data; s2, carrying out data preprocessing and resampling; s3, multi-temporal data of Sentinel-2 are analyzed, a key phenological period is determined, normalized spectral features are developed, and vegetation is classified; s4, constructing a response canopy height of an explanatory variable and ICESat-2 ATL08, and training a forest canopy height regression model based on forest features based on a random forest regression algorithm; and S5, performing spatial inversion on the continuous interpretation variables by using the trained model, and generating a large-scale forest canopy height map. According to the method, discrete height data of ICESat-2, Sentinel-1 / 2 and auxiliary data are combined, and the forest canopy height is accurately inverted by adopting a random forest regression algorithm.

Description

Technical Field

[0001] The invention relates to the technical field of forest resource monitoring, and in particular to a large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data. Background Art

[0002] Forest canopy height (FCH) is the distance from the top of the canopy to the ground. It is a key parameter for accurately quantifying forest biomass and carbon storage, and plays a vital role in estimating forest-related carbon emissions, forest management, and decision-making. Although traditional field survey methods can provide very accurate FCH information over a large area, this method requires a lot of manpower and financial resources, is time-consuming, and is difficult to carry out quickly and smoothly in complex mountainous areas.

[0003] In recent years, airborne laser scanning (ALS) has been widely used as a technology for accurately measuring FCH, but the application of ALS technology in large-scale FCH mapping faces challenges due to the high cost of data collection on a large scale. The emergence of satellite-borne laser radar (LiDAR) technology has provided a possible solution to this problem. Its short revisit time and wide coverage make it possible to obtain regional or even global LiDAR data in a relatively short time.

[0004] Among them, ICESat-2 (Ice Cloud and Land Elevation Satellite-2) has received widespread attention since its launch in 2018 due to its onboard Advanced Topographic Laser Altimetry System (ATLAS). ICESat-2 uses photon counting LiDAR technology with a high repetition rate, which significantly enhances the detection capability of three-dimensional information. However, the ATL08 product of ICESat-2 can only provide sparse FCH data along the orbital direction, making it difficult to directly obtain spatially continuous forest canopy heights over a large area. Therefore, a method is needed to associate sparse LiDAR point clouds with spatially continuous variables to achieve accurate mapping of FCH over a large area.

[0005] To address the problem of spatial discontinuity of ICESat-2 data, researchers have tried to combine ICESat-2 data with multi-sensor data to draw more comprehensive forest canopy height maps. However, the application of this approach in large-scale canopy height modeling is still limited. A major problem in the existing technology is that multivariate quality control is insufficiently applied in canopy height estimation, inversion inaccuracies caused by data collection anomalies are rarely addressed, and the use of fixed response variables and a single inversion model may limit the accuracy of the estimation results due to the complexity of environmental variables and differences in forest conditions.

[0006] In addition, accurate vegetation distribution maps provide important information on the spatial distribution of vegetation types, canopy height distribution, and canopy structural characteristics, and can subdivide large areas of vegetation into smaller, more uniform areas. Among all remote sensing images, Sentinel-2 images are widely used in vegetation classification due to their high spatial resolution and rich spectral information. However, image processing of a single time series may mask important phenological events and affect classification accuracy. How to optimize temporal features and spectral data sets to achieve more accurate and robust vegetation classification remains a difficult point in current research. Summary of the invention

[0007] In view of the shortcomings of the prior art, the present invention provides a large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data, which solves the problems of spatial discontinuity and low precision caused by data sparsity and surface diversity when accurately mapping forest canopy height over a large area.

[0008] To achieve the above objectives, the present invention is implemented by the following technical solutions: A large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data, comprising the following steps:

[0009] Collect ICESat-2 ATL08 product data, Sentinel-1 and Sentinel-2 multi-temporal image data, auxiliary data, and sample plot field survey data;

[0010] Preprocess the collected data and resample all data to a uniform spatial resolution;

[0011] Analyze multi-temporal data from Sentinel-2 to identify key phenological periods and develop normalized spectral signatures to classify vegetation types in the study area;

[0012] Construct explanatory variables and the response canopy height of ICESat-2 ATL08. Based on the discrete response canopy height and the corresponding explanatory variables, train a forest canopy height regression model based on forest characteristics using the random forest regression algorithm.

[0013] The trained forest canopy height regression model is used to perform spatial inversion of continuous explanatory variables to generate a large-scale forest canopy height map.

[0014] Preferably, the data preprocessing step includes:

[0015] For the ATL08 data product of ICESat-2, the weak beam track samples were screened and seven indicators including the maximum canopy height, RH98, RH95, canopy height uncertainty, terrain uncertainty, sample point longitude, and sample point latitude were extracted for the strong beam track samples.

[0016] For Sentinel-1 images, two bands, VV and VH, are extracted;

[0017] For Sentinel-2 images, cloud removal was performed and visible bands 2 to 4, red edge and near-infrared bands 5 to 8, and short-wave infrared bands 11 to 12 were extracted;

[0018] For the sample plot survey data, the vegetation type of the sample plot and the accurate forest canopy height are mainly collected;

[0019] The collected multispectral data of Sentinel-2, SAR data of Sentinel-1, auxiliary data and sample survey data were uniformly resampled in spatial resolution.

[0020] Preferably, the step of determining the key phenological period includes:

[0021] Calculate vegetation indices for band data extracted by Sentinel-2, including normalized difference vegetation index, enhanced vegetation index and land water index;

[0022] Savitzky-Golay filtering algorithm is used to smooth the time series data to reduce the impact of noise;

[0023] The spectral coefficient of variation and change rate in the time series were analyzed, and three key phenological periods were determined in combination with the peak and valley characteristics of vegetation growth, namely the T1 period from the 90th to the 120th day, the T2 period from the 150th to the 180th day, and the T3 period from the 270th to the 300th day.

[0024] Preferably, the step of developing the normalized spectral features comprises:

[0025] Normalize the band data extracted from Sentinel-2 multi-temporal images to generate normalized spectral features;

[0026] Improved vegetation classification and spectral feature selection techniques were applied to screen out normalized spectral features with high spectral separation;

[0027] By analyzing the performance of normalized spectral features on different vegetation combinations, the spectral separation index and correlation coefficient of each normalized spectral feature were calculated;

[0028] According to the spectral separability index and correlation coefficient, the normalized spectral features that can maximize the separability of vegetation types are selected for the subsequent forest canopy height inversion.

[0029] Preferably, the selection of the normalized spectral feature is optimized based on the spectral separation index and correlation coefficient of the normalized spectral feature, and the specific steps include:

[0030] The spectral separation index is calculated for all possible combinations of normalized spectral features between each pair of vegetation types and is determined by the following formula:

[0031]

[0032] in, and Represents the average values ​​of different vegetation types a and b on feature m, σ a and σ b is the corresponding standard deviation;

[0033] The correlation between normalized spectral features was calculated using the coefficient of determination To measure the feature m i and m j The correlation between them, the coefficient of determination is determined by the following formula:

[0034]

[0035] Among them, cov n (m i ,m j ) represents feature m i and m j The covariance between var n (m i ) and m n (m j ) are their variances respectively;

[0036] Combining the spectral separability index and correlation coefficient, normalized spectral features with high separability and low correlation were selected to maximize the separability of vegetation types.

[0037] Preferably, the step of training a forest canopy height inversion regression model based on forest characteristics by a random forest regression algorithm comprises:

[0038] The scattering and optical characteristics of Sentinel-1 / 2 and auxiliary data were used to generate explanatory variables, and the response canopy height provided by ICESat-2 data was integrated to form the training and validation datasets of the model;

[0039] Apply the random forest regression algorithm to train the input data and optimize the model's predictive performance by setting model parameters, including the number and depth of trees;

[0040] The cross-validation method is used to adjust the hyperparameters of the random forest algorithm to ensure that the model has the best generalization ability on the validation dataset;

[0041] After completing the model training, an independent validation dataset was used to test the model’s prediction accuracy to ensure the accuracy and reliability of the forest canopy height inversion.

[0042] Preferably, the step of performing spatial inversion on the continuous explanatory variables using the trained forest canopy height regression model comprises:

[0043] Substitute the explanatory variables of different forest characteristic areas into the trained regression model based on forest characteristics, and obtain the forest canopy height maps of different forest characteristic areas respectively;

[0044] The forest canopy height maps of different forest characteristic areas are merged to obtain a complete and continuous forest canopy height map;

[0045] Preferably, the explanatory variables include:

[0046] Variables from the Sentinel-1 dataset, including:

[0047] VH polarization signal

[0048] VV polarization signal

[0049] VV1 signal, where VV1 = VV*VH

[0050] VV2 signal, where VV2 = VV / VH

[0051] VV3 signal, where VV3 = VH / VV

[0052] VV4 signal, where VV4 = (VV + VH) / 2

[0053] VV5 signal, where VV5 = square(VV*VH);

[0054] Variables from the Sentinel-2 dataset, including:

[0055] Original bands, including B2 blue light, B3 green light, B4 red light, B5 red edge 1, B6 red edge 2, B7 red edge 3, B8 near infrared, B11 short wave infrared 1, B12 short wave infrared 2

[0056] DVI, where DVI = B8-B4

[0057] EVI, the formula is: EVI = 2.5 * (B8-B4) / (B8 = 6 * B4-7.5 * B2 + 1)

[0058] RVI, where RVI = B8 / B4

[0059] SVAI, the formula is: SVAI = 1.5 * (B8-B4) / (B8 + B4 + 0.5)

[0060] SVAI1=(B8–B4) / (B8+B4+0.5)

[0061] LAI, formula: LAI = (2.5*(B8-B4 / B8+6*B4-7.5*B2+1))*3.618-0.118

[0062] ARVI, the formula is: ARVI = B8-2*(B4+B2") / (B8")+2*B4-B2)

[0063] RGVI, the formula is: RGVI = (B4-B3) / (B4+B3)

[0064] NDWI, the formula is: NDWI = (B8-B11) / (B8+B11)

[0065] NDVI, the formula is: NDVI = (B8-B4) / (B8+B4)

[0066] NDVI_max, where NDVI_max is the maximum NDVI value from June to August

[0067] NDVI_min, where NDVI_min is the minimum NDVI value from June to August

[0068] NDVI_stdDes, where NDVI_stdDes is the standard deviation of NDVI from June to August

[0069] NDVI1, formula: NDVI1 = (B8-B5) / (B8+B5)

[0070] NDVI2, formula: NDVI2 = (B8-B6) / (B8+B6)

[0071] NDVI3, formula: NDVI3 = (B8-B7) / (B8+B7)

[0072] Auxiliary datasets, including:

[0073] aspect, elevation, slope and VCF

[0074] Texture feature datasets, including:

[0075] Texture indices derived from NDVI, VV, and VH, which are calculated based on a 3×3 kernel, including the gray-level co-occurrence matrix texture index:

[0076] Contrast: NDVIcontrast, VVcontrast, VHcontrast

[0077] Variance: NDVIvar, VVvar, VHvar

[0078] Entropy: NDVIent, VVent, VHent

[0079] Correlation: NDVIcorr, VVcorr, VHcorr

[0080] Inertia: NDVIinertia, VVinertia, VHinertia

[0081] Difference variance: NDVIdvar, VVdvar, VHdvar.

[0082] Preferably, the forest canopy height regression model based on forest characteristics includes:

[0083] Based on the vegetation type map, resampling was used to generate forest type maps with 30m resolution and 30*100m resolution, and large forest areas were subdivided into smaller and more uniform areas with different forest characteristics. First, different response canopy heights were set for the forest canopy height regression model according to the forest height characteristics. The high-stem forest used ATL08 RH98 as the response canopy height, and the low-stem forest used ATL08 RH95 as the response canopy height, which was recorded as the preset response canopy height. Then, according to the forest canopy closure characteristics, different forest canopy height regression models were trained in different forest characteristic areas, which were recorded as forest canopy height regression models based on forest characteristics.

[0084] Preferably, the method further comprises performing multivariate quality control on the ICESat-2 data to remove low-quality samples, wherein the multivariate quality control comprises:

[0085] Sample point screening: remove non-forest type sample points and those sample points with maximum height deviation greater than a preset threshold from the ICESat-2 data. Specifically, remove sample points with ICESat-2 maximum canopy height index exceeding 30 meters or less than 3 meters.

[0086] Uncertainty assessment: Sample points with high uncertainty were excluded, including those with canopy height uncertainty or terrain uncertainty exceeding 30;

[0087] Geographic location correction: Using the information of ICESat-2’s geolocation error of 5 meters and the spatial resolution of the sample points of 11*100m, sample points that are inconsistent with the forest type map with a resolution of 30m and 30*100m are eliminated;

[0088] Statistical analysis methods: Cook's distance was applied to evaluate the impact of each sample point on the model prediction, which was calculated as follows: The Cook's distance of each sample measured the degree to which the prediction of the generalized linear model fitting the response canopy height to the explanatory variables would change after excluding that sample.

[0089] Preferably, the method further comprises verifying the generated forest canopy height map, specifically comprising:

[0090] Compare the resulting forest canopy height maps with ground truth data, including canopy height data acquired using aerial laser scanning technology;

[0091] Statistical comparison indicators, including root mean square error, absolute mean error and correlation coefficient, were used to evaluate the consistency and accuracy between the inverted forest canopy height map and the measured data;

[0092] Perform spatial consistency analysis to ensure that the canopy height map fully corresponds to the measured data in terms of geographical location;

[0093] Based on the validation results, necessary corrections and parameter optimizations are made to the forest canopy height map to improve the accuracy and reliability of future model predictions.

[0094] The present invention provides a large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data. It has the following beneficial effects:

[0095] 1. The present invention combines the discrete height data of ICESat-2 and the continuous multi-temporal image data of Sentinel-1 and Sentinel-2 as well as auxiliary data, and uses the random forest regression algorithm to train a forest canopy height regression model based on forest characteristics to accurately invert the forest canopy height. This method effectively improves the accuracy and spatial resolution of canopy height mapping, making the canopy height map more detailed, which is helpful for accurate monitoring and management of forest resources.

[0096] 2. The present invention significantly improves the classification accuracy of vegetation types by identifying key phenological periods and developing normalized spectral features. Through time series analysis, the growth stages and physiological states of different vegetation types are effectively distinguished, thereby improving the accuracy of classification results and providing more reliable data support for ecological research and land management.

[0097] 3. The present invention uses a multivariate quality control method to effectively identify and remove low-quality data points that affect the accuracy of model prediction during the data preprocessing stage. This process reduces the uncertainty in canopy height inversion and improves the overall stability and reliability of the model.

[0098] 4. The method of the present invention is not only applicable to specific forest types, but also can be adapted to different regions and different types of forests by adjusting different response variables and feature sets. This flexibility makes the model widely applicable to canopy height estimation in different forest regions around the world. BRIEF DESCRIPTION OF THE DRAWINGS

[0099] Figure 1 It is a schematic diagram of the method flow of the present invention;

[0100] Figure 2 It is a time series diagram of different vegetation types in the original band and vegetation index according to an embodiment of the present invention;

[0101] Figure 3 A coefficient of variation diagram between vegetation types in the original band and the vegetation index according to an embodiment of the present invention;

[0102] Figure 4 A graph showing the rate of change of different vegetation types in the original band and vegetation index according to an embodiment of the present invention;

[0103] Figure 5 Schematic diagram of separation (SI) of the best normalized spectral features (best NSF) in 28 vegetation combinations (VC) according to an embodiment of the present invention;

[0104] Figure 6 The vegetation type distribution map based on Scheme 7 of the embodiment of the present invention is shown in the illustrations, and the areas of different vegetation types in the classification results are shown in the illustrations. (ag) respectively represent the comparison between the vegetation type map drawn from high-resolution images of Google Earth or Sentinel-2 in different regions and the classification results of the four schemes. Among them, points, surfaces and columns of the same color are the same vegetation type;

[0105] Figure 7 The scatter plots are of the predicted forest canopy height and the response canopy height of the model constructed with ATL08 RH95, ATL08 RH98 and the preset response canopy height as the response canopy height in the independent validation set in the embodiment of the present invention. The (a) row and the (b) row represent the results generated before and after the multivariate quality control, respectively, where the black dotted line and the red solid line represent the 1:1 and regression lines, respectively, and the point density is represented by the color bar;

[0106] Figure 8 The scatter plots of the predicted forest canopy height and the response canopy height in the independent validation sets of different forest areas by constructing a single model and a forest canopy height regression model based on forest characteristics in all forest areas when the preset response canopy height is used as the response canopy height in the embodiment of the present invention. Rows (a) and (b) represent the results generated before and after multivariate quality control, respectively, where the black dashed line and the red solid line represent the 1:1 and regression lines, respectively, and the point density is represented by the color bar;

[0107] Fig. 9 A predicted FCH distribution map of a study area according to an embodiment of the present invention;

[0108] Fig.10 It is a scatter plot between ALS RH90 and predicted FCH according to an embodiment of the present invention. DETAILED DESCRIPTION

[0109] The following will be combined with the drawings of the present invention specification to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0110] Please see attached Figure 1 The embodiment of the present invention provides a large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data, comprising the following steps:

[0111] S1. Collect ICESat-2 ATL08 product data, Sentinel-1 and Sentinel-2 multi-temporal image data, auxiliary data (topography and vegetation coverage data), and sample survey data (field survey of vegetation types and UAV lidar data);

[0112] First of all, data collection is required, which is the basis and one of the key steps of the method of the present invention. Through comprehensive and accurate data collection, sufficient input data can be provided for the subsequent forest canopy height inversion to ensure the accuracy and reliability of the model. In the present invention, data collection covers a variety of data sources, including satellite lidar data, optical imaging data, and auxiliary terrain and climate data. The comprehensive use of these data helps to overcome the limitations of insufficient or inconsistent data in traditional methods, and lays a solid foundation for large-scale forest canopy height mapping.

[0113] In this embodiment, ATL08 product data is first obtained from the ICESat-2 satellite of the National Aeronautics and Space Administration (NASA). ICESat-2 is equipped with the Advanced Terrain Laser Altimetry System (ATLAS), which can accurately measure the three-dimensional structure of the earth's surface and vegetation. ATL08 product data contains forest canopy height information worldwide and is the core data source for the inversion of forest canopy height in the present invention.

[0114] In addition, it is necessary to collect multi-temporal image data from Sentinel-1 and Sentinel-2. The synthetic aperture radar (SAR) data provided by the Sentinel-1 satellite has all-weather and all-time observation capabilities, and can penetrate clouds to ensure that high-quality surface data can be obtained under all weather conditions. The Sentinel-2 satellite provides high-spatial resolution optical image data, which contains information from multiple bands and is suitable for analyzing the growth status and type distribution of vegetation. By acquiring Sentinel-2 image data of multiple phases in the study area, the dynamic changes in vegetation growth can be effectively captured, supporting subsequent key phenological period analysis and forest canopy height inversion.

[0115] In terms of auxiliary data, this embodiment also collects 30-meter resolution SRTM (Shuttle Radar Topography Mission) digital elevation model (DEM) data and Landsat Vegetation Continuum Facility (vegetation coverage) products, providing detailed terrain information and vegetation coverage information in the study area. These auxiliary data can help analyze the impact of environmental factors on forest canopy height.

[0116] In terms of sample survey data, this embodiment collects ground survey data of major vegetation types for training samples and verification samples for vegetation type mapping. In addition, this embodiment collects drone lidar data for verification of forest canopy height maps. By comprehensively utilizing these multi-source data, the present invention can overcome the limitations that may exist in a single data source, achieve accurate inversion of forest canopy height, and provide high-quality data support for large-scale forest management and ecological monitoring.

[0117] S2, preprocessing the collected data and uniformly resampling the spatial resolution of all data;

[0118] After completing the collection of multi-source data, data preprocessing is an important step to ensure data quality and consistency. This process is crucial for subsequent vegetation classification and forest canopy height inversion, because only high-quality and spatially consistent data can ensure the accuracy and reliability of the inversion results. Data preprocessing includes removing noise in the data (such as cloud effects) and unifying the spatial resolution of each data set to enable effective integration and analysis between different data sources.

[0119] In this embodiment, the collected Sentinel-2 image data is first de-clouded. Since clouds and cloud shadows can seriously affect the quality of optical image data, especially in multi-temporal analysis, cloud occlusion can cause data loss or distortion. Therefore, this embodiment uses the CF masking algorithm (Cloud-Free Masking Algorithm) in combination with a cloud probability image to process each Sentinel-2 image one by one, automatically identifying and removing cloud and cloud shadow areas. The CF masking algorithm can effectively detect and mask clouds, thereby retaining clear surface and vegetation information. This process ensures that the continuity and integrity of the data are maintained in images of all phases.

[0120] Next, we extract the visible light bands 2 to 4, the red edge and near-infrared bands 5 to 8, and the short-wave infrared bands 11 to 12 of the Sentinel-2 image, and the VV and VH bands of the Sentinel-1 image.

[0121] Furthermore, the collected multispectral data of Sentinel-2, SAR data of Sentinel-1, SRTM terrain data and vegetation coverage data are uniformly resampled for spatial resolution. The original spatial resolutions of different data sources may differ. For example, the multispectral band resolution of Sentinel-2 is 10 meters or 20 meters, the SAR data of Sentinel-1 is 10 meters, and the SRTM terrain data and vegetation coverage data are 30 meters. In order to ensure that these data can be seamlessly integrated in subsequent analysis, this embodiment resamples the data to a unified spatial resolution of 10 meters and 30 meters for vegetation classification and canopy height inversion, respectively. During the resampling process, the nearest neighbor interpolation method is used to retain the spatial characteristics and accuracy of the original data.

[0122] Then, for the ATL08 data product of ICESat-2, the weak beam track samples were screened and seven indicators of the strong beam track samples, including the maximum canopy height, RH98, RH95, canopy height uncertainty, terrain uncertainty, sample point longitude, and sample point latitude, were extracted;

[0123] Finally, for the UAV LiDAR data, denoising and classification techniques were applied to the point cloud data to remove noise and to classify the point cloud into ground points and canopy points. A digital elevation model (DEM) and digital surface model (DSM) with a resolution of 1 meter were generated from the classified point cloud data using linear interpolation, and the DEM was subtracted from the DSM to obtain the canopy height model (CHM). To verify the accuracy of the forest canopy height inversion results, the 30-meter percentage aggregation technique was applied to generate the ALSRH90.

[0124] S3, analyze the multi-temporal data of Sentinel-2, identify key phenological periods, develop normalized spectral signatures, and classify vegetation types in the study area;

[0125] Accurate classification of vegetation types is a key step in forest canopy height inversion. In order to improve classification accuracy, this paper makes full use of Sentinel-2's multi-phase data, analyzes the spectral characteristics of vegetation at different periods, determines key phenological periods, and develops normalized spectral features. Through in-depth analysis and optimization of these features, different types of vegetation can be effectively distinguished, thereby providing accurate input data for subsequent forest canopy height inversion. This process combines advanced remote sensing technology and spectral analysis methods, and is an important step in achieving high-precision forest canopy height maps.

[0126] In this example, the original bands (bands 2-8 and 11-12) with a resolution of 10 m were first extracted and resampled using Sentinel-2 images from April to October, and three vegetation indices, namely the normalized difference vegetation index (NDVI), the enhanced vegetation index (EVI) and the land water index (LSWI), were calculated.

[0127] On this basis, in order to generate a continuous time series graph, this embodiment adopts the 10-day average synthesis method to generate 21 synthetic images. These images are processed by the Savitzky-Golay (SG) algorithm, which is used to remove noise and smooth data to ensure the stability and reliability of the time series.

[0128] By analyzing these processed original bands and time series images of vegetation indices, this embodiment further identifies the key phenological periods of vegetation in the study area. Specifically, the characteristics of the original spectral curve and the peak-to-valley changes in the vegetation growth cycle were analyzed, and three key phenological periods were determined, namely T1 (90 to 120 days), T2 (150 to 180 days) and T3 (270 to 300 days). During these key phenological periods, the spectral coefficient of variation (CV) and spectral change rate (V) were calculated. These two indicators are used to quantify the changes in spectral characteristics of different vegetation types at different times.

[0129] Next, this embodiment develops a normalized spectral feature set (NSFS) by combining the original bands in pairs. These normalized spectral features are intended to reduce the spectral overlap between different vegetation types, thereby improving classification accuracy. To this end, this embodiment calculates the separation (SI) of the normalized spectral feature set on different vegetation combinations, the mean separation (SIM) of the original bands and vegetation index on different vegetation combinations, and combines the correlation (R) between the features of the normalized spectral feature set (NSFS) to obtain the best classification accuracy. 2), and finally the optimal normalized spectral feature (NSF) that can maximize the separability of vegetation types is screened out.

[0130] In the process of vegetation type mapping, this embodiment designs multiple classification schemes to evaluate the impact of key phenological periods and normalized spectral features on classification accuracy. Finally, this embodiment uses the random forest (RF) algorithm to classify vegetation types, combining the average value of the original band, the spectral time measurement (STM) of the vegetation index, and the spectral time measurement (STM) of the best normalized spectral features selected during the key phenological period to generate a vegetation type map in the study area.

[0131] Through the above steps, this embodiment achieves accurate classification of vegetation types in the study area, and provides high-quality input data for subsequent forest canopy height inversion.

[0132] S4, based on the discrete ICESat-2 ATL08 product data and the corresponding explanatory variables, the forest canopy height regression model based on forest characteristics was trained by the random forest regression algorithm;

[0133] On the basis of determining the vegetation type and key phenological period, constructing an accurate forest canopy height inversion model is one of the core steps of the present invention. To achieve this goal, the present invention adopts a random forest regression algorithm, combines normalized spectral features and vegetation index of key phenological periods, and constructs a canopy height inversion model that can reflect the characteristics of different forest types through the fusion of multi-source data. This model aims to effectively solve the shortcomings of traditional methods in dealing with complex forest structures, especially in applications in large-scale, highly heterogeneous areas, and significantly improves the accuracy and reliability of predictions.

[0134] In this embodiment, the best vegetation type distribution map obtained in step S3 is first used as the forest type map. In order to ensure that the forest type map is consistent with the explanatory variables in the random forest regression model, the grassland data is deleted in the embodiment, and the maximum resampling technique is used to unify all forest type maps to 30 meters resolution and 30*100m resolution. This is done to reduce the data error caused by inconsistent resolution and reduce the impact of ICESat-2 sample point positioning error and sample point resolution, thereby improving the stability and accuracy of the model.

[0135] Next, according to the different characteristics of forest types, the forests in the study area were divided into high-stem forests and low-stem forests, as well as high canopy closure forests and low canopy closure forests. High-stem forests include tree species with canopy heights generally greater than 8 meters, such as Quercus variabilis, Robinia pseudoacacia, and plantation forests; while low-stem forests include tree species with canopy heights generally less than 8 meters, such as Platycladus orientalis, Pinus tabulaeformis, Gleditsia sinensis, and Cotinus coggygria. High canopy closure forests include Quercus variabilis, Platycladus orientalis, Pinus tabulaeformis, Cotinus coggygria, and plantation forests; while low canopy closure forests include Robinia pseudoacacia and Gleditsia sinensis.

[0136] To further improve the accuracy of the model, this embodiment performs multivariate quality control on ICESat-2 data. First, sample points that do not belong to the forest type and sample points whose canopy height deviates seriously from the normal range (i.e., the canopy height is greater than 30 meters or less than 3 meters) are eliminated. Then, sample points with high uncertainty are excluded, including those sample points whose canopy height uncertainty or terrain uncertainty exceeds 30. Secondly, using the information that the geolocation error of ICESat-2 is about 5 meters and the spatial resolution of the sample point is 11*100m, sample points whose forest type map with a resolution of 30m is inconsistent with the forest type map with a resolution of 30*100m are eliminated. Finally, the Cook distance (D) method is used to evaluate and exclude sample points that produce large prediction errors in the model to reduce the uncertainty in the data. Cook distance is a technique widely used in regression analysis, which can effectively identify and exclude outliers that have a greater impact on the model prediction results. In this embodiment, Cook distance is used to quantify the influence of each sample point in the generalized linear model, and the threshold of the outlier is determined by the IQR rule to optimize the performance of the model. In this process, the ATL08 product data of ICESat-2 was used, and ATL08 RH98 was used as the canopy height of the tall-stem forest, while ATL08 RH95 was used as the canopy height of the low-stem forest, which were recorded as the preset response canopy heights.

[0137] Subsequently, this embodiment uses the random forest regression algorithm (RF) to construct forest canopy height regression models in low canopy closure forests and high canopy closure forests, respectively, which are recorded as regression models based on forest characteristics. The random forest algorithm is a powerful machine learning method that has the ability to process a large number of input variables and is particularly suitable for classification and regression tasks of remote sensing data. In the embodiment, the explanatory variables of the RF algorithm include 30-meter resolution Sentinel-1 and Sentinel-2 image data from June to August, SRTM terrain data, and vegetation coverage data. Specifically, the explanatory variables include:

[0138]

[0139]

[0140] In order to evaluate and optimize the performance of the model, this embodiment uses multiple quantitative indicators, including root mean square error (RMSE), relative root mean square error (rRMSE) and Pearson correlation coefficient (R). These indicators are used to measure the difference between the model prediction value and the reference value, the impact of data scale, and the correlation between the prediction value and the reference value. Through these indicators, the performance of the model under different forest types and canopy closure is evaluated, so as to further optimize the model parameters.

[0141] Finally, through these steps, this embodiment constructs a forest canopy height inversion model based on the random forest regression algorithm. This model can achieve high-precision prediction of forest canopy height under different forest types and complex terrain conditions, providing important technical support for the monitoring and management of large-scale forest resources.

[0142] S5. Use machine learning algorithms to spatially interpolate discrete ICESat-2 altitude data to generate large-scale forest canopy height maps.

[0143] After completing the construction of the forest canopy height inversion model, the next step is to generate a continuous large-scale forest canopy height map through spatial interpolation of discrete ICESat-2 height data. This step is crucial because although ICESat-2 data is highly accurate, its observation data is spatially discrete, and what is needed in practical applications is continuous forest canopy height information covering the entire study area. Spatial inversion of continuous explanatory variables through machine learning algorithms can effectively generate high-precision continuous canopy height maps, providing important data support for regional and even global forest monitoring and management.

[0144] In this embodiment, the explanatory variables of different canopy closure areas are first substituted into the trained regression model based on forest characteristics, and the forest canopy height map of high canopy closure and the forest canopy height map of low canopy closure are obtained respectively. Then, the two forest canopy height maps are merged to obtain a continuous forest canopy height map. Finally, the continuous canopy height map is compared with ALS RH90 to test the accuracy of the forest canopy height map.

[0145] Through the above steps, this embodiment successfully realizes the conversion from discrete observation data to continuous forest canopy height map, ensures the accurate prediction and spatial continuity of large-scale forest canopy height, and provides strong data support for large-scale forest resource management and ecological assessment.

[0146] The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data provided by the present invention effectively constructs a forest canopy height prediction model applicable to a wide range of regions by integrating ICESat-2 high-precision lidar data, multi-phase Sentinel image data, auxiliary data (topography and vegetation coverage data) and sample plot survey data (vegetation type field survey and drone lidar data). This method improves the accuracy of vegetation classification by using normalized spectral features and key phenological period information, thereby improving the accuracy of forest canopy height prediction. At the same time, a machine learning algorithm is introduced to construct a regression model for discrete ICESat-2 data, generating a continuous canopy height map, which greatly improves the application value and spatial resolution capability of canopy data. In addition, the implementation of this method contributes to more effective forest resource management, ecological monitoring and environmental protection, and is of great significance to ecological research and policy formulation.

[0147] Example:

[0148] This embodiment provides a large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data. The method combines satellite lidar data and multi-source remote sensing data and uses a random forest regression model to achieve high-precision large-scale forest canopy height mapping. The implementation steps of the present invention are described in detail below and further illustrated in conjunction with the accompanying drawings.

[0149] 1. Data collection and preprocessing

[0150] In this embodiment, this embodiment collects ICESat-2's ATL08 product data, Sentinel-1 and Sentinel-2 multi-temporal image data, auxiliary data (terrain and vegetation coverage data), and sample survey data (vegetation type field survey and drone lidar data). Specifically, the collected Sentinel-2 image data is first de-clouded. Since clouds and cloud shadows can seriously affect the quality of optical image data, especially in multi-temporal analysis, cloud occlusion can cause data loss or distortion. Therefore, this embodiment uses the CF masking algorithm (Cloud-Free Masking Algorithm) combined with the cloud probability image to process each Sentinel-2 image one by one, automatically identifying and removing cloud and cloud shadow areas. The CF masking algorithm can effectively detect and mask clouds, thereby retaining clear surface and vegetation information. This process ensures that the continuity and integrity of the data are maintained in images of all phases. Next, we extract the visible light bands 2 to 4, the red edge and near-infrared bands 5 to 8, and the short-wave infrared bands 11 to 12 of the Sentinel-2 image, and the VV and VH bands of the Sentinel-1 image.

[0151] In order to ensure data quality, this embodiment performs a unified resampling process on the spatial resolution of the collected Sentinel-2 multispectral data, Sentinel-1 SAR data, SRTM terrain data, and vegetation coverage data. The original spatial resolutions of different data sources may differ. For example, the multispectral band resolution of Sentinel-2 is 10 meters or 20 meters, the SAR data of Sentinel-1 is 10 meters, and the SRTM terrain data and vegetation coverage data are 30 meters. In order to ensure that these data can be seamlessly integrated in subsequent analysis, this embodiment resamples the data to a unified spatial resolution of 10 meters and 30 meters for vegetation classification and canopy height inversion, respectively. During the resampling process, the nearest neighbor interpolation method is used to retain the spatial characteristics and accuracy of the original data.

[0152] In addition, for the ATL08 data product of ICESat-2, the weak beam track samples were screened and the maximum canopy height, RH98, RH95, canopy height uncertainty, terrain uncertainty, sample point longitude, and sample point latitude of the strong beam track samples were extracted; finally, for the UAV lidar data, denoising and classification techniques were applied to the point cloud data to remove noise, and the point cloud was divided into ground points and canopy points. The digital elevation model (DEM) and digital surface model (DSM) with a resolution of 1 meter were generated from the classified point cloud data using linear interpolation, and the DEM was subtracted from the DSM to obtain the canopy height model (CHM). In order to verify the accuracy of the forest canopy height inversion results, the 30-meter percentage aggregation technique was applied to generate ALS RH90.

[0153] 2. Vegetation type mapping

[0154] In this example, the key phenological periods were determined by analyzing the multi-temporal data of Sentinel-2, and the normalized spectral signature was developed. The spectral-temporal signature of vegetation shows that the spectral responses of different vegetation types during the growth cycle are significantly different (e.g. Figure 2 In order to further improve the accuracy of vegetation type classification, this embodiment calculates the coefficient of variation of each spectral band and vegetation index (as shown in Figure 3 ) and the rate of change (as shown Figure 4 Based on this, three key phenological periods were identified: T1 (90-120 days), T2 (150-180 days) and T3 (270-300 days).

[0155] In each key phenological period, this example further constructed a normalized spectral feature set (NSFS), and optimized the feature set by calculating the separation index (SI) and correlation (R2), thereby selecting 18 optimal normalized spectral features for the final vegetation type mapping (e.g. Figure 5By comparing different classification schemes (as shown in Table 1), the scheme 7, which combines the average value of the original bands, the spectral time measurement (STM) of the vegetation index and the spectral time measurement (STM) of the best normalized spectral features screened out during the key phenological period, performed best in terms of classification accuracy, with an overall accuracy (OA) of 90.84% ​​and an average F1 score of 89.64.

[0156] Table 1 Accuracy of vegetation classification according to different schemes

[0157] Scheme OA (%) PA(%) UA(%) Avg.F1 Scheme1 80.10 81.10 76.37 77.79 Scheme2 74.72 78.02 72.90 72.95 Scheme3 79.84 80.93 75.23 77.05 Scheme4 88.15 88.77 83.38 84.86 Scheme5 88.96 90.20 86.12 87.20 Scheme6 89.45 90.80 86.96 88.38 Scheme7 90.84 91.97 88.23 89.64

[0158] The resulting vegetation type map (such as Figure 6 The figure shows the spatial distribution of the main vegetation types in the study area, including Platycladus orientalis, Pinus tabulaeformis, Quercus variabilis, Robinia pseudoacacia, Gleditsia sinensis, Cotinus coggygria, plantations and grasslands.

[0159] 3. Construction of canopy height inversion model

[0160] After the vegetation type mapping is completed, this embodiment constructs explanatory variables based on 30-meter resolution Sentinel-1 and Sentinel-2 image data from June to August, SRTM terrain data, and vegetation coverage data. In order to improve the accuracy of the model, this embodiment divides different forest types into high-stem forests and low-stem forests, as well as high canopy closure forests and low canopy closure forests, and uses ICESat-2's ATL08 RH98 and ATL08 RH95 as canopy height response variables for high-stem forests and low-stem forests, respectively, denoted as preset response canopy height. Finally, the forest canopy height regression model was trained using the random forest regression algorithm in high canopy closure forests and low canopy closure forests using the explanatory variables and the preset response canopy height, respectively, and denoted as a regression model based on forest characteristics.

[0161] In order to demonstrate the superiority of the present invention, this embodiment uses all ICESat-2 sample points with ATL08RH98 as the response canopy height, all ICESat-2 sample points with ATL08 RH95 as the response canopy height, and all ICESat-2 sample points use the preset response canopy height to train a regression model in all forest areas for comparison. The regression model based on forest characteristics invented in this embodiment shows significant advantages (such as Figure 7 , Figure 8 shown).

[0162] In addition, in order to further reduce the uncertainty in the data, this example introduces a multivariate quality control method and uses the Cook distance to screen the sample points, thereby eliminating outliers that have a significant impact on the model prediction results. In this process, D* quality control shows significant advantages in improving model performance (such as Figure 7 shown).

[0163] 4. Generation of canopy height map and accuracy verification

[0164] After the model training is completed, this embodiment uses the regression model based on forest characteristics to invert the continuous explanatory variables. First, the explanatory variables of different canopy closure areas are respectively substituted into the trained regression model based on forest characteristics to obtain the forest canopy height map of high canopy closure and the forest canopy height map of low canopy closure. Then, the two forest canopy height maps are merged to obtain a continuous forest canopy height map (such as Fig. 9 The inversion results show that the FCH values ​​in the study area range from 1.42 to 20.51 meters, with an average of 8.39 meters and a standard deviation of 3.28 meters.

[0165] In order to verify the accuracy of the interpolation results, this embodiment cross-compares the generated FCH map with the airborne laser scanning (ALS) data (e.g. Fig.10 As shown), the results show that there is a good correlation between the predicted FCH and ALSRH90, the R value is 0.73, and the root mean square error is 2.99 meters, which verifies the reliability of the method of this embodiment.

[0166] 5. Discussion of results

[0167] The results of this example show that by combining discrete observation data from ICESat-2 and multi-source remote sensing data from Sentinel-1 and Sentinel-2, a large-scale, high-precision forest canopy height map can be effectively generated using the random forest regression algorithm. This method shows good adaptability when dealing with different vegetation types and complex terrains, especially in both high-stem and low-stem forests in the study area. High prediction accuracy is achieved.

[0168] Although embodiments of the present invention have been shown and described, it will be appreciated by those skilled in the art that various changes, modifications, substitutions and variations may be made to the embodiments without departing from the principles and spirit of the present invention, and that the scope of the present invention is defined by the appended claims and their equivalents.

Claims

1. A large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data, characterized in that: The following steps are involved: Collect ICESat-2 ATL08 product data, Sentinel-1 and Sentinel-2 multi-temporal image data, auxiliary data, and sample survey data; Preprocess the collected data and resample all data to a uniform spatial resolution; Analyze multi-temporal data from Sentinel-2 to identify key phenological periods and develop normalized spectral signatures to classify vegetation types in the study area; Construct explanatory variables and the response canopy height of ICESat-2 ATL08. Based on the discrete response canopy height and the corresponding explanatory variables, train a forest canopy height regression model based on forest characteristics using the random forest regression algorithm. The trained forest canopy height regression model is used to perform spatial inversion of continuous explanatory variables to generate a large-scale forest canopy height map.

2. A large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 1, characterized in that: The data preprocessing steps include: For the ATL08 data product of ICESat-2, the weak beam track samples are removed and the seven indicators of the strong beam track samples, namely, maximum canopy height, RH98, RH95, canopy height uncertainty, terrain uncertainty, sample point longitude, and sample point latitude, are extracted; For Sentinel-1 images, two bands, VV and VH, are extracted; For Sentinel-2 images, cloud removal was performed and visible bands 2 to 4, red edge and near-infrared bands 5 to 8, and short-wave infrared bands 11 to 12 were extracted; For the sample plot survey data, the vegetation type of the sample plot and the accurate forest canopy height are mainly collected; The collected multispectral data of Sentinel-2, SAR data of Sentinel-1, auxiliary data and sample survey data were uniformly resampled in spatial resolution.

3. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 1, characterized in that: The steps for determining the key phenological period include: Calculate vegetation indices for band data extracted by Sentinel-2, including normalized difference vegetation index, enhanced vegetation index and land water index; Savitzky-Golay filtering algorithm is used to smooth the time series data to reduce the impact of noise; The spectral coefficient of variation and change rate in the time series were analyzed, and three key phenological periods were determined in combination with the peak and valley characteristics of vegetation growth, namely the T1 period from the 90th to the 120th day, the T2 period from the 150th to the 180th day, and the T3 period from the 270th to the 300th day.

4. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 1, characterized in that: The steps of developing the normalized spectral signature include: Normalize the band data extracted from Sentinel-2 multi-temporal images to generate normalized spectral features; Improved vegetation classification and spectral feature selection techniques were applied to screen out normalized spectral features with high spectral separation; By analyzing the performance of normalized spectral features on different vegetation combinations, the spectral separation index and correlation coefficient of each normalized spectral feature were calculated; According to the spectral separability index and correlation coefficient, the normalized spectral features that can maximize the separability of vegetation types are selected for subsequent vegetation type classification.

5. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 1, characterized in that: The selection of the normalized spectral features is optimized based on the spectral separation index and correlation coefficient of the normalized spectral features, and the specific steps include: The spectral separation index is calculated for all possible combinations of normalized spectral features between each pair of vegetation types and is determined by the following formula: in, and Represents the average values ​​of different vegetation types a and b on feature m, σ a and σ b is the corresponding standard deviation; The correlation between normalized spectral features was calculated using the coefficient of determination To measure the feature m i and m j The correlation between them, the coefficient of determination is determined by the following formula: Among them, cov n (m i ,m j ) represents feature m i and m j The covariance between n (m i ) and var n (m j ) are their variances respectively; Combining the spectral separability index and correlation coefficient, normalized spectral features with high separability and low correlation were selected to maximize the separability of vegetation types.

6. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 1, characterized in that: The step of training the forest canopy height inversion regression model based on forest characteristics by using the random forest regression algorithm comprises: The scattering and optical characteristics of Sentinel-1 / 2 and auxiliary data were used to generate explanatory variables, and the response canopy height provided by ICESat-2 data was integrated to form the training and validation datasets of the model; Apply the random forest regression algorithm to train the input data and optimize the model's predictive performance by setting model parameters, including the number and depth of trees; The cross-validation method is used to adjust the hyperparameters of the random forest algorithm to ensure that the model has the best generalization ability on the validation dataset; After completing the model training, an independent validation dataset was used to test the model’s prediction accuracy to ensure the accuracy and reliability of the forest canopy height inversion.

7. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 6, characterized in that: The step of using the trained forest canopy height regression model to perform spatial inversion on the continuous explanatory variables includes: Substitute the explanatory variables of different forest characteristic areas into the trained regression model based on forest characteristics, and obtain the forest canopy height maps of different forest characteristic areas respectively; The forest canopy height maps of different forest characteristic areas are merged to obtain a complete and continuous forest canopy height map; The explanatory variables include: Variables from the Sentinel-1 dataset, including: VH polarization signal VV polarization signal VV1 signal, where VV1 = VV*VH VV2 signal, where VV2 = VV / VH VV3 signal, where VV3 = VH / VV VV4 signal, where VV4 = (VV + VH) / 2 VV5 signal, where VV5 = square(VV*VH); Variables from the Sentinel-2 dataset, including: Original bands, including B2 blue light, B3 green light, B4 red light, B5 red edge 1, B6 red edge 2, B7 red edge 3, B8 near infrared, B11 short wave infrared 1, B12 short wave infrared 2 DVI, where DVI = B8-B4 EVI, the formula is: EVI = 2.5 * (B8-B4) / (B8 = 6 * B4-7.5 * B2 + 1) RVI, where RVI = B8 / B4 SVAI, the formula is: SVAI = 1.5 * (B8-B4) / (B8 + B4 + 0.5) SVAI1=(B8–B4) / (B8+B4+0.5) LAI, formula: LAI = (2.5*(B8-B4 / B8+6*B4-7.5*B2+1))*3.618-0.118 ARVI, the formula is: ARVI = B8-2*(B4+B2") / (B8")+2*B4-B2) RGVI, the formula is: RGVI = (B4-B3) / (B4+B3) NDWI, the formula is: NDWI = (B8-B11) / (B8+B11) NDVI, the formula is: NDVI = (B8-B4) / (B8+B4) NDVI_max, where NDVI_max is the maximum NDVI value from June to August NDVI_min, where NDVI_min is the minimum NDVI value from June to August NDVI_stdDes, where NDVI_stdDes is the standard deviation of NDVI from June to August NDVI1, formula: NDVI1 = (B8-B5) / (B8+B5) NDVI2, formula: NDVI2 = (B8-B6) / (B8+B6) NDVI3, formula: NDVI3 = (B8-B7) / (B8+B7) Auxiliary datasets, including: aspect, elevation, slope and VCF Texture feature datasets, including: Texture indices derived from NDVI, VV, and VH, which are calculated based on a 3×3 kernel, including the gray-level co-occurrence matrix texture index: Contrast: NDVIcontrast, VVcontrast, VHcontrast Variance: NDVIvar, VVvar, VHvar Entropy: NDVIent, VVent, VHent Correlation: NDVIcorr, VVcorr, VHcorr Inertia: NDVIinertia, VVinertia, VHinertia Difference variance: NDVIdvar, VVdvar, VHdvar.

8. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 6, characterized in that: The forest canopy height regression model based on forest characteristics includes: Based on the vegetation type map, resampling was performed to generate forest type maps with a resolution of 30m and 30*100m, subdividing large forest areas into smaller and more uniform areas with different forest characteristics; Firstly, different response canopy heights were set for the forest canopy height regression model according to the forest height characteristics. The high-stem forest used ATL08 RH98 as the response canopy height, and the low-stem forest used ATL08 RH95 as the response canopy height, which were recorded as the preset response canopy heights. Then, according to the characteristics of forest canopy closure, different forest canopy height regression models are trained in different forest characteristic areas, which are recorded as forest canopy height regression models based on forest characteristics.

9. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 1, characterized in that: The method further includes performing multivariate quality control on the ICESat-2 data to remove low-quality samples, wherein the multivariate quality control includes: Sample point screening: remove non-forest type sample points and those sample points with maximum height deviation greater than a preset threshold from the ICESat-2 data. Specifically, remove sample points with ICESat-2 maximum canopy height index exceeding 30 meters or less than 3 meters. Uncertainty assessment: Sample points with high uncertainty were excluded, including those with canopy height uncertainty or terrain uncertainty exceeding 30; Geographic location correction: Using the information of ICESat-2’s geolocation error of 5 meters and the spatial resolution of the sample points of 11*100m, sample points that are inconsistent with the forest type map with a resolution of 30m and 30*100m are eliminated; Statistical analysis methods: Cook's distance was applied to evaluate the impact of each sample point on the model prediction, which was calculated as follows: The Cook's distance for each sample measures the degree to which the prediction of the generalized linear model fitting the response canopy height to the explanatory variables changes after excluding that sample.

10. The large-scale forest canopy height inversion method based on ICESat-2 and multi-source remote sensing data according to claim 1, characterized in that: The method also includes verifying the generated forest canopy height map, specifically comprising: Compare the resulting forest canopy height maps with ground truth data, including canopy height data acquired using aerial laser scanning technology; Statistical comparison indicators, including root mean square error, absolute mean error and correlation coefficient, were used to evaluate the consistency and accuracy between the inverted forest canopy height map and the measured data; Perform spatial consistency analysis to ensure that the canopy height map fully corresponds to the measured data in terms of geographical location; Based on the validation results, necessary corrections and parameter optimizations are made to the forest canopy height map to improve the accuracy and reliability of future model predictions.

Images