A method for improving the precision of human footprint index

Abstract

This invention relates to the field of ecological protection zone identification technology, and particularly to a method for improving the accuracy of the Human Footprint Index. The method comprises the following steps: S1, Population Grid Redistribution: A construction land mask is extracted from land cover data. Population statistics data from the construction land mask and population grid data are input into the iURBAN algorithm to redistribute the population grid data. Population statistics data are then used to verify the accuracy of the population grid data before and after redistribution. S2, Fragmentation Index Scoring: A movement speed grid is calculated. City tiers are classified based on land cover data. Based on the city tier classification and movement speed grid, different levels of urban patches representing urban edges are obtained. Landscape pattern indicators representing fragmentation are calculated, and fragmentation index scores are assigned based on the landscape pattern indicator calculation results. S3, Human Footprint Index Calculation. This invention redistributes population data from the population grid to optimize the population distribution and improve the accuracy of the population grid data.

Description

Technical Field

[0001] This invention relates to the field of ecological protection zone identification technology, and in particular to a method for improving the accuracy of the human footprint index. Background Technology

[0002] Human activity intensity has become an important reference for scientific research and policy-making, providing guidance for spatial planning, land use management, infectious disease transmission, and habitat protection. Calculating human activity intensity is crucial for public administration and maintaining biodiversity. Human activities not only exacerbate ecosystem stress but also lead to problems such as forest degradation, habitat loss, and species extinction, thus reducing biodiversity. We have now entered a new geological era dominated by humans, where human activities have become a key factor affecting the survival of all species on Earth. Stress-based assessments focus on measuring the stress that various human activities bring to ecosystems, including indices such as human footprint, human alteration index, and ecological footprint. Among studies related to human activity intensity, the human footprint index is currently the most frequently used indicator.

[0003] Currently, research on human footprints focuses more on ecologically important areas, aiming to measure the intensity of human activity in these areas through human footprint indices, thereby providing a reference for ecological protection. However, population grid data commonly used to calculate human footprints, including Worldpop and LandScan, often suffer from inaccurate population distribution in underdeveloped areas such as rural areas. These data can also represent large populations in uninhabited areas of the real world, such as primary forests. This affects the accuracy of subsequent human footprint calculations. Existing human footprint indices do not adequately consider human stress, failing to account for human-induced fragmentation of farmland and forest land. Therefore, they often underestimate the intensity of human activity in urban fringe areas, affecting the accuracy of human activity intensity mapping, and urgently need improvement. Summary of the Invention

[0004] The features and advantages of the present invention are set forth in part in the description which follows, or may be apparent from the description, or may be learned by practicing the invention.

[0005] To overcome the problems of existing technologies, this invention provides a method for improving the accuracy of the human footprint index, specifically including the following steps: S1. Population grid redistribution: The construction land mask is extracted from the land cover data. The population statistics data of the construction land mask and the population grid data are input into the iURBAN algorithm to realize the redistribution of population grid data. The accuracy of the population grid data before and after redistribution is verified by the population statistics data. S2. Fragmentation Index Scoring: Navigable waterways are calculated using river flow data, and movement speed grids are calculated using road data, railway data, navigable waterway data, and land cover data; impervious patches are obtained based on land cover data, and these patches are classified into city tiers based on the redistributed population; urban fringe areas of different city tiers are obtained based on city tier classification and movement speed grids, and landscape pattern indicators representing fragmentation are calculated within the urban fringe areas of the urban patches; fragmentation index scoring is then assigned based on the landscape pattern indicator calculation results. S3. Human Footprint Index Calculation: Select land cover data, redistributed population grid data, railway data, highway data, nighttime light data, and navigable waterway data to calculate land use indicators, population density indicators, road accessibility indicators, power infrastructure indicators, and navigable waterway indicators. Add these scores to obtain the original human footprint score. Then, add the fragmentation indicator score to the original human footprint score to obtain the human footprint index.

[0006] Preferably, the population grid data redistribution mainly includes the following steps: S101. Based on land cover data and population grid data, obtain the area size and population of each built-up area patch, and determine the importance of the built-up area patch through the two. After comparing with the urban population in the statistical yearbook, the built-up area patches are divided into urban built-up area patches and rural built-up area patches. S102. Based on the size of the rural built-up area patches and the population of the rural built-up area patches, allocate the rural population data in the statistical yearbook to each rural built-up area patch. S103. Based on the proportion of built-up area and population of the pixels within each rural built-up area patch, the statistical yearbook population allocated to each rural built-up area patch is assigned to the pixels within the patch.

