Post-earthquake recovery and its driving forces of ecological environment quality using remote sensing and GIScience, a case study of 2015 Ms8.1 Nepal earthquake

Abstract

Natural hazards always have significant impacts on the regional ecological environment, causing huge damages to biodiversity, soil erosion, and impoverishment. This study aimed to analyze the variation of ecological environment quality in Nepal from 2014 to 2018, with a specific focus on the impact and recovery of the 2015 Ms8.1 Nepal earthquake. We utilized Google Earth Engine (GEE) and remote sensing based ecological index (RSEI) to evaluate the ecological environment quality. Moreover, Moran’s I and Geodetector model were employed to explore spatial clustering and driving forces of post-earthquake recovery. The results indicate that after the 2015 Nepal earthquake, there was a decrease of 0.0225 in the mean RSEI in Nepal. The monitoring of the severe-affected areas showed that the regional environmental quality gradually improved within three years after the earthquake. After a span of three years, the ecological environment managed to recover to its pre-earthquake level. In terms of spatial distribution, there is a notable correlation in Nepal, with the ecological environment exhibiting higher quality in the southern regions and lower quality in the northern regions. By examining the driving forces behind the changes in ecological environment quality during the post-earthquake recovery period, it has been determined that altitude, precipitation, and their interaction serve as influential factors. These findings provide valuable insights into the state of ecological environment recovery after the 2015 Nepal earthquake. Furthermore, this study establishes a robust theoretical foundation for local decision-makers to develop ecological protection measures and enhance local disaster risk management.

Keywords:

• Remote sensing
• recovery
• RSEI
• geodetector
• Moran’s I index
• the 2015 Ms8.1 Nepal earthquake

1. Introduction

To this day, a wide range of natural hazards continues to occur frequently, including storms, floods, tsunamis, earthquakes, mudslides, and more (Liu, Liu, Lin, et al. 2023; Liu, Zhang, et al. 2023). These events often lead to ecological deterioration, building damage, economic losses, and psychological issues for the affected population (Kane et al. 2018; Liu et al. 2021). The process of post-disaster recovery is both complex and lengthy (Gao, Chen, et al. 2020). Among these natural hazards, earthquakes stand out as one of the most destructive geological phenomena, inflicting significant and enduring negative effects (Liu et al. 2019, 2020; Li et al. 2020). The effective recovery after such natural hazards holds paramount importance for local human survival and sustainable social development. Despite its crucial significance, the post-disaster recovery stage is often overlooked or neglected within the various stages of disaster risk management (Ghaffarian et al. 2021). This, in turn, hampers the provision of crucial information to decision-makers and planners regarding the extent of damage and the recovery process after a natural hazard.

With the advent of remote sensing technology, the availability and analysis of large-scale, real-time ground information have greatly improved, leading to its widespread utilization in disaster risk management. For instance, Lin et al. (2004) conducted an evaluation of vegetation recovery in landslide-prone areas resulting from the 921 earthquake in central Taiwan, employing SPOT satellite images. Contreras et al. (2016) utilized Italy as a case study, effectively integrating remote sensing, geographic information system (GIS), and ground observations to assess the progress of post-earthquake recovery. Additionally, Cooner, Shao, and Campbell (2016) employed remote sensing technology in conjunction with machine learning algorithms to detect damage to urban areas after the 2010 Haiti earthquake. These studies showcase the significant role that remote sensing plays significance in monitoring and analyzing the impacts of earthquakes and guiding post-disaster recovery efforts.

In ecological environment assessment studies, the development of remote sensing technology has also made satellite imagery a valuable tool for collecting data in the post-disaster phase (Lin et al. 2004; Ghaffarian, Kerle, et al. 2018; Ghaffarian, Rezaie Farhadabad, et al. 2020; Ghaffarian et al. 2021). So far, various remote sensing evaluation indexes have been developed to characterise or quantify ecological quality status (Yang et al. 2021). Such as enhanced vegetation index (EVI) (Wang et al. 2008), land surface temperature (LST) (Liao and Jiang 2020), normalized difference build-up index (NDBI) (Luo and Li 2014), and hyperspectral flower index (HFI) (Chen et al. 2009). Some researchers employ the above or other single indexes for evaluation, but this often cannot truly reflect the complex ecological environment of the region. Xu (2013) proposed the RSEI integrating greenness, humidity, dryness, and temperature, which can quickly and objectively reflect the regional ecological status. Compared with the evaluation of a single index, RSEI can more accurately quantify ecological quality (Xiong et al. 2021) and has been widely used (Hu et al. 2020; Yan et al. 2021; Cui, Zhu, Mi, et al. 2022; Yang, Meng, Fu, Zhang, et al. 2022; Yang and Li 2023). Meanwhile, the emergence of the Google Earth Engine (GEE) platform has brought a lot of convenience to remote sensing image processing (Geng et al. 2022; Yang, Tian, et al. 2022). It has a geographic database that stores massive images (Zhao et al. 2021; Cui, Zhu, Liang, et al. 2022; Yang, Meng, Fu, Wang, et al. 2022; Campos et al. 2023), including Landsat, MODIS, meteorological and socioeconomic data sets, etc. The GEE platform overcomes the low efficiency of traditional remote sensing image processing by writing simple programs online, directly accessing the database, and quickly processing and analyzing geospatial data (Kumar and Mutanga 2018; Amani et al. 2020; Waleed et al. 2023). For example, based on the GEE platform and the RSEI, Zhang, She, et al. (2022) conducted a 20-year long-term analysis of the Chang-Zhu-Tan metropolitan circle. The results found that the ecological environment level in this area is mainly good and excellent, and it is on the rise. Huang et al. (2021) also used the GEE cloud computing platform to build the RSEI to assess and monitor the changes in the quality of the ecological environment in the Lhasa metropolitan area. In addition to discovering ecological changes, the results also confirmed that GEE is an efficient platform that can be used for large-scale environmental monitoring.