[0007] Preferably, the root mean square error is used to evaluate the accuracy of the population grid data before and after redistribution using demographic data.

[0008] Preferably, step S2 further includes quantifying the city's reachability within one hour, specifically including the following steps: S20101, City Classification: Statistically determine the population within impermeable patches and classify impermeable patches into large cities, medium cities, and small cities based on the population size. S20102. Create a movement speed grid: Create corresponding movement speed grids based on roads, railways, navigable waterways, land cover, altitude, and slope. The resolution of each movement speed grid is 1000 m, and the cell value represents the movement speed in each cell. S20103. The cumulative cost algorithm is used to obtain the reachable range of each city patch within 1 hour.

[0009] Preferably, the scoring of the fragmentation index further includes the following steps: S20201. Generate 8-directional buffer zones: Generate a 1 km concentric circular buffer zone within the patch city edge area of ​​each city, and dequantize the fragmentation of forests and farmland within the buffer zone; Since the gradient of the fragmentation of forests and farmland varies in different directions, the concentric circular buffer zone is divided into 8 sectors, each sector expanding radially outward at 45°. S20202, Calculation of landscape indicators for farmland-forest fragmentation; the landscape indicators include patch density, maximum patch index, average patch area, aggregation index, patch connectivity index, and average fractal dimension. S20203. Fragmentation score: The degree of fragmentation is scored based on the landscape index calculation results of cultivated land and forest in each concentric buffer zone.

[0010] Preferably, the formula for the patch density index is as follows: , In the formula, For patch density, This represents the number of patches of a specific land cover class in the landscape. The total area of ​​the landscape; The formula for the maximum plaque index is as follows: , In the formula, The maximum plaque index, This indicates the area of ​​the largest patch in the landscape; this indicator is used to determine the dominant patch type in the landscape. The formula for the average patch area is as follows: , In the formula, The average patch area, This represents the area of ​​a specific land cover category within the landscape. The number of patches in a landscape for a specific land cover category; The formula for the aggregation index is as follows: , In the formula, The aggregation index, This represents the number of adjacent patches of the same land cover category for a given patch. This represents the maximum number of adjacent patches of the same land cover category within a given land cover category. The formula for the patch connectivity index is as follows: , In the formula, This is the plaque connectivity index. The perimeter of the patch, expressed in pixels. The area of ​​the patch is expressed in pixels. This represents the total number of pixels in the landscape. The formula for the average fractal dimension is as follows: , In the formula, The average fractal dimension, The first in a certain land cover category The perimeter of each patch, The first in a certain land cover category The area of ​​each patch, This represents the total number of patches within a given land cover category.

[0011] Preferably, when the urban edge areas of different urban patches overlap, the fragmentation score of the overlapping parts is processed according to the following rules: (1) If the urban fringe areas of high-level cities and low-level cities overlap, the fragmentation score of the overlapping area shall be the fragmentation score of the urban fringe area of ​​the high-level city. (2) If the urban fringe areas of cities of the same level overlap, the fragmentation score of the overlapping area shall be the maximum value of the fragmentation score of that area.

[0012] The beneficial effects of this invention are: This invention measures the development of construction land using land use intensity, calculates the potential human activity pressure on the surrounding environment from construction land with different land use intensities, and assigns a pressure score to this potential human activity pressure. Simultaneously, it redistributes population data from the population grid to optimize population distribution and improve the accuracy of the population grid data. Finally, it uses indicators such as power infrastructure, land use, road accessibility, population density, and navigable waterways, along with the pressure score assigned to construction land on the surrounding environment from potential human activity, to calculate the human footprint in each region. The aim is to obtain high-precision human activity intensity data, fill the gaps in high-precision human footprint maps, and identify human activity hotspots and colds, providing data support for protected area planning, biodiversity conservation, and other related work. Attached Figure Description