On the basis of RSEI, it can reflect the changing trend of regional ecological quality, and in addition, the driving forces of image RSEI can be explored. Geodetector is a model that combines GIS spatial overlay technology with set theory (Wang and Xu 2017; Cui, Zhu, Liang, et al. 2022), which can reveal the driving factors behind their changes, and can also quantitatively describe the interaction between driving factors (Gao, Kasimu, et al. 2020; An et al. 2022; Wang et al. 2023; Yang and Su 2023), so it is widely used in the ecological quality assessment field. For instance, Cui, Zhu, Liang, et al. (2022) delved into the driving forces behind land use/land cover changes within the Shandong section of the Yellow River Basin over the period from 2000 to 2020. Zhu et al. (2020) utilized Geodetector to analyze the influential effects of both natural and human factors on NDVI alterations in the middle reaches of the Heihe River Basin. Likewise, Huo and Sun (2021) employed Geodetector as a foundation to explore the dynamic spatiotemporal variations in vegetation and the factors influencing these changes, selecting the northwest Yunnan Plateau as their case study area. The implementation of Geodetector in these studies showcased their efficacy in unraveling the complexities of geographical phenomena and provided valuable insights for understanding environmental dynamics and decision-making processes.

Situated at the foothills of the Himalayas, Nepal stands among the countries frequently struck by natural hazards (Aryal 2014). On April 25, 2015 (11:56 am local time), the epicenter was located in the Gorkha District, and a Ms8.1 earthquake occurred (Okamura et al. 2015). Nepal’s complex terrain exacerbates the impact of earthquakes and their secondary disasters on the regional ecological environment. The preservation of a healthy ecological environment is of utmost significance to humanity. Consequently, meticulous monitoring of ecological changes and post-earthquake recovery becomes imperative. However, there is still a lack of understanding regarding the changes in the quality of the ecological environment after the earthquake.

In this study, we first employed the GEE platform as the foundation for utilizing remote sensing data to develop the RSEI, facilitating a swift and accurate assessment of the ecological quality status before and after the 2015 Ms8.1 Nepal earthquake. Moreover, conducting raster analysis on the severe-affected areas enhances our comprehension of post-earthquake recovery dynamics. After that, Moran’s I was employed to delve into the spatial autocorrelation analysis of the ecological environment in Nepal. Finally, the Geodetector model was introduced, considering both natural and social factors, to investigate the driving forces behind changes in ecological quality. Additionally, this study possessed the potential to make significant contributions to the advancement of more comprehensive and effective methods for ecological quality assessment and disaster risk management in Nepal. By combining remote sensing technology and GIScience, this study significantly strengthened the ongoing efforts towards enhancing resilience and promoting sustainable development in the aftermath of seismic events.

2. Study area and data

2.1. Study area

Nepal, located in South Asia (26°22′N-30°27′N, 80°04′E-88°12′E) (Pangali Sharma et al. 2019), is a landlocked country bordered by the Tibet Autonomous Region of China in the north and India in the east, south, and west (shown in Figure 1). The country has a total area of 147,181 square kilometers, spanning 800 kilometers from east to west and 200 kilometers from north to south. With altitudes ranging from 45 meters in the south to 8,729 meters in the north, Nepal’s terrain varies widely, with the south characterized by low-lying areas and the north marked by high peaks. Nepal is divided into five climatic zones that roughly correspond to altitude. These include the tropical and subtropical zones below 1,200 m, the temperate zone between 1,200 and 2,400 m, the cold zone between 2,400 and 3,600 m, the subarctic zone between 3,600 and 4,400 m, and the Arctic zone above 4,400 m. And, Nepal experiences five seasons: summer, monsoon, autumn, winter, and spring. Annual precipitation varies from 200 mm in parts of the northern Himalayas to 5,500 mm in the south, with a national average of 1,800 mm. Meanwhile, average temperatures in the south can reach as high as 30 °C, whereas in the north they can drop as low as −10 °C. In summary, Nepal is a diverse country with varying topography and climate. Its climatic zones and seasons are largely determined by altitude, and the country experiences a wide range of temperatures and precipitation levels across its different regions.

Figure 1. The location of Nepal.

On April 25, 2015, a devastating Ms8.1 earthquake struck Nepal, known as the Gurkha Earthquake (Liu, Fang, et al. 2021). Strong aftershocks with Ms7.0, Ms7.1, and Ms7.5 occurred at 34 min, 24 h, and 17 days, respectively (presented in Figure 2). The 2015 Nepal earthquake resulted from the low-angle subduction of the Indian plate towards the Eurasian plate along the Main Himalayan Thrust (MHT) fault zone. The Himalayan tectonic belt is mainly composed of three giant thrust fault systems: the Main Central Thrust (MCT) Zone, the Main Boundary Thrust (MBT) Zone and the Main Front Thrust (MFT) Zone from north to south, and the fault zone has a high frequency of seismic activity. Among them, there have been additional earthquakes in Nepal, such as the 1505 Ms8.2 Lo Mustang earthquake, the 1833 Ms8.0 Kathmandu earthquake, and the 1934 Ms8.1 Bihar earthquake. Of note, the 2015 Nepal earthquake stands as the largest recorded earthquake in Nepal since 1934 (Okamura et al. 2015) and left an indelible impact on the affected regions. More than 8,790 people died in the earthquake and 22,000 were injured (Earthquake 2015). The total economic loss is around $10 billion, or about half of Nepal’s gross domestic product, the report said (Goda et al. 2015). The earthquake’s epicenter was located in Gurkha District, approximately 80 kilometers northwest of the capital city, Kathmandu. The earthquake had a focal depth of 15 kilometers and reached a high seismic intensity level of Ⅹ in the severe-affected areas. The rupture of the earthquake extended from west to east (Goda et al. 2015), resulting in powerful tremors that were felt in Kathmandu and the surrounding cities. The Central Development Region and Western Development Region of Nepal, which were former administrative divisions, were among the severe-affected areas by the earthquake.

Figure 2. Seismic intensity map of 2015 Nepal earthquake

(Sun and Yan, 2015).

2.2. Data collection

Multi-source data was employed in this study, including remote sensing data, digital elevation models (DEM), precipitation data, population density data, and nighttime lights data. Remote sensing data played a crucial role in this study, and the MODIS dataset was selected due to its high temporal resolution, making it well-suited for large-scale monitoring purposes (Xu et al. 2019; Lu et al. 2022). Specifically, the MODIS dataset’s MOD09A1, MOD11A2, and MOD13A1 subsets were selected for preprocessing on the GEE platform, with July to September being used to generate four indicators and the RSEI. Within this timeframe, T0 represents the pre-earthquake period in 2014, T1 represents the year 2015 just after the earthquake, and T2 to T4 represents the subsequent post-disaster recovery stage from 2016 to 2018. Apart from remote sensing data, this study incorporated DEM and precipitation data, which were also processed using the GEE platform. The elevation data was utilized to calculate slope and aspect information. As a result, six factors—elevation, slope, aspect, precipitation, nighttime lights, and population density—were employed to explore the driving forces behind the recovery of the RSEI. Table 1 presents comprehensive details regarding the study data, while Figure 3 illustrates the spatial distribution of the selected impact factors.

Figure 3. Spatial distribution map of impact factors—taking T2 time as an example.

Table 1. Details of datasets.

Data	Time	resolution	Time resolution	Source
Administrative division	/	/	/	https://gadm.org/data.html
MOD09A1 v061	2014-2018 (7.01-9.25)	500m	8-day	https://www.earthdata.nasa.gov/
MOD11A2 v061		500m	8-day
MOD13A1 v061		500m	16-day
DEM	2016-2018	30m	/	https://www.earthdata.nasa.gov/
Precipitation		0.1 degree	30min	https://www.earthdata.nasa.gov/
Population density		1000m	/	https://www.worldpop.org/
Nighttime lights		500m	year	https://engine-aiearth.aliyun.com/#/

3. Methodology

The flowchart of this study is shown in Figure 4. This study used MODIS datasets as the basis for calculating four component indexes—greenness, dryness, humidity, and temperature—using the GEE cloud computing platform. Principal component analysis (PCA) was then conducted to generate the RSEI to monitor the spatiotemporal changes of the ecological environment quality before and after the 2015 Nepal earthquake. In order to explore the recovery of the ecological environment quality after the earthquake, a grid analysis was carried out for the RSEI in the severe-affected areas. Additionally, we analyzed the spatial distribution of RSEI from T0 to T4 with the help of Moran’s I index. Finally, the driving factors behind the changes in RESI during the recovery period were explored based on the Geodetector model.

Figure 4. Workflow of this study.

3.1. Methods of RSEI

RSEI is a quantitative evaluation method of the ecological environment that integrates greenness, humidity, dryness, and temperature and is completely based on remote sensing information (Xu et al. 2019). To represent the greenness index, we chose NDVI, which is closely related to plant biomass, leaf area index, and vegetation coverage (Goward et al. 2002) and is the most widely used vegetation index (de Araujo Barbosa et al. 2015; Mishra et al. 2015). The humidity index in this study is represented by the third component WET in the tasseled cap transformation, which reflects the water content of water, soil, and vegetation (Crist 1985). In urban areas, normalized difference built-up and soil index (NDBSI) constructed by bare soil index and building index is used to represent dryness since building land can also have an impact on the dryness of the surface. Lastly, LST is the most common index of urban heat and is employed as the temperature index in this study (Nichol 2005). The RSEI can be expressed as a function of these four indicators as follows: (1) $$\mathit{RSEI}=f(\mathit{NDVI},\mathit{WET},\mathit{NDBSI},\mathit{LST})$$(1)

We calculate the four indicators mentioned above using MODIS data and GEE platform. NDVI is calculated from the MOD13A1 dataset, while LST is calculated from the MOD11A2 dataset. WET and NDBSI are calculated by MOD09A1. The calculation formula of WET (Lobser and Cohen 2007) and NDBSI (Yuan et al. 2021) is as follows: (2) $$\mathit{WET}=0.1147{\rho }_1+0.2489{\rho }_2+0.2408{\rho }_3+0.3132{\rho }_4-0.3122{\rho }_5-0.6416{\rho }_6-0.5087{\rho }_7$$(2) (3) $$\mathit{NDBSI}=\frac{1}{2}\frac{\frac{2{\rho }_6}{{\rho }_6+{\rho }_2}-\left[\frac{{\rho }_2}{{\rho }_2+{\rho }_1}+\frac{{\rho }_4}{{\rho }_4+{\rho }_6}\right]}{\frac{2{\rho }_6}{{\rho }_6+{\rho }_2}+\left[\frac{{\rho }_2}{{\rho }_2+{\rho }_1}+\frac{{\rho }_4}{{\rho }_4+{\rho }_6}\right]}+\frac{1}{2}\frac{\left[\right({\rho }_6+{\rho }_1)-({\rho }_2+{\rho }_3\left)\right]}{\left[\right({\rho }_6+{\rho }_1)+({\rho }_2+{\rho }_3\left)\right]}$$(3)

Where ${\rho }_1$ to ${\rho }_7$ denote bands 1 to 7 corresponding to the red, NIR1, blue, green, NIR2, short-wavelength infrared 1 (SWIR1), and SWIR2 bands of the MODIS image, respectively.

To avoid any weight imbalances among the four indicators mentioned above, which have different dimensions, normalization was necessary before conducting the PCA. This ensured that the range of all indexes was unified between 0 and 1 (Jing et al. 2020). After normalization, a PCA was performed, and the first principal component (PC1) was extracted to calculate the RSEI. For comparative analysis, the RSEI is normalized again as follows (Xiong et al. 2021). After normalization, RSEI is evenly divided into five levels at intervals of 0.2, including Poor: 0-0.2; Fair: 0.2-0.4; Moderate: 0.4-0.6; Good: 0.6-0.8 and Excellent: 0.8-1. (4) $${NI}_i=(I_i-l_{\min })/(l_{\max }-l_{\min })$$(4) (5) $${\mathit{RSEI}}_0=PC1\left[f\right(\mathit{NDVI},\mathit{WET},\mathit{NDBSI},\mathit{LST}\left)\right]$$(5) (6) $$\mathit{RSEI}=\left({\mathit{RSEI}}_0-{\mathit{RSEI}}_{0\_\min }\right)/\left({\mathit{RSEI}}_{0\_\max }-{\mathit{RSEI}}_{0\_\min }\right)$$(6)