[0013] The present invention will be described in detail below with reference to the accompanying drawings and examples. The advantages and implementation methods of the present invention will become more apparent from this description. The accompanying drawings are for illustrative purposes only and do not constitute any limitation on the present invention. In the accompanying drawings: Figure 1This is a flowchart illustrating a method for improving the accuracy of the human footprint index in a specific embodiment of the present invention. Figure 2 This is a schematic diagram of an 8-directional buffer zone for a method to improve the accuracy of the human footprint index in a specific embodiment of the present invention. Figure 3 This is a map showing the improved human footprint index mapping in the Indochina Peninsula from 2015 to 2022, as part of a method for improving the accuracy of the human footprint index in a specific embodiment of the present invention. Figure 4 This is a comparison chart of the distribution of human footprint values ​​in the Indochina Peninsula before and after the improvement in a method for improving the accuracy of the human footprint index in a specific embodiment of the present invention; Figure 5 for Figure 4 Enlarged images of points A, B, and C after the improvements were made, along with high-resolution Google Earth images. Detailed Implementation

[0014] The embodiments of the present invention will be described in further detail below with reference to the accompanying drawings and examples. The following examples are for illustrative purposes only and should not be construed as limiting the scope of the invention.

[0015] like Figure 1 As shown, the present invention provides a method for improving the accuracy of the human footprint index, specifically including the following steps: The datasets used in this embodiment include: NPP / VIIRS (Suomi National Polar-orbiting Partnership / Visible Infrared Imaging Radiometer) data, GLC_FCS30D data, OSM (OpenStreet Map) road data, LandScan data, and demographic data. All the above raster data were resampled to a uniform resolution of 30″, with data acquired between 2015 and 2022. Among them, the GLC_FCS30D data is the first 30-m resolution global fine-grained land cover dynamic product using continuous change detection technology. The system includes data on 10 major land cover types: cultivated land, forest, shrubland, grassland, tundra, wetland, impervious surfaces, bare land, water bodies, and permanent snow and ice. The Visible Infrared Imaging Radiometer (VIIRS) collects radiometric images of land, atmosphere, ice, and ocean in the visible and infrared bands. It is an extension and improvement of the High Spatial Resolution Radiometer (VHRR) and the Medium Resolution Imaging Spectroradiometer (MODIS) series of Earth observation instruments. Nighttime light remote sensing images acquired by the Visible Infrared Imaging Radiometer aboard the new-generation Suomi NPP satellite accurately record the intensity of nighttime light radiation. Annual composite version 2.1 NPP / VIIRS images from 2015–2021 and annual composite version 2.2 from 2022 are used. The NPP / VIIRS imagery comprises the NPP / VIIRS imagery dataset from 2015 to 2022. LandScan data provides more realistic population distribution data. Demographic statistics are derived from national or governmental census data, statistics published by the United Nations Department of Economic and Social Affairs, and publicly available data from the World Bank. HydroSHEDS acc provides centralized river flow data within the dataset.

[0016] S1. Population grid redistribution: The construction land mask is extracted from the land cover data. The population statistics data of the construction land mask and the population grid data are input into the iURBAN algorithm to realize the redistribution of population grid data. The accuracy of the population grid data before and after redistribution is verified by the population statistics data. Furthermore, in this embodiment, the iURBAN model is a model for redistributing the rural population based on the size of built-up area patches and the rural population figures in statistical yearbooks. Since rural populations are often distributed in uninhabited areas such as forests, it is necessary to redistribute the rural population. The redistribution of rural population grid data based on the iURBAN model mainly includes the following steps: S101. Based on land cover data and population grid data, the area size and population of each built-up area patch are obtained, and the importance of the built-up area patch is determined by the two. After comparing with the urban population in the statistical yearbook, the built-up area patches are divided into urban built-up area patches and rural built-up area patches. S102, Based on the size of the rural built-up area patch and the population of the rural built-up area patch, the rural population in the statistical yearbook is allocated to each rural built-up area patch; S103, based on the proportion of built-up area and population of the pixels inside the rural built-up area patch, the statistical yearbook population allocated to each rural built-up area patch is assigned to the pixels inside the patch.

[0017] Furthermore, in this embodiment, the GLC_FCS30D data is first resampled to a resolution of 1000 m. If the proportion of impermeable surface area in the 1000 m grid is greater than 50%, the reclassified result is an impermeable surface. An impermeable surface mask is generated based on the resampled result. Then, based on the impermeable surface mask, the population count in the LandScan population grid is restricted to within the impermeable surface according to Formula 1, as shown below: Formula 1, In the formula, To limit the population to the built-up area after the first Population of each built-up area patch For the first The original population of each built-up area patch This represents the total population of the grid within the study area. This represents the total original population of all built-up area patches. The iURBAN model is then used to redistribute the LandScan population grid.