Where ${NI}_i$ is the image normalization result; $I_i$ is the numerical pixel value of the image; $l_{\max }$ and $l_{\min }$ are the maximum and minimum values of the image pixel, respectively; $PC1$ is the first principal component of the four indicators; ${\mathit{RSEI}}_0$ indicates the initial ecological index value; ${\mathit{RSEI}}_{0\_\max }$ and ${\mathit{RSEI}}_{0\_\min }$ are the maximum and minimum ${\mathit{RSEI}}_0$ value, respectively.

3.2. Spatial correlation analysis

3.2.1. Global Moran’s I

Global Moran’s I is a widely used statistic in global correlation analysis. It is primarily utilized to describe the average degree of correlation between all spatial units and their surrounding areas within the entire study area (Zhang, Zhao, et al. 2022), and the formula is as follows: (7) $$\mathit{Global}\mathrm{ }\mathit{Moran}\mathrm{’}s\mathrm{ }I=\frac{n{\sum }_{i=1}^n{\sum }_{j=1}^n{W}_{ij}(x_i-\overline{x})(x_j-\overline{x})}{{\sum }_{i=1}^n{\sum }_{j=1}^n{W}_{ij}{\sum }_{i=1}^n{{(x}_i-\overline{x})}^2}$$(7)

Where $n$ is the number of spatial elements, $x$ presents the variable of the study object, $W_{ij}$ donates the binary weight matrix of the general cross-production statistic, $\overline{x}$ means the average of the property values of the study object.

The value range of $\mathit{Global}\mathrm{ }\mathit{Moran}\mathrm{’}s\mathrm{ }I$ is $[-1,\mathrm{ }1].$ When $I>0,$ means that all regional attribute values have a positive correlation in space; when $I=0,$ indicates that all regional attribute values are randomly distributed without spatial correlation; when $I<0,$ presents that all regional attribute values have negative spatial correlation (Luo et al. 2022).

3.2.2. Local Moran’s I

Local Indicators of Spatial Association (LISA) can further determine the degree of correlation between the RSEI of each spatial unit and the RSEI of its neighbors (Anselin 1995). Meanwhile, it can indicate the location and degree of its spatial aggregation or dispersion. The formula of the Local Moran’s I is as follows: (8) $$\mathit{Local}\mathrm{ }\mathit{Mora}{n}^{\prime }s\mathrm{ }I=\frac{n(x_i-\overline{x}){\sum }_{j=1}^n{W}_{ij}(x_j-\overline{x})}{{\sum }_{i=1}^n{(x_i-\overline{x})}^2}$$(8)

Where $x_i,$ $x_j,$ $x_j$ , and $W_{ij}$ involved in the calculation formula all have the same parameters as formula (7(7) $$\mathit{Global}\mathrm{ }\mathit{Moran}\mathrm{’}s\mathrm{ }I=\frac{n{\sum }_{i=1}^n{\sum }_{j=1}^n{W}_{ij}(x_i-\overline{x})(x_j-\overline{x})}{{\sum }_{i=1}^n{\sum }_{j=1}^n{W}_{ij}{\sum }_{i=1}^n{{(x}_i-\overline{x})}^2}$$(7) ).

3.3. Geodetector model

The Geodetector is a statistical method that is used to detect spatial heterogeneity and reveal its driving factors by detecting spatial differentiation (Wang et al. 2010; Wu et al. 2015). This model has been extensively employed in various fields, such as geography and ecology. Here, we focus on the detection of the factors influencing ecological environment changes and their interactions using the Geodetector model.

1. Factor detector

   This module is capable of quantitatively detecting the spatial heterogeneity of the dependent variable and the explanatory power of the influencing factor $X_i$ on the spatially heterogeneous characteristics of RSEI changes (Zhang, Jiang, et al. 2022); $q$ can be expressed mathematically as follows (Wang et al. 2021; Cai et al. 2023; Liu, Liu, Qiu, et al. 2023):

(9) $$q=1-\frac{\sum _{h=1}^L{N}_h{\sigma }_h^2}{N{\sigma }^2}$$(9) where $L$ is the number of layers of the independent variable; $h$ is the stratification of the category number of driving factor $X;$ $N_h$ and $N$ are the numbers of units in layer $h$ and the whole region, respectively; ${\sigma }^2$ is the overall variance of the sample. The value of $q$ is consistent with the explanatory power of the indicated driving factor on the RSEI.

• Interaction detector

  The module can identify the interactive effect of two driving factors on RSEI (Zhu et al. 2020). Assuming that $q(x_1)$ and $q(x_2)$ are the explanatory powers of the influence factors $x_1$ and $x_2,$ respectively, $q(x_1\mathrm{ }\cap \mathrm{ }{x}_2)$ is the explanatory power of the two factors when they interact with each other, and five patterns of influence exist (as presented in Table 2).

Table 2. Types of interaction between two covariates.

Criterion	Interactive Forms
$q(x_1\cap x_2)<\mathrm{Min}(q(x_1),\mathrm{ }q(x_2))$	Weakened, nonlinear
$\mathrm{Min}(q(x_1),\mathrm{ }q(x_2))<q(x_1\cap x_2)<\mathrm{Max}(q(x_1),\mathrm{ }q(x_2))$	Weakened, single factor nonlinear
$q(x_1\cap x_2)>\mathrm{Max}(q(x_1),\mathrm{ }q(x_2))$	Enhanced, double factors
$q(x_1\cap x_2)=q(x_1)+q(x_2)$	Independent
$q(x_1\cap x_2)>q(x_1)+q(x_2)$	Enhanced, nonlinear

4. Results

4.1. Ecological quality assessment before and after the 2015 Nepal earthquake

4.1.1. Principal component analysis of the ecological environment index

Based on the aforementioned formula principles, four indicators that were used to construct RSEI in remote sensing image were extracted, and their corresponding mean values are shown in Table 3. There are differences in the trend of changes in the four indicators from T0 to T4. Among them, the changes in NDVI indicate that vegetation continued to decline after the earthquake, until the NDVI decreased by 8.18% at T2, and vegetation began to recover and grow after T2. The average WET representing humidity increased by 8.13% at T1 after the earthquake and showed a downward trend during the recovery period. NDBSI represents dryness, and there is no obvious pattern of dryness fluctuation from T0 to T4. The average LST increased slightly at T1 after the earthquake, with an increase rate of 0.7761%, but gradually decreased during the recovery period. Comparing the changes in the four indicators before the earthquake (T0) and in the third year after recovery (T4), it was found that the mean NDVI increased by 7.4599%, indicating a positive change in vegetation. On the other hand, the mean WET decreased by 10.7575%, implying a decrease in humidity. Furthermore, NDBSI and LST, which are negatively associated with the ecological environment, decreased by 2.5150% and 6.6774%, respectively.

Table 3. Statistics of four indicators from T0 to T1.

Indicators	NDVI	WET	NDBSI	LST
T0	0.6917	0.2891	0.6004	0.6829
T1	0.6832	0.3126	0.5323	0.6882
T2	0.6351	0.2993	0.6693	0.7107
T3	0.6880	0.2699	0.3829	0.6693
T4	0.7433	0.2580	0.5853	0.6373
Rate T0-T4(%)	+7.4599	−10.7575	−2.5150	−6.6774

Table 4 shows the results of the principal component analysis. The principal component 1 (PC1) of each year exhibits a consistent sign pattern among the variables NDVI, WET, NDBSI, and LST. As shown in Table 4, the loadings of NDVI and WET in PC1 exhibit positive values, indicating their positive contribution to the ecology. Conversely, the loadings of LST and NDBSI have negative values, suggesting their negative impact on the ecology, as expected (Xu and Deng 2022). In contrast, PC2, PC3, and PC4 do not demonstrate a regular pattern among the symbols of the four indicators. Moreover, PC1 has a contribution rate exceeding 70% for each year, implying that it captures the majority of the information. Consequently, utilizing PC1 to construct the RSEI is deemed more appropriate for this study.

Table 4. PCA results of RSEI from T0 to T1.

Year	Indicators	PC1	PC2	PC3	PC4
T0	NDVI	0.9682	0.2114	0.1011	0.0872
	WET	0.0594	0.0262	−0.9318	0.3572
	LST	−0.2074	−0.9766	0.0075	0.0568
	NDBSI	−0.1265	0.0298	0.3486	0.9282
	Percent eigenvalue	79.72%	14.60%	5.01%	0.67%
T1	NDVI	0.9588	0.2686	−0.0911	−0.0175
	WET	0.0778	0.0555	0.9937	−0.0580
	LST	−0.2724	−0.9616	0.0323	−0.0015
	NDBSI	−0.0217	−0.0065	−0.0562	−0.9982
	Percent eigenvalue	80.97%	15.36%	3.64%	0.02%
T2	NDVI	0.9362	0.3364	0.0804	0.0620
	WET	0.0361	0.0437	−0.9136	0.4026
	NDBSI	−0.3350	−0.9406	−0.0076	0.0547
	LST	−0.0997	0.0143	0.3985	0.9116
	Percent eigenvalue	83.78%	11.14%	4.23%	0.85%
T3	NDVI	0.9015	0.4095	0.0980	0.0999
	WET	0.0411	0.0286	−0.9124	0.4063
	NDBSI	−0.4115	0.9112	−0.0003	0.0218
	LST	−0.1275	−0.0360	0.3975	0.9080
	Percent eigenvalue	78.90%	16.02%	4.54%	0.53%
T4	NDVI	0.8395	0.4389	−0.2808	−0.1543
	WET	0.0200	0.0405	0.5728	−0.8185
	LST	−0.4284	−0.8965	−0.0736	−0.0854
	NDBSI	−0.3337	−0.0445	−0.7666	−0.0854
	Percent eigenvalue	70.09%	20.93%	7.73%	1.25%

4.1.2. The spatiotemporal variation characteristics of RSEI

Figure 5 shows the spatial distribution of RSEI grades from T0 to T4. The RSEI grades are indicators of ecological environment quality, with green indicating high quality, dark yellow representing medium quality, and light colors indicating poor quality. As can be seen from Figure 5 the distribution of RSEI in Nepal is high in the south, low in the north region, and transitional in the middle region. Notably, in T1 (shown in Figure 5b), there is a noticeable decrease in green colors and an increase in light colors near the epicenter. This indicates that the quality of the ecological environment in this area is disturbed by the earthquake. To enhance the visual representation of the changes in the ecological environment before and after the earthquake, we calculated the annual mean RSEI and performed area statistics based on different grades of RSEI. These findings are presented in Figure 6. The average RSEI values from T0 to T1 were 0.6932, 0.6707, 0.6480, 0.6862, and 0.7232, respectively. The results demonstrate that the overall ecological environment in Nepal showed a pattern of decreasing from T0 to T2 and then increasing from T2 to T4 over time. The decline of the RSEI value from T0 to T1 indicates that the ecological environment in Nepal has been affected by the earthquake and deteriorated. The intensity of the earthquake was high, and after the earthquake, there were many aftershocks and secondary disasters, which increased the ecological vulnerability of the area involved. In addition, the heavy rain at T2 in Nepal caused large-scale floods and debris flow disasters, which caused the average RSEI to continue to decline. It was only in the second year after the earthquake that a rebound in the ecological environment was observed. Through the calculation of the area occupied by each level of RSEI, it was observed that from T0 to T4, the Excellent grade consistently accounted for over 40% of the total area, maintaining the highest proportion each year. Conversely, the Poor grade consistently had the lowest area, never exceeding 5% of the total area. However, after the occurrence of the earthquake in T1, there was a decrease in the proportion of areas in the Excellent grade from 48% to 40%, accompanied by an increase in the proportion of areas in the Poor grade. One year after the earthquake, there was a simultaneous increase in the proportions of areas in the Poor and Fair grades, exceeding a total of 30%. At T3 and T4, the proportion of areas in the Excellent grade increased to 46% and 51%, respectively, while the proportion of areas with moderate or lower ecological environment levels decreased. Comparing the numerical values of the overall regional statistical data, it becomes evident that the ecological quality during the post-earthquake recovery period consistently improves, eventually reaching the level observed before the earthquake. These findings offer valuable insights into the ecological environment quality in Nepal from T0 to T4, as well as the changes observed during the post-earthquake recovery phase.

Figure 5. Spatial distribution of RSEI in Nepal from T0 to T4: (a) RSEI in T0; (b) RSEI in T1; (c) RSEI in T2; (d) RSEI in T3; (e) RSEI in T4.