[0018] Furthermore, the steps for verifying the accuracy of population grid data before and after population redistribution are as follows: The total population of each country / region in the LandScan population grid data before and after redistribution is compared with the total population in official statistics. The accuracy is evaluated using the root mean square error (RMSE). Official statistics include, but are not limited to, population census data from various levels of government, as well as UN population statistics. The specific formula is shown in Formula 2: Formula 2, In the formula, For the first Total population grid of a country / region For the first The total population of a country / region as recorded in the United Nations population statistics. Number of countries / regions.

[0019] S2. Fragmentation Index Scoring: Navigable waterways are calculated using river flow data; movement velocity grids are calculated using road, railway, navigable waterway, and land cover data; impervious patches are obtained based on land cover data, and these patches are classified into city tiers according to the redistributed population; urban fringe areas of different city tiers are obtained based on city tier classification and movement velocity grids; landscape pattern indicators representing fragmentation are calculated within the urban fringe areas of these city patches; and fragmentation index scoring is assigned based on the landscape pattern indicator calculation results. The specific steps include: S201, Quantification of the city's reachable range within 1 hour

[0020] The urban fringe areas of cities of different tiers are determined by calculating the reachable range within one hour for cities of different tiers. This involves the following steps: S20101. City Classification: Population within impermeable patches is statistically analyzed, and these patches are classified into large cities, medium-sized cities, and small cities based on population size. The specific classification is shown in Table 1. Table 1. City Ranking Based on Population , S20102. Create a movement speed grid: Create corresponding movement speed grids based on roads, railways, navigable waterways, land cover, altitude, and slope. The resolution of each movement speed grid is 1000 m, and the cell value represents the movement speed in each cell.

[0021] The road speed raster is derived from the `maxspeed` attribute in OSM road vector data. This attribute records the maximum speed of a road vector. Some road vectors have empty `maxspeed` values; these are filled with the average `maxspeed` value of roads of the same level in the corresponding country, based on the road category and country of the road vector. Road vectors from 2015–2022 are converted to raster data and assigned the `maxspeed` attribute. If multiple road vectors pass through a cell, the cell value is the `maxspeed` attribute of the highest-level road.

[0022] The railway speed raster is obtained based on OSM railway vector data. Unlike road vector data, railway vector data is not classified. A uniform speed of 24.3 km / h is assigned to the railway vector data. The railway vector data from 2015 to 2022 is converted into a raster to obtain the railway speed raster.

[0023] The navigable waterway velocity grid is derived from navigable waterway data calculated based on HydroSHEDS river flow data from 2015–2022. River depth is calculated using the HydroSHEDS river flow data according to the following formula 3-6: Formula 3, Formula 4, Formula 5, Formula 6, In the formula, For river flow, For the width of the river, For river flow velocity, The cross-sectional area of ​​the river. The criteria for determining navigable waterways are: (1) the river depth is greater than 2 m; (2) the distance from the luminous pixels in the nighttime light data is within 4 km.

[0024] The land cover movement speed raster was obtained based on GLC_FCS30D land cover data, with cell values ​​representing the speed of walking through different map cover types. Referring to existing methods for assigning movement speed values ​​to different land cover types, the 30 land cover categories in the GLC_FCS30D data of the study area were reclassified and then their movement speed values ​​were assigned, as shown in Table 2. Table 2. Assignment of Land Cover Movement Speed , The movement speed grid is obtained by overlaying the railway movement speed grid, the navigable waterway movement speed grid, and the land cover movement speed grid. If the land cover movement speed grid and other movement speed grids have values ​​greater than 0 in a certain cell, then the value of that cell is the value of the other movement speed grid.

[0025] S20103. The cumulative cost algorithm is used to obtain the 1-hour reachability range for each urban patch: The cumulative cost algorithm calculates the minimum cumulative cost required to reach each point in space from the target location. The raster image cell values ​​obtained through this algorithm represent the minimum cumulative cost to reach the cell from the target location. The 2015–2022 motion speed raster is converted into a cost raster using GEE according to the following formula 7: Formula 7, In the formula, The converted pixel value, The values ​​are the pixel values ​​before conversion. Each city patch after classification and the cost grid are sequentially input into the cumulative cost algorithm to obtain the cumulative cost grid for each city patch. The pixel value represents the time required to travel from the city patch to each point in space. A mask is applied with a 1-hour threshold to obtain the reachable range of each city patch within 1 hour.