Figure 6. Areas of different RSEI quality levels in Nepal from T0 to T4.

4.1.3. Analysis of RSEI transfer

The transfer matrix provides a quantitative representation of the relationship between different grades of the RSEI in different periods, offering a comprehensive assessment of the impact of earthquakes on Nepal’s ecological environment and the changes during the post-earthquake recovery phase. Table 5 and Figure 7 present the statistical results of RSEI level transfers from T0 to T4 in Nepal. The findings indicate that there was mutual conversion between different RSEI levels, with overall frequent conversions observed. From the data displayed, it is evident that the area where the ecological environment changed from Excellent to Good was 19,556.046 km² from T0 to T1, representing the highest value for this type of transformation. Moreover, the amount of transformation area for this type decreased gradually in subsequent years. Additionally, the results reveal a decreasing trend in the area transitioning to the Poor ecological grade from T0 to T4, reaching a minimum transfer amount of 0 km². This decrease can be attributed to the significant damage caused by the earthquake to areas with initially Excellent ecological conditions. During the recovery period, the transfer of areas with Excellent ecological environments gradually decreased, and the overall ecological environment started to stabilize.

Figure 7. The trajectory map of RSEI level changes in Nepal from T0 to T4.

Table 5. Statistics of RSEI level transfer types and area in Nepal from T0 to T4.

Year	Grade	Poor	Fair	Moderate	Good	Excellent
T0-T1 (Km^2)	Poor	2798.381	801.435	55.439	9.032	1.102
	Fair	1880.013	11051.304	1786.844	1084.996	820.390
	Moderate	284.986	3193.101	5987.430	2878.750	2857.839
	Good	126.953	1278.465	5180.334	18317.013	9322.317
	Excellent	19.109	822.926	2467.115	19556.046	41448.9836
T1-T2 (Km^2)	Poor	3730.304	1537.007	103.579	26.881	1.102
	Fair	2452.437	11644.627	2672.875	657.724	233.434
	Moderate	567.614	3639.767	6889.692	3489.072	1278.811
	Good	158.443	2707.334	5290.0389	17534.810	16147.121
	Excellent	17.196	1575.361	4848.462	11337.691	36683.673
T2-T3 (Km^2)	Poor	3346.598	2767.704	498.532	161.308	12.351
	Fair	1022.366	10568.248	4429.551	3196.929	2257.748
	Moderate	36.406	1845.134	7299.145	5149.307	5592.793
	Good	3.075	732.848	3348.092	16263.507	12737.857
	Excellent	0.000	260.339	1510.755	11271.968	41286.446
T3-T4 (Km^2)	Poor	2200.2131	1790.882	299.152	141.305	2.422
	Fair	1265.263	9366.348	2227.185	2927.169	651.415
	Moderate	13.771	1718.321	7106.509	5806.512	2872.315
	Good	0.000	221.178	1985.706	20244.260	14437.704
	Excellent	0.000	132.150	1070.998	8233.261	53730.754

aa

4.2. Characteristics of ecological environment recovery after the 2015 Nepal earthquake

To conduct a detailed analysis of the changes and recovery of the ecological environment after the Nepal earthquake, we performed the grid analysis specifically focusing on the RSEI in the severe-affected areas. Figures 8 and 9 present the results of post-earthquake damage and recovery at the administrative district level in Nepal. The degree of damage was derived from the analysis of RSEI changes at T0 and T1, and it is classified using the natural break point method. In Figure 8, larger values indicate a smaller degree of damage, while smaller values indicate a greater degree of damage. On the other hand, the degree of recovery is generated based on the RSEI values at T0 and T2, T3, and T4, respectively. Values closer to 1 indicate a higher post-earthquake recovery level, as they closely approach the RSEI value at T0. Conversely, smaller values in Figure 9 indicate a lower recovery level, signifying that the ecological environment has not yet fully recovered from the impact of the earthquake.

Figure 8. Map of the degree of post-earthquake damage in Nepal from T0 to T1.

Figure 9. Map of post-earthquake recovery degree in Nepal. (a) T2; (b) T3; (c) T4.

It is evident that the ecological environment of Gorkha District, where the earthquake’s epicenter was located, experienced severe damage. Additionally, districts surrounding Gorkha were also affected to varying degrees, with most of the damage falling into grades IV and V. As we move outward from the epicenter, the degree of damage tends to decrease. In general, the southwest and northeast of the development region suffered relatively minor damage, while the central region near the epicenter faced more serious damage. Figure 9(a) and (b) illustrates variations in the level of recovery across different districts. Notably, Kaski District and Manang District experienced the same degree of damage, yet Kaski District demonstrated a higher recovery level during the same recovery period. This observation highlights that Kaski District has a faster recovery speed compared to Manang District. In fact, the degree and speed of recovery not only reflect the changes in the ecological environment of a region but also indirectly indicate the region’s ability to withstand hazards. This highlights the importance of decision-makers considering the specific conditions of each region when formulating ecological and environmental protection measures. When looking at the two development regions as a whole, it becomes apparent that areas with less damage return a higher degree of recovery. Referring to Figure 9c, it is evident that after three years of recovery, both the Central and West development regions have achieved a relatively considerable recovery.

5. Discussion

5.1. Spatial autocorrelation analysis

5.1.1. Global spatial autocorrelation

We established a 2 km × 2km grid within the study area to investigate the potential spatial correlation of the RSEI in Nepal. Subsequently, Moran’s I index and LISA were utilized to analyze the spatial autocorrelation of RSEI. The results of the analysis, including Moran’s I index and scatter diagram, are presented in Table 6 and Figure 10. Over the period from T0 to T4, the p-value of 0.0010 confirms the scientific validity of the spatial autocorrelation analysis. The Global Moran’s I values for T0 to T4 are 0.8742, 0.8865, 0.8824, 0.8971, and 0.9217, respectively, indicating a clustering phenomenon in the spatial distribution of RSEI during these periods. Notably, the Global Moran’s I value exhibits an increasing trend, with the spatial correlation of RSEI in Nepal being the strongest at T4. As shown in Figure 10, the annual scatter points are mainly clustered in the first and third quadrants, indicating a positive spatial correlation in Nepal’s RSEI from T0 to T4.