[0026] S202. Fragmentation index scoring, specifically including the following steps: S20201. Generating an 8-directional buffer zone: In each city, the 1-hour reachable area of ​​a city patch is used as the city edge zone. A 1 km concentric circular buffer zone is generated from the city patch boundary to the city edge zone boundary. As the distance from the city patch boundary gradually increases, the fragmentation of forests and farmland tends to gradually improve, exhibiting a spatial gradient. Therefore, a 1 km concentric circular buffer zone is generated within the city edge zone of each city patch, and the fragmentation of forests and farmland is quantified within this buffer zone. Simultaneously, because the gradient of forest and farmland fragmentation changes differently in different directions, the concentric circular buffer zone is divided into 8 sectors, each sector expanding radially outward at 45°, as detailed below. Figure 2 As shown.

[0027] S20202. Calculation of Landscape Indicators for Farmland and Forest Fragmentation: Due to urban expansion, the landscape pattern of farmland and forest around cities is gradually becoming fragmented. Landscape indices are quantitative indicators that condense information about landscape patterns, reflecting the structural composition and spatial configuration of the landscape. Therefore, landscape indices are used to quantify the degree of farmland and forest fragmentation in urban fringe areas. The GLC_FCS30D data of the target area is input into Fragstas software to calculate the landscape indices within each concentric buffer zone. Specifically, the following indicators are used at the categorical scale: Patch Density (PD), Maximum Patch Index (LPI), Average Patch Area (AREA_MN), Aggregation Index (AI), Patch Connectivity Index (COHESION), and Average Fractal Dimension (FRAC_MN). 1) The formula for the patch density index is shown in Formula 8: Formula 8, In the formula, For patch density, This represents the number of patches of a specific land cover class in the landscape. This represents the total area of ​​the landscape. The higher the value of this indicator, the greater the patch density, and the more fragmented the landscape.

[0028] 2) The formula for the maximum plaque index is shown in Formula 9: Formula 9, In the formula, The maximum plaque index, This represents the area of ​​the largest patch in the landscape. This index is used to determine the dominant patch type in the landscape. The larger the index, the greater the proportion of the patch in the landscape and the less fragmented it is.

[0029] 3) The formula for the average patch area is shown in Formula 10: Formula 10, In the formula, The average patch area, This represents the area of ​​a specific land cover category within the landscape. This represents the number of patches of a specific land cover category within a landscape. A higher value indicates a less fragmented landscape type.

[0030] 4) The formula for the aggregation index is shown in Formula 11: Formula 11, In the formula, The aggregation index, This represents the number of adjacent patches of the same land cover category for a given patch. This represents the maximum number of adjacent patches of the same land cover category within a given land cover category. This index is used to calculate the degree of clustering of patches within a specific land cover category. A higher value indicates greater clustering of patches within that category and a less fragmented landscape.

[0031] 5) The formula for the plaque connectivity index is shown in Formula 12: Formula 12, In the formula, This is the plaque connectivity index. The perimeter of the patch, expressed in pixels. The area of ​​the patch is expressed in pixels. This represents the total number of pixels in the landscape. This index is used to calculate the connectivity of patches in the landscape; a higher value indicates better connectivity and less fragmentation.

[0032] 6) The formula for the average fractal dimension is shown in Formula 13: Formula 13, In the formula, The average fractal dimension, The first in a certain land cover category The perimeter of each patch, The first in a certain land cover category The area of ​​each patch, This represents the total number of patches within a specific land cover category. This index is used to calculate the shape complexity of patches in a landscape; a higher value indicates a more complex patch shape.

[0033] S20203. Fragmentation Score: The degree of fragmentation is scored based on the landscape index calculation results of cultivated land and forest in each concentric buffer zone. The specific formula is shown in Formula 14-17. Formula 14, Formula 15, Formula 16, Formula 17, In the formula, , For the first The density of forest and farmland patches within concentric circular buffer zones , For the first The maximum patch index of forests and cultivated land within concentric circular buffer zones , For the first The average patch size of forest and cultivated land within a concentric circular buffer zone , For the first Forest and arable land aggregation indices within concentric buffer zones , For the first The connectivity index of forest and cultivated land patches in a concentric circular buffer zone , For the first The average fractal dimension of forests and farmland within a concentric circular buffer zone. , For the first The scores for the degree of forest and farmland fragmentation within concentric circular buffer zones. For the first The overall fragmentation score in a concentric circular buffer zone.

[0034] Since the urban fringe areas of different city patches overlap, the fragmentation score of the overlapping parts is handled according to the following rules: (1) If the urban fringe areas of high-level cities and low-level cities overlap, the fragmentation score of the overlapping area is the fragmentation score of the area in the urban fringe area of ​​the high-level city; (2) If the urban fringe areas of cities of the same level overlap, the fragmentation score of the overlapping area is the maximum value of the fragmentation score of the area.

[0035] S3. Human Footprint Index Calculation: Select land cover data, redistributed population grid data, railway data, highway data, nighttime light data, and navigable waterway data to calculate land use indicators, population density indicators, road accessibility indicators, power infrastructure indicators, and navigable waterway indicators. Add these scores to obtain the original human footprint score. Then, add the fragmentation indicator score to the original human footprint score to obtain the human footprint index.

[0036] Furthermore, this embodiment selects GLC_FCS30D land cover data, redistributed LandScan population grid data, OSM railway data, OSM highway data, NPP / VIIRS data, and HydroSHEDS river flow data. The Human Footprint Index is calculated using an assignment method, named the Pre-Improved Human Footprint Index. The Improved Human Footprint Index is obtained by adding a fragmentation score to the Pre-Improved Human Footprint Index. The specific indicators of the Pre-Improved Human Footprint Index are shown below: (1) Land use indicators

[0037] The land use index includes different land cover types. Among them, the built-up area is most affected by human activities and is assigned 10 points; cultivated land has a strong impact on the natural surface and soil physical and chemical properties, so it is assigned 7 points; the remaining land cover types are all assigned 0 points.

[0038] (2) Population density index

[0039] Population density is represented using optimized population grid data. Existing population density assignment methods are used, assigning values ​​greater than 1000 people / km². 2 The pixel value is assigned 10 points, which is less than 1000 people / km. 2 The pixels are assigned according to the following formula 18: Formula 18, In the formula, For the first Population density is assigned to each pixel. For the first Population density per pixel.

[0040] (3) Road accessibility index

[0041] Road data is used to represent road accessibility. The pressure exerted by roads can be divided into pressure exerted by railways and pressure exerted by highways, and scored according to formulas 19-20: Formula 19, Formula 20, In the formula, For the first Each pixel represents a railway pressure score. For the first The distance of one pixel from the railway. For the first Each pixel represents a highway pressure score. For the first The distance of each pixel from the highway.

[0042] (4) Power infrastructure indicators

[0043] The nighttime light data for each year is normalized and then assigned a score of 0–10. See Formula 21 for details: Formula 21, In the formula, For the first The pressure on power infrastructure is assigned to each pixel. For the first The pixel value of a pixel. This represents the maximum pixel value.

[0044] (5) Navigation waterway indicators

[0045] The navigable waterway calculated as described above is used to calculate this index, and the navigable waterway pressure is scored according to the following formula 22: Formula 22, In the formula, Assigning pressure scores to navigable waterways For the first Distance of each pixel from the navigable waterway.

[0046] The improved Human Footprint Index calculation involves summing land use scores, population density scores, road accessibility scores, nighttime light scores, and navigable waterway pressure scores to obtain the original human footprint score. A fragmentation score is then added to this original score to arrive at the improved Human Footprint Index, as shown in formulas 23-24. Formula 23, Formula 24, In the formula, For the first The improved human footprint assignment per pixel before the pixel was used. For the first Improved human footprint scoring per pixel For the first Land use allocation per pixel For the first Population density is assigned to each pixel. For the first Highway pressure score assigned to each pixel For the first The railway pressure score is assigned to each pixel. For the first The nighttime lighting assignment of each pixel For the first The pressure assignment of each pixel in the navigable waterway For the first The degree of fragmentation of each pixel.

[0047] Human footprint accuracy verification: The required sample size is calculated using the minimum sample size calculation formula. After calculating the minimum sample size, random sampling is then performed. 1 km 2 Sampling points were used to verify the accuracy of human footprints. Visual interpretation was employed to assess the intensity of human activity at the sampling points, which served as the validation dataset. The visual interpretation rules were as follows: within a 1 km... 2 The landscape within Google's high-resolution imagery was manually visually interpreted. Visual interpretation stress scores were assigned based on the landscape composition. The landscapes requiring scoring mainly included urban areas, farmland, roads, deforestation areas, non-urban residential areas, and infrastructure. A single image may contain multiple landscape types. The final visual interpretation stress score is the sum of the scores for each landscape type. The landscape scoring rules are shown in Table 3. Table 3 Landscape Visual Interpretation Stress Assignment Rules