Figure 10. Moran’s I scatter plots of the RSEI in Nepal from T0 to T4.

Table 6. Global Moran’s I of RSEI in Nepal from T0 to T4.

Year	T0	T1	T2	T3	T4
I	0.8742	0.8865	0.8824	0.8971	0.9217
p-value	0.0010	0.0010	0.0010	0.0010	0.0010
z-score	94.3341	94.8547	95.0967	97.3857	101.3949

5.1.2. Local spatial autocorrelation

In order to visually depict the spatial distribution pattern of the RSEI in Nepal and further analyze the spatial correlation between various geographical units, a LISA cluster diagram was created using a scatter plot, as illustrated in Figure 11. The High-High clustering areas from T0 to T4 are primarily concentrated in low-altitude regions in southern Nepal. On the other hand, the Low-Low clustering areas are mainly distributed in high-altitude regions in northern Nepal. The areas with insignificant distribution are mainly found in central Nepal.

Figure 11. LISA cluster map of the RSEI in Nepal from T0 to T4. (a) T0; (b)T1; (c) T2; (d) T3; (e) T4.

5.2. Analysis of driving mechanism of ecological environment quality

Understanding the driving forces behind the recovery process is crucial for developing effective strategies to mitigate the negative impacts of future hazards. This contributes significantly to developing a more comprehensive comprehension of the intricate relationship between natural and social factors in the ecological environment’s recovery after a major hazard. Ultimately, this can be harnessed to develop more efficient disaster risk management strategies in Nepal.

5.2.1. Single-factor detection analysis

The ecological environment is characterized by complexity and spatial heterogeneity, and its distribution and changes are influenced by a multitude of factors. To investigate the driving forces behind the post-earthquake ecological environment recovery, we employed six key factors encompassing both natural and socioeconomic aspects. These factors, including elevation, slope, aspect, precipitation, population density, and nighttime lights, were selected to examine their impact on RSEI changes, using the Geodetector model. These factors were classified using the natural breakpoint method. Table 7 presents the results of the single-factor detection analysis. From Table 7, it is evident that only four factors consistently passed the significance test, and the explanatory power of natural factors is greater than that of social factors. These factors, namely population density, elevation, slope, and precipitation, were identified as having significant driving forces on the change of RSEI after the 2015 Nepal earthquake. Among these four factors, elevation exhibited the most substantial contribution rate, indicating that post-earthquake RSEI in Nepal is primarily influenced by changes in elevation. The corresponding q-values for elevation at T2, T3, and T4 were 0.5742, 0.6886, and 0.7835, respectively. Additionally, the precipitation factor also exerted a notable impact, with q-values of 0.1566, 0.2721, and 0.2734 at T2, T3, and T4, respectively. The trend suggests that the influence of elevation and precipitation on the quality of the post-earthquake ecological environment has been increasing over time. The slope factor, as the third largest driving factor, has q-values of 0.0881, 0.0778, and 0.0729 at T2, T3, and T4, respectively, with a contribution rate of less than 10%. Lastly, the impact of population density on the changes in RSEI during the recovery period is very low, with no more than 5%.

Table 7. Explanatory power of the factors to the RSEI detection for T2 to T4.

Year	T2		T3		T4
	q	P	q	p	q	p
Population density	0.0470	0.0000	0.0043	0.0355	0.0124	0.0000
Nighttime lights	0.0007	0.9902	0.0018	0.8325	0.0020	0.8108
Precipitation	0.1566	0.0000	0.2721	0.0000	0.2734	0.0000
Elevation	0.5742	0.0000	0.6886	0.0000	0.7835	0.0000
Slope	0.0881	0.0000	0.0778	0.0000	0.0729	0.0000
Aspect	0.0003	0.5326	0.0008	0.0454	0.0006	0.1212

5.2.2. Interaction detection analysis

Changes in the ecological environment result from the combined influence of various factors, including natural elements and human activities. The interaction among these driving factors plays a crucial role in shaping RSEI in Nepal. Therefore, it is of utmost importance to understand the interactive influence and the mode of action of these factors behind RSEI changes. Interaction detection helps to reveal individual interactions between two factors.

As depicted in Figure 12, the results illustrate that the interaction effect of any two factors is more significant than that of a single factor. The types of factor interaction are mainly double synergy and nonlinear synergy. These outcomes highlight that the changes in RSEI after the earthquake are not solely attributed to the individual effects of population density, elevation, slope, and precipitation, but also to the interactions between these factors. Particularly noteworthy is the interaction between elevation and precipitation, which has the most substantial impact on RSEI, accounting for 63.41% to 78.63% of the observed changes. Moreover, the interaction between elevation and slope also exerts a significant influence on RSEI, contributing to over 57% of the changes. However, the single factor of the slope has a comparatively minimal effect on RSEI variations. In summary, the interaction between elevation and precipitation has the most substantial impact, followed by the interactions involving elevation and slope, and elevation and population density.

Figure 12. Results of the interaction detection.

According to the results of the Geodetector analysis, it is suggested that the difference in recovery rates between Kaski District and Manang District can be attributed to variations in natural conditions, specifically the elevation and precipitation. These findings underscore the crucial significance of taking into account the interactions among multiple factors when studying the changes in RSEI during the recovery process. In light of these insights, the local government should enhance the environmental protection policy system and strive for a harmonious balance between ecological governance and social development. Implementing appropriate and moderate measures to develop ecological civilization will lead to the achievement of a win-win scenario, wherein the regional ecological and economic environments are both improved.

5.3. Strengths and limitations

This study employed the RSEI as an ecological environment quality assessment indicator and utilized remote sensing and GIScience methods for post-earthquake recovery monitoring. In contrast to previous studies (Jiao et al. 2014; Yang and Qi 2017), we utilized RSEI instead of a single remote sensing index to detect changes in the ecological environment. RSEI is a fast evaluation method entirely based on remote sensing, and the PCA used in its construction can avoid human subjectivity. Thus, the application of RSEI in this study to monitor the changes and recovery of the ecological environment in Nepal after hazards is scientifically and logically sound. However, the application of RSEI holds limitations in terms of regional environmental monitoring. Optical remote sensing images are often affected by non-biological factors that can negatively impact the inversion of RSEI and ecological monitoring performance. For large study areas, it is always a challenging task that obtaining the complete high-resolution Landsat image data without cloud or snow. Therefore, we used MODIS data products with lower resolution, which applied methods such as quality control and perspective field angle to improve data quality to some extent. Additionally, previous studies revealed that RSEI is sensitive to different image acquisition times, so RSEI was combined with the Continuous Change Detection and Classification (CCDC) algorithm and a time-series RSEI was plotted (Sun et al. 2022). If more accurate quantification of changes in the regional ecological environment is needed, future research should introduce high-resolution images and algorithms to improve the accuracy of RSEI inversion.