[0048] The visual interpretation results and the improved human footprint indices obtained in this application are normalized. The root mean square error (RMSE) and goodness-of-fit (R²) are used to calculate the normalized improved human footprint indices. The Kappa coefficient is used to calculate their consistency. When calculating the Kappa coefficient, if the human footprint score differs from the visual interpretation result by less than 20% in the 0–1 interval, it is considered a match to the visual interpretation result. Specific formulas are shown in Formulas 25-28. Formula 25, Formula 26, Formula 27, Formula 28, In the formula, For the first Human footprint score for each verification patch For the first Visual interpretation of human activity intensity scores for each verification patch. For the first The number of samples that are correctly classified into each category. For the first The number of true samples in each category For the first The number of samples in each category For the sample size, This represents the number of categories.

[0049] In this embodiment, the method of the invention is applied to the Indochina Peninsula (located between 0° and 29°N, and 92° and 110°E), with data collection conducted from 2015 to 2022. The mapping results are as follows: Figure 3 As shown, the map was created using the improved method for obtaining the Human Footprint Index, i.e., the method of this invention. Human footprint values ​​are higher in the southern, western, and eastern parts of the Indochina Peninsula, and lower in the northern part. The areas with the highest human footprint values ​​(red areas) are mainly distributed around large cities and coastlines. Figure 4 As shown, taking 2022 as an example, Figure 4 The left-hand map shows the map created using the method before the improvement, while the right-hand map shows the map created using the method after the improvement. Regions A, B, and C in both maps represent Hanoi, Bangkok, and Ho Chi Minh City, respectively, all major cities in the Indochina Peninsula. Figure 5 As shown in the high-resolution imagery from Google Earth, these cities are surrounded by relatively fragmented farmland and woodland (green areas). In the mapping based on the original method of obtaining the Human Footprint Index, the human footprint values ​​in these areas were relatively low. However, in the mapping based on the improved method, the human footprint values ​​in these areas increased compared to the former. The improved method incorporates indicators such as the degree of fragmentation to characterize actual human activity, achieving accurate quantification of the intensity of human activity in these areas. The increase in values ​​is highly consistent with multi-source verification data such as visual interpretation and nighttime light, and it also makes the human footprint index around the city show a reasonable spatial pattern of gradient decay from the built-up area outwards, effectively correcting the identification bias of the original method and significantly improving the accuracy of the index.

[0050] In summary, this application addresses the issues of inaccurate population data distribution and more fragmented farmland and forest land in urban fringe areas by improving the calculation process of the human footprint index. It delineates urban fringe areas by the reachability of urban patches within one hour, quantifies the degree of fragmentation of farmland and forest land in urban fringe areas, and incorporates a fragmentation index into the human footprint index, thereby improving the accuracy of human footprint remote sensing mapping.

[0051] The preferred embodiments of the present invention have been described above with reference to the accompanying drawings. Those skilled in the art can implement the present invention in various modifications without departing from its scope and spirit. For example, a feature shown or described in one embodiment can be used in another embodiment to obtain yet another embodiment. The above are merely preferred embodiments of the present invention and are not intended to limit the scope of the present invention. All equivalent changes made based on the description and drawings of the present invention are included within the scope of the present invention.

Claims

1. A method for improving the accuracy of the human footprint index, characterized in that, Specifically, the following steps are included: S1. Population grid redistribution: Construction land mask is extracted from land cover data, and population statistics data and population grid data of construction land mask are input into the iURBAN algorithm to realize population grid data redistribution. Population statistics data are used to verify the accuracy of population grid data before and after redistribution. S2. Fragmentation Index Scoring: Navigable waterways are calculated using river flow data, and movement speed grids are calculated using road data, railway data, navigable waterway data, and land cover data; impervious patches are obtained based on land cover data, and these patches are classified into city tiers based on the redistributed population; urban fringe areas of different city tiers are obtained based on city tier classification and movement speed grids, and landscape pattern indicators representing fragmentation are calculated within the urban fringe areas of the urban patches; fragmentation index scoring is then assigned based on the landscape pattern indicator calculation results. S3. Human Footprint Index Calculation: Select land cover data, redistributed population grid data, railway data, highway data, nighttime light data, and navigable waterway data to calculate land use indicators, population density indicators, road accessibility indicators, power infrastructure indicators, and navigable waterway indicators. Add these scores to obtain the original human footprint score. Then, add the fragmentation indicator score to the original human footprint score to obtain the human footprint index.

2. The method for improving the accuracy of the human footprint index according to claim 1, characterized in that, The population grid data reallocation mainly includes the following steps: S101. Based on land cover data and population grid data, obtain the area size and population of each built-up area patch, and determine the importance of the built-up area patch through the two. After comparing with the urban population in the statistical yearbook, the built-up area patches are divided into urban built-up area patches and rural built-up area patches. S102. Based on the size of the rural built-up area patches and the population of the rural built-up area patches, allocate the rural population data in the statistical yearbook to each rural built-up area patch. S103. Based on the proportion of built-up area and population of the pixels within each rural built-up area patch, the statistical yearbook population allocated to each rural built-up area patch is assigned to the pixels within the patch.

3. The method for improving the accuracy of the human footprint index according to claim 1, characterized in that, The accuracy of the population grid data before and after redistribution was evaluated using root mean square error (RMSE) based on demographic data.

4. The method for improving the accuracy of the Human Footprint Index according to claim 1, characterized in that, Step S2 also includes quantifying the city's reachability within one hour, specifically including the following steps: S20101, City Classification: Statistically determine the population within impermeable patches and classify impermeable patches into large cities, medium cities, and small cities based on the population size. S20102. Create a movement speed grid: Create corresponding movement speed grids based on roads, railways, navigable waterways, land cover, altitude, and slope. The resolution of each movement speed grid is 1000 m, and the cell value represents the movement speed in each cell. S20103. The cumulative cost algorithm is used to obtain the reachable range of each city patch within 1 hour.

5. The method for improving the accuracy of the Human Footprint Index according to claim 1, characterized in that, The scoring of the fragmentation index also includes the following steps: S20201. Generate 8-directional buffer zones: Generate a 1 km concentric circular buffer zone within the patch city edge area of ​​each city, and dequantize the fragmentation of forests and farmland within the buffer zone; Since the gradient of the fragmentation of forests and farmland varies in different directions, the concentric circular buffer zone is divided into 8 sectors, each sector expanding radially outward at 45°. S20202, Calculation of landscape indicators for farmland-forest fragmentation; the landscape indicators include patch density, maximum patch index, average patch area, aggregation index, patch connectivity index, and average fractal dimension. S20203. Fragmentation score: The degree of fragmentation is scored based on the landscape index calculation results of cultivated land and forest in each concentric buffer zone.

6. The method for improving the accuracy of the human footprint index according to claim 5, characterized in that, The formula for the patch density index is as follows: In the formula, For patch density, This represents the number of patches of a specific land cover class in the landscape. Total landscape area; The formula for the maximum plaque index is as follows: In the formula, The maximum plaque index, This indicates the area of ​​the largest patch in the landscape; this indicator is used to determine the dominant patch type in the landscape. The formula for the average patch area is as follows: In the formula, The average patch area, This represents the area of ​​a specific land cover category within the landscape. The number of patches in a landscape for a specific land cover category; The formula for the aggregation index is as follows: In the formula, The aggregation index, This represents the number of adjacent patches of the same land cover category for a given patch. This represents the maximum number of adjacent patches of the same land cover category within a given land cover category. The formula for the patch connectivity index is as follows: In the formula, This is the plaque connectivity index. The perimeter of the patch, expressed in pixels. The area of ​​the patch is expressed in pixels. This represents the total number of pixels in the landscape. The formula for the average fractal dimension is as follows: In the formula, The average fractal dimension, The first in a certain land cover category The perimeter of each patch, The first in a certain land cover category The area of ​​each patch, This represents the total number of patches within a given land cover category.

7. The method for improving the accuracy of the human footprint index according to claim 5, characterized in that, When the urban edge areas of different urban patches overlap, the fragmentation score of the overlapping parts is processed according to the following rules: (1) If the urban fringe areas of high-level cities and low-level cities overlap, the fragmentation score of the overlapping area shall be the fragmentation score of the urban fringe area of ​​the high-level city. (2) If the urban fringe areas of cities of the same level overlap, the fragmentation score of the overlapping area shall be the maximum value of the fragmentation score of that area.

Images