To ensure the reliability of earthquake recovery monitoring, this study focuses on the three years after the earthquake, taking into account the severe flooding in Nepal in 2019, which marks the fourth year post-earthquake. Moreover, this study conducted a more comprehensive analysis and statistical assessment of the RSEI based on the transfer matrix. We further performed grid operations and displayed administrative regions to gain valuable insights into the post-earthquake damage and recovery degree in Nepal. Finally, Moran’s I index was introduced to explore the strong spatial correlation as evidence for non-random patterns in Nepal’s RSEI, and the driving forces behind the RSEI changes after the Nepal earthquake were investigated using the Geodetector model. While we acknowledged data limitations in the selection of influencing factors, future studies should pay attention to the diversity of factor selection to more comprehensively reflect the driving forces of post-disaster recovery. This will provide a solid theoretical basis for relevant government departments and decision-makers to effectively build and develop the ecological environment. In conclusion, we focused on the three-year recovery period, in-depth analysis, and spatial autocorrelation exploration, as well as the application of the Geodetector model, which contributes valuable insights to post-earthquake recovery monitoring and offers crucial guidance for ecological environment development in Nepal.

6. Conclusion

In this study, we analyzed the changes in the quality of the ecological environment in Nepal from 2014 to 2018, with a particular focus on investigating the impact and extent of post-recovery after the 2015 Nepal earthquake. Initially, we employed the RSEI generated through the integration of GEE to assess the overall changes in Nepal’s ecological environment. Subsequently, we utilized a transfer matrix to examine the shifts observed among different RSEI grades. Moreover, we conducted a detailed analysis of the extent of damage and subsequent recovery in the severely affected areas following the earthquake. Finally, we employed spatial autocorrelation analysis and utilized GIScience to explore the driving forces responsible for the observed changes in Nepal’s RSEI. The findings of this study can be concluded as follows:

1. The analysis of mean RSEI values reveals a noticeable delay in the recovery of the regional ecological environment after the earthquake, with signs of an improvement in the ecological environment becoming evident from the second year of recovery. Notably, the mean RSEI value of Nepal surpassed the pre-earthquake mean value in the third year of recovery, indicating a positive trend in ecological recovery. Furthermore, during the study period, frequent conversions between different ecological grades were observed. However, the conversion from high to low ecological grades exhibited a declining trend, while the area occupied by high ecological grades displayed an increasing trend. This suggests an overall improvement in the ecological conditions over time, with fewer regions experiencing degradation and a growing portion of the area demonstrating higher ecological quality.

2. From the district-level perspective, the districts within the Central Development Region and Western Development Region of Nepal experienced varying degrees of damage. However, the analysis of the recovery degree shows that the ecological environment of each district consistently improved throughout the recovery period. By the third year of recovery, most areas in the severe-affected areas had better recovered their ecological conditions to or even surpassed the levels observed before the earthquake. This positive trend in the overall recovery process indicates significant progress and recovery in the severe-affected areas.

3. Based on the analysis results using GIScience, a robust positive correlation was observed in the spatial distribution of RSEI in Nepal. The high-high agglomeration primarily occurred at high altitudes, while the low-low agglomeration was mainly concentrated at lower elevations in Nepal. This suggests that areas with similar ecological conditions tend to cluster together based on altitude. Additionally, this study delved further into the driving forces behind the RSEI changes by conducting single-factor and interaction tests. The results indicated that elevation, precipitation, and their interaction were identified as the primary factors influencing RSEI. These findings emphasize the significance of both elevation and precipitation in shaping the ecological environment’s recovery, and the interaction between these two factors plays a crucial role in influencing RSEI changes in Nepal.

4. To further enhance the achievements of ecological recovery and disaster risk management after the 2015 Nepal earthquake, relevant departments should optimize ecological protection measures on an existing basis. To begin with, implement land recovery and vegetation recovery initiatives to stabilize the soil, mitigate soil erosion, and facilitate ecosystem rehabilitation. Furthermore, it is suggested to establish and reinforce natural reserves and sanctuaries to safeguard the affected ecosystems and vulnerable species post-earthquake. Lastly, it is imperative to enhance public awareness and understanding of disaster risks, environmental preservation, and sustainable development through informative campaigns, educational programs, and training activities aimed at fostering environmental consciousness.

This study conducted a comprehensive monitoring of the ecological environment quality both before and after the 2015 Nepal earthquake, thereby generating valuable insights into the changes that transpired during the subsequent recovery phase. The obtained findings provide a profound understanding of the recovery process, thereby establishing a scientific foundation for relevant governmental departments to enhance their efforts in safeguarding and advancing ecological civilization and disaster risk management in the future. Furthermore, the outcomes of this study can serve as a theoretical framework for localizing Sustainable Development Goal (SDG) 11, which specifically targets sustainable cities and communities. By promoting the realization of SDG 11, this study contributes to the broader agenda of achieving sustainable development, consequently fostering a more resilient and ecologically balanced future for Nepal.

Acknowledgments

We sincerely thank the anonymous reviewers and associated editors for their valuable and constructive comments and suggestions that helped to improve this paper.

Disclosure statement

The authors report there are no competing interests to declare.

Data availability statement

The data that support the findings of this study are available upon reasonable request.

Additional information

Funding

This research was supported in part by the National Natural Science Foundation of China, grant numbers 42201077, 42177453 and 42101343; the Natural Science Foundation of Shandong Province, grant number ZR2021QD074; the China Postdoctoral Science Foundation, grant number 2023M732105; the Lhasa National Geophysical Observation and Research Station, grant number NORSLS22-05