An Assessment of Landscape and Land Use/Cover Change and Its Implications for Sustainable Landscape Management in the Chittagong Hill Tracts, Bangladesh

: Human‑caused environmental change has profoundly impacted resource management and land use patterns in Bangladesh’s Chittagong Hill Tracts. This study used multi‑temporal Land‑ sat images from 1998, 2008, and 2018 to analyze land use and land cover changes, particularly those associated with forest cover changes, in Bangladesh’s Chittagong Hill Tracts. Using object‑based image classification, Landsat images from 1998, 2008


Introduction
The recognition of the invaluable services provided by tropical forests to both host nations and the global community is widely acknowledged [1]. These unique ecosystems, characterized by their unparalleled array of plant and animal species, possess irreplaceable biodiversity and genetic resources [2,3]. Since the Earth Summit, forest policies in numerous countries have prioritized the objective of sustainable forest management, regardless of the extent of human interventions within forested areas [4]. Forests and woodlands are vital to keeping the environment carbon-free [5][6][7]. However, ecological systems, such as forests, experience constant changes due to natural biological processes, resulting in continued instability [4].
Bangladesh has approximately 2.52 million hectares of tropical forest, accounting for 10% of its area [8]. Rural areas house over two-thirds of Bangladesh's population, and their livelihoods are closely linked to the forests, either directly or indirectly [9]. For several decades, the national forests in Bangladesh have faced significant and rapid depletion, reaching a critical stage of concern [1]. The natural forests of Bangladesh experienced a consistent annual decline of 2.1% over 20 years until the early 1980s, which further accelerated to a rate of 2.7% between 1984 and 1990 [10]. Between 1964 and 1985, there was a decline in the growing stock within Chittagong Hill Tracts reserve forests, decreasing from 23.8 million m 3 to below 19.8 million m 3 [11]. Between 2000 and 2005, around 2000 hectares of forest cover were lost annually [12].
As a result of degradation, forest communities are motivated to utilize their traditional knowledge and practices to engage in activities such as conservation, reforestation, bushfire control, and the prevention of illegal forest exploitation and encroachment [13]. The combined efforts of these local communities have led to the regeneration of forested lands and the enhancement of biodiversity levels [14]. Several policies and practices in South and Southeast Asia have emerged centered around community-based forest management [15,16]. The global community places significant emphasis on preserving biodiversity, promoting forest health, ensuring sufficient forest productivity, and safeguarding the socio-economic functions associated with forest resources [17].
In Bangladesh, traditional forest management techniques have historically aimed to achieve economic and ecological objectives [18]. However, rapid deforestation occurred due to various socio-economic and socio-political factors [19,20], diminishing the effectiveness of traditional forest planning and management approaches. Unplanned human activities and unforeseen pressures exceeded planned conservation efforts, resulting in extensive deforestation and fragmentation of forest resources [19]. Given the country's dense population and limited land area, policymakers had to explore alternative management practices. In the late 1970s, social forestry was introduced as a successful alternative that transitioned the Forest Department's role from a custodial to a participatory model, involving local communities in forest protection, reforestation activities, and benefit-sharing arrangements [21].
During Bangladesh's current period of sovereignty, the Forest Act underwent its first amendment in 1989, which aimed to enhance forest protection by imposing stricter penalties and limiting the discretionary powers of forest officials and local magistrates [18]. Although this amendment primarily focused on strengthening traditional forest protection measures, it did not introduce the concept of social forestry until 2000, when another amendment was introduced, leading to the emergence of social forestry in Bangladesh [22]. The Forest (Amendment) Act of 2000 marked a significant milestone as it paved the way for the formulation of the groundbreaking 2004 Social Forestry Rules (SFR) by the government [18]. Bangladesh's fundamental principle of social forestry revolves around integrating local communities in reforestation activities, aiming to achieve multiple objectives encompassing ecological, economic, and social benefits [23].
One of the main requirements for global change research is to evaluate and track the condition of the earth's surface [24,25]. As the foundation for all living things and a key factor in global climate change, vegetation classification and mapping are crucial tech- nological undertakings for managing natural resources [25,26]. Land use/cover change (LUCC) is most frequently associated with logging, globalization, and agricultural expansion that alter the natural vegetation [27,28]. At both the local and global levels, LUCC causes several environmental issues, such as biodiversity loss brought on by greenhouse gas emissions [28,29], variations in land surface temperature (LST), and shifts in precipitation [30]. Urbanization's adverse environmental effects, which include population increase, extensive infrastructure development, and constantly shifting landscapes, are a worldwide issue [31].
Bangladesh, known for its high population density, faces increasing land pressure, particularly in forested areas, due to food production, urban settlements, and industrial development [1]. Approximately 60% of the total land is utilized for agriculture, which serves as a primary source of livelihood for over two-thirds of the rural population [32]. However, limited land availability per person, constrained by geographical factors and inadequate farming practices, hinders the country's food production capacity [9]. The forests and agricultural lands play a crucial role in the lives of people residing in Bangladesh's Chittagong Hill Tracts (CHT) region [33]. In the past, these forested landscapes provided various local and regional benefits, including food, energy, timber, water, and healthcare, while also contributing to national revenue generation. Nonetheless, the exploitation and degradation of forests, which started in the previous century [34] and continue to this day, have significant implications for the sustainable livelihoods of forest-dependent communities in terms of both direct and indirect resources [35].
Bangladesh is not exceptional in these environmental changes. Urbanization has also affected the local environment to certain degrees and made it more susceptible to land degradation. Thus, timely and accurate information about local spatial coverage, distributions of LULC categories, and their dynamics are prerequisites for the country's planning, socio-economic development, and sustainable land management. No evidence-based studies have been conducted in the Chittagong region to understand how implementing a community-based forestry policy affects forest cover change and regional land dynamics. This is the first study of its kind to examine the impact of community-based forestry on land use and land cover change in the Chittagong Hill Tracks and whether or not there are factors that influence this change. This study will highlight the aspects responsible for a land cover change in the Chittagong Hill Tracks. The study aimed to access the dynamic LULC change detection in the Chittagong Hill Tracks using remote sensing data and Geographic Information Systems technologies for 1998, 2008, and 2018. The study also used the cellular automata-Markov model (CA-Markov) to predict future land use changes under a simulated 2048 scenario. The main objectives of this study were: (i) to detect and identify the dynamics change in LULC from 1998 to 2018 in the Chittagong Hill Tracts, Bangladesh; (ii) to identify the driving forces behind these changes; and (iii) to predict the LULC map for the year 2048 using CA-Markov.

Study Area
The Chittagong Hill Tracts are part of the Chittagong Division, located in the southeastern hilly area of Bangladesh, and play a vital role (Figure 1) [9]. It geographically lies between 21.025 ′ N to 23.045 ′ N latitude and 91.045 ′ E to 92.050 ′ E longitude and covers an area of approximately 13,198.16 Km 2 , almost 40% of the total forest land area, and 8% of the total land area in Bangladesh [36]. The study area borders Myanmar in the southeast, the Indian state of Tripura to the north, Mizoram to the east, and the Chittagong district in the west. The Chittagong Hill Tracts are divided into three hill districts: Khagrachhari, Bandarban, and Rangamati. The area comprises forest, grassland, watershed, bare land, and agricultural land [9]. Paths and rows cover this area: 136/44, 136/45, and 135/45. In this region, almost twelve ethnic groups are the principal inhabitants, highly dependent on natural resources for their livelihood [37]. The major sources of income are agriculture, livestock, and the harvesting of forest products. The landscape and local people are unique parts of the Chittagong Hill Tracts. The traditional landscape of Chittagong Hill Tracts combines a mix of land uses: natural forests and plantation forests cover more than 70% of the land area in three districts. The primary plant species in the region are tropical wet evergreen/semi-evergreen and deciduous, classified as 'hill forests'. The large area covers mixed natural and planted forests [34]. This area is characterized by a tropical monsoon climate with three dominant seasons: the dry winter season (November to March), the premonsoon season (April to May), and the monsoon season (June to October). The rainfall usually starts to increase from March to June. Annual rainfall ranges from 2540 mm in the north to 2450 mm in the east, reaching nearly 3810 mm in the south and west. 136/44, 136/45, and 135/45. In this region, almost twelve ethnic groups are the principal inhabitants, highly dependent on natural resources for their livelihood [38]. The major sources of income are agriculture, livestock, and the harvesting of forest products. The landscape and local people are unique parts of the Chittagong Hill Tracts. The traditional landscape of Chittagong Hill Tracts combines a mix of land uses: natural forests and plantation forests cover more than 70% of the land area in three districts. The primary plant species in the region are tropical wet evergreen/semi-evergreen and deciduous, classified as 'hill forests.' The large area covers mixed natural and planted forests [35]. This area is characterized by a tropical monsoon climate with three dominant seasons: the dry winter season (November to March), the pre-monsoon season (April to May), and the monsoon season (June to October). The rainfall usually starts to increase from March to June. Annual rainfall ranges from 2540 mm in the north to 2450 mm in the east, reaching nearly 3810 mm in the south and west.

Remote Sensing Data
The data used in this study were obtained from various Landsat sensors, including Thematic Mapper (TM), Enhanced Thematic Mapper Plus (ETM+), and Landsat 8 OLI for 1998, 2008, and 2018 (Table 1). A total of 3 senses were acquired between November and December. All images were obtained from the United States Geological Survey (USGShttps://earthexplorer.usgs.gov/ (accessed on 12 January 2023)). All raw satellite images had a cloud cover of less than 10%. Each image cloud was removed by supervised classification extraction by the mask tool in ENVI 5.1. For the Thematic Mapper (TM), the

Remote Sensing Data
The data used in this study were obtained from various Landsat sensors, including Thematic Mapper (TM), Enhanced Thematic Mapper Plus (ETM+), and Landsat 8 OLI for 1998, 2008, and 2018 (Table 1). A total of 3 senses were acquired between November and December. All images were obtained from the United States Geological Survey (USGS-https://earthexplorer.usgs.gov/ (accessed on 12 January 2023)). All raw satellite images had a cloud cover of less than 10%. Each image cloud was removed by supervised classification extraction by the mask tool in ENVI 5.1. For the Thematic Mapper (TM), the Enhanced Thematic Mapper Plus (ETM+) image band combination was red-green-blue (RGB) [4-3-2], and the Landsat 8 OLI image band combination was red-green-blue (RGB) [5-3-2].

Satellite Image Pre-Processing
Geometric rectification is essential for producing corrected LULC maps ( Figure 2). For change detection, various pre-processing requirements, for example, geometric correction, radiometric collaboration, and atmospheric corrections, are the most important to avoid errors in the results. In this study, Landsat Collection Level-1 downloaded images were rectified and corrected-geometrically and topographically. Moreover, one of the more significant preconditions for remote sensing data analysis is an atmospheric correction, which has been done by applying the Model FLAASH (Fast Line-of-Sight Atmospheric Analysis of Spectral Hypercubes) to improve image information by transforming radiance as a sensor into surface reflectance values. Mosaic was applied for the paths and rows: 136/44, 136/45, and 135/45. The image pre-processing was done using ENVI 5.1 and ArcGIS 10.5 software.

Image Classification
Image classification was performed using a supervised classification method based on the Maximum Likelihood Classifier of given classes ( Table 2). Training samples were selected with sufficient homogeneity to symbolize the spectral and spatial characteristics of each LULC class and maximize classification accuracy. This step is the most significant element of supervised classification because spectral signatures extracted from training samples will define the overall accuracy of the classification and, therefore, the final LULC map. Table 2. Classification scheme/categories.

Code LULC Classes Description
1 Water Area covered with ponds, lakes, rivers, and seasonal streams. 2 Bare land Areas described with non-vegetative cover, rocks, and some limited settlement areas. 3 Forest Areas with natural or artificial woody vegetation and canopy cover more than <10%. 4 Grassland Areas covered by grasses, agriculture, and some sparse shrubs in between.

Accuracy Assessment
Accuracy assessment is a post-classification step accomplished through the correspondence evaluation of classified LULC maps against assumed factual geographical reference data [38]. The image classification accuracy was assessed with ground truth data in the form of reference data points obtained from Google Earth images of 1998, 2008, and 2018, used to obtain overall accuracy, user accuracy, producer accuracy, and Kappa coefficient. The classifier, carried out using 58, 90, and 211 points, was generated randomly for 1998, 2008, and 2018 supervised images, respectively. Each point has a specific pixel value based on ground-truth data and visual interpretation. Identified points were considered classified values, and pixel values were considered reference values.
Classification results were carried out on thematic maps using different truth reference data, and accuracy assessments were presented for all classifications. Overall accuracy, user accuracy, producer's accuracy, and Kappa coefficient accuracy derived from the error matrices were used for the accuracy assessment of the final maps. An error matrix and Kappa statistics were generated from reference and classified data.

Image Classification
Image classification was performed using a supervised classification method based on the Maximum Likelihood Classifier of given classes ( Table 2). Training samples were selected with sufficient homogeneity to symbolize the spectral and spatial characteristics of each LULC class and maximize classification accuracy. This step is the most significant element of supervised classification because spectral signatures extracted from training samples will define the overall accuracy of the classification and, therefore, the final LULC map.  In this study, satellite images were accurately assessed with the Landsat images Landsat 5-TM 1988, Landsat 7-ETM+ 2008, and Landsat 8 OLI 2018, which the ground truth data likely equates to. Overall accuracy was calculated by dividing the sum of correctly classified pixels (diagonal) [39] by the total number of reference pixels.
Overall accuracy = (∑_(e = 1)^Lncc)/N The overall classification accuracy was calculated for all classifications, as well as classspecific accuracy, which can be created on the user and producer levels. Users' accuracy, i.e., the error of commission [40] and user accuracy assessment, was calculated by dividing Users' accuracy = ncc/(n + c) (2) Producer accuracy is the measure of the error of omission [40], and producer accuracy was calculated by dividing the number of correctly classified pixels in each category by the total number of reference pixels in that category (total column).
Producer accuracy = ncc/(nc+) The Kappa coefficient was calculated by using the following equations: where r = the number of rows in the error matrix; P_ii = the proportion of pixels in row 'r' and column 'i'; P_(i+) = the proportion of the marginal total of row 'i'; P_(+i) = the proportion of the marginal total of column 'i'.

Cellular Automata Model
The cellular automata model mainly consists of cell, time, and rule neighbor. The filter of the cellular model determines the neighbors. Cellular automata (CA) are a spatiotemporal calculation of a dynamic process model used for LULC change. The commonly used neighborhoods are Moore, the extended Moore, and von Neumann. The expression of the Cellular Automata equation is: where S is the set of states of the finite cells, the (t, t + 1) are different times, N is the neighborhood of cells, and f is the transformation rule of local space [41].

CA-Markov Model
The CA-Markov model is valuable for modeling LULC changes and can simulate prediction changes. The spatial prediction accuracy can be effectively simulated simultaneously [42]. CA-Markov model effectively combines the Markov model and the CA model. This approach is based on a Markov stochastic probability matrix for predicting the transition from one condition to another [30]. The Markov chain model is most commonly used to simulate transitions, parameters, and trends. It created probability transition matrices to anticipate and categorize probable land use/cover change (LUCC) and urban development scenarios and investigate land change simulation trends [30]. Equations ((6)-(8)) based on the conditional probability formula were used to estimate trends: where S(t) is the state of the system at time t, S (t + 1) is the state of the system at time (t + 1), and Pij is the matrix of the transition probability in a state. By projecting 2018 and 2048, the cellular automata (CA) and Markov chain models are employed to generate the LUCC

Landscape Metrics
The FRAGSTATES software version 4.2.1 [43] was applied for computed landscape metrics in this study area. Landscape metrics characterized the spatial fragmentation and heterogeneity for 1998, 2008, and 2018, covering three decades. The selected landscape patterns were divided into three levels, i.e., patch, class, and landscape. FRAGSTATS is a spatial pattern analysis program implemented by the decision maker, forest manager, and ecologists to analyze landscape fragmentation or describe characteristics of the landscape and elements of those landscapes [44]. FRAGSTATS version 4.

LULC Change-Driving Forces Model
Previous studies have suggested that the main LULC of the Chittagong Hill Tracts area driving factor is climate change, while anthropogenic contributed to LULC slightly. Therefore, this study only considered climate change as a leading driving factor.
The monthly meteorological data from 1998 to 2018 were collected from the meteorological station of Bangladesh. Each year of the CHT area was obtained by substituting time into the linear regression model. The annual average temperature (AAT), annual average relative humidity (AARH), and annual average rainfall (AAR) were used in the analysis. Climate change factors such as the average yearly temperature and annual average relative humidity were used to establish corresponding data sets. All the meteorological data processing steps were carried out in linear regression.
Climate change factors such as the annual average temperature and annual average relative humidity were used to establish corresponding data sets. The relationship between dependent variables and independent variables is given below: where Y is a dependent variable, a is a constant turn, x 1 , x 2 , and x 3 are independent variables, and b 1 , b 2 , and b 3 are the coefficients of independent variables. The constant coefficient illustrates the value of Y in all dependent variables at a zero time zone. In addition, the parameter coefficients state a change in Y for one unit increase in the dependent variable [45]. T-tests were applied to verify the statistical regression coefficients' significance and constant in the linear regression model. The determination sample coefficient (R 2 ) was applied to explain the contribution to the dependent variable in the linear regression model. R 2 means how closely dependent variables are related to independent variables [46].

LULC Classes and Their Distribution
Land use cover changes (LULCCs) were computed for 1998, 2008, and 2018, focusing on vegetation, grassland, water bodies, and the bare land area in the study area ( Figure 3). In 1998, the maximum land area covered by grassland accounted for 44.71% (5901.11 Km 2 ),  Figure 4).

LULC Change Detection
The post-classification comparison of LULC change for each class within the study period from 1998 to 2018 is shown in Figure 5, whereas the result analysis of LULC change is indicated in Table 5. From 1998 to 2008, only bare land increased by 665.04 km 2 . In contrast, water, forest, and grassland decreased significantly by −30.40 Km 2 , −528.04 Km 2 , and −106.60 Km 2 , respectively, with a significant decrease in bare land and grassland, i.e., −692.07 Km 2 and −2614.17 Km 2 , between 2008 and 2018 (Table 5). During the three decades, the grassland decreased by −2720.77 km 2 , while the forest increased by an area of about 2708.36 km 2 .
Regarding LULC conversions and transformations, the most dynamic change in the study area is the conversion of grassland to forest by 13.34% (Table 6), followed by forest into grassland by nearly 58.03%. The transformation of bare land to forest was 1.68%, and the transformation of bare land into grassland was calculated to be 6.28%. However, forest land conversion to bare land in the region was estimated to be 12.10%, while the transformation of grassland into bare land was observed to be 29.50%. Therefore, between 1998 and 2018, the region experienced many changes in its LULC.

LULC Change Detection
The post-classification comparison of LULC change for each class within the study period from 1998 to 2018 is shown in Figure 5, whereas the result analysis of LULC change is indicated in Regarding LULC conversions and transformations, the most dynamic change in the study area is the conversion of grassland to forest by 13.34% (Table 6), followed by forest into grassland by nearly 58.03%. The transformation of bare land to forest was 1.68%, and the transformation of bare land into grassland was calculated to be 6.28%. However, forest land conversion to bare land in the region was estimated to be 12.10%, while the transformation of grassland into bare land was observed to be 29.50%. Therefore, between 1998 and 2018, the region experienced many changes in its LULC.

Prediction of LULC Change Based on the Markov Model
The state transition area map was created according to LULC maps from 1998 to 2018, which can be used to predict, using the CA-Markov model in IDIRISI software version 17.0, the land requirements for the different LULC types in 2048. The predictive results map for 2018 is obtained with a 5 × 5 contiguity filter, whose running cycle is 30 years. The combination of cellular automata (CA) and the stochastic transition matrix of the Markov chain model resulted in LUCC for the projected period of 2048 ( Figure 6). Map accuracy for the projected land use/cover change for predictive years was classified by a sufficient Kappa coefficient value of 0.97. A 2048 map predicted that the maximum area covered by forest accounts for 9129 Km 2 (69.17%) and the minimum area covered by water accounts for 665 Km 2 (5.05%). A decrease of 1.51% of bare land and 5.61% of grassland cover areas were estimated during 2018-2048, respectively (Table 7). However, forest and water area cover will expand by 6.40% and 0.72%, respectively, during 2018-2048.

Landscape Metrics Analysis of Land Use and Land Cover Structure
From LULC maps from 1998 to 2018, the most changing classes, such as water, bare land, forest, and grassland, were chosen to analyze spatial landscape metrics at class and landscape levels. The most significant change among land use and land cover classes is increased forest and decreased grassland. The statistical results of the landscape metrics in the Chittagong Hill Tracts area are shown below (Table 8).

Landscape Metrics Analysis of Land Use and Land Cover Structure
From LULC maps from 1998 to 2018, the most changing classes, such as water, bare land, forest, and grassland, were chosen to analyze spatial landscape metrics at class and landscape levels. The most significant change among land use and land cover classes is increased forest and decreased grassland. The statistical results of the landscape metrics in the Chittagong Hill Tracts area are shown below (Table 8).  (Table 9). There is no significant change in the heterogeneity of the landscape.

Analysis of the Driving Forces behind LULC Change
The linear regression model result indicated that all driving factors had significant values in adjusted R square for the model of change in the forest, bare land, grassland, and water of 0.619, 0.559, 0.718, and 0.752, respectively (Table 10). The study's results explained that AT was highly significant and AR substantially impacted forests (p < 0.01, p < 0.05; Tables 8 and 11). AR was highly influential, and ARH significantly impacted bare land (p < 0.01, p < 0.05; Tables 8 and 11). AT had a highly significant impact on grassland (p < 0.01; Tables 8 and 11). AT had a highly substantial effect, and ARH, AR had a considerable effect on water (p < 0.01, p < 0.05; Tables 8 and 11) and the change of CHT in Bangladesh.  Note: Y1 = forest land, Y2 = bare land, Y3 = grassland, and Y4 = water. AR = average rainfall, AT = average temperature, and ARH = average relative humidity. ** = highly significant; * = significant; NA = non-significant.

Accuracy Assessment
After calibration, land cover data were extracted from photographs and corrected for atmospheric effects. Landsat 5 TM 1998 had six spectral bands plotted upon them, but Landsat 7 ETM+ 2008 and Landsat 8 OLI 2018 had seven bands, respectively. In three different periods, a summary of supervised classification accuracy was recorded as the highest in 1998 (89.65%) and 2008 (84.44%). This result is consistent with the study findings by the authors [44]. The study indicated that all classes of producers were classified over 85% accurately. Similarly, a study by Dewan and Yamaguchi [47] using MSS found that the lowest overall accuracy was 85.6%. Another study by Kayiranga et al. [48] revealed that from 1986 to 2015, all images' overall accuracy and Kappa coefficient results were greater than 75%. In this study, among the three periods of 1998, 2008, and 2018, Kappa coefficients were above 0.8% in all classes except 2018. Very similar results were found in the study conducted by del Castillo [49], who investigated that all forest classes' Kappa coefficients were above 0.8 except Q. pyrenaica in 2010.

Trends in Forest Cover Changes
According to estimates of forest cover change (loss and gain) derived from Landsat images during the past twenty years, forest and tree cover increased from 1998 to 2018. The national assessment led by FAO [50] showed that the net growth of tree cover has remained stable since 1990. Different spatial and temporal extents may, in part, explain these discrepancies. Furthermore, while our results show a net gain in forest cover over time, details on the visibility and psychological impact [51] on forest loss vs. gain may also shape perspectives on forest dynamics [52].
The most notable LULC changes were the transformation of grassland to forest. The analysis of the Landsat images reveals some trends that complement the LULC from 1998 to 2018. From 1998 to 2018, the change in LULC classes shows that the Chittagong Hill Tracts region has been impacted and degraded with forest (2708.36 Km 2 ) and bare land (−27.03 Km 2 ), respectively. A similar result was contrary to the range reported elsewhere [53,54]. A similar result was proposed by Baral et al. [55], who reported that forest cover increased with a similar loss in the extent of all other land cover types. The human population was the most important factor in increasing or decreasing land use and land cover change.
Here is the evidence that variations in forest cover changes may be associated with the types of management regimes [56]. Community-based management or secure tenure can result in positive social interactions, better forest conditions, and opportunities to increase tree cover and distribute economic benefits [57]. State-led plantation initiatives and the National Forest Policy of 1994 have been primarily responsible for expanding forests in the CHT region, focusing on plantations made possible by the relocation of native peoples [58]. Bangladesh has demonstrated encouraging progress in curbing the annual deforestation rate. It has reduced deforestation from 2.1% from 1960 to 1980 [59] to approximately 0.2% between 1990 and 2010 [18].

Future LULC Change Prediction
In contrast, between 2018 and 2048, the predicted LULC in Chittagong Hill Tracts grassland significantly decreased after 1998, the forest increased, and the remaining area increased by a certain amount in 2048. Similar predictive results were proposed by Liping et al. [41]; they reported that water and bare land increased and the woodland area decreased significantly. Predictive results were revealed by [60]: bare land will increase and potentially lead to the loss of ecosystem services. The Markov predictive results have shown that the probability of each class changing in LULC in the future is high.

Landscape Analysis
During the study period  in Chittagong Hill Tracts, changes occurred in four major land classes: water, bare land, forest, and grassland. The results of the landscape indices (forest) revealed that the values of NP decreased proportionally with increased CA, ED, IJI, PSCOV, and MPS. Similar results were found in the study conducted by Kayiranga et al. [46], who reported that during the whole study period , the values of MPS increased with a decrease of NP inside the Nyungwe-Kibira Park. del Castillo et al. [49] also reported that the forest inside the Moncayo Natural Park (Spain) was slightly fragmented. This study claims landscape pattern metrics at heterogeneity with a similar identic landscape configuration. According to research correlations on landscape metrics, Shannon's Evenness Index (SEI) and Shannon's Diversity Index (SDI) were also selected as indices for the landscape pattern analysis [61]. During the study period, Shannon's Diversity Index and Shannon's Evenness Index showed a slightly decreased landscape. The landscape metrics analysis provided valuable information about land use and land cover change, especially fragmentation and heterogeneity weighed at the landscape class and landscape level. The landscape structure became more fragmented and heterogeneous [62]. Moreover, during the study period, landscape fragmentation increased and forests became a major land area [63,64]. Landscape (fragmentation and heterogeneity) and land use/landcover results indicated that forest management significantly increases, decreases, and maintains landscape patterns [49].

LULC Change-Climate Factors
The linear regression results assumed that the climate change factors model variables indicated that three impact factors were selected for the forest: bare land, grassland, and water. These three factors are less than significantly related to forest, bare land, grassland, and water. According to Salman-Mahini and Kamyab [65], there was a linear relationship between dependent and independent variables in multiple linear regression. The result of linear regression indicated that the climate change factors (average rainfall, average humidity, and average temperature) model as independent variables estimated the influence of drivers behind changes in four kinds of LULC change: forest, bare land, grassland, and water. Grassland was most affected by rainfall, the forest was significantly affected by rainfall and temperature, bare land was significantly affected by rainfall and humidity, and water was significantly affected by rainfall, humidity, and temperature. Rainfall was the main factor driving LULC changes. A similar result was found in a Northern Iranian study by Jahanifar et al. [46]. The linear alternative of income per capita, rain, and temperature with a definite coefficient of 0.4 as independent variables qualified for estimating forest area reduction.
Although research into the spatial-temporal trends and connections between forest cover and local land use practices will aid in identifying regions where changes have taken place and making predictions about future changes, this research alone will not be sufficient to explain the neighborhood's driving factors of change, which are necessary for making trained, evidence-based decisions. Understanding the mechanisms of forest cover changes, like deforestation and forest regeneration, in the Chittagong Hill Tracts at various geographical scales requires knowledge of the local driving variables at work. Therefore, future research should investigate the causes and consequences of land use and forest cover changes within this landscape's socioeconomic, policy, institutional, and past histories and how these shifts relate to larger shifts happening at the regional and global levels.

Conclusions
Knowledge of past land use trends and land cover change is vital to comprehending the connection between landscape dynamics and ecological responses. Our method, which integrates multi-temporal remote sensing data with GIS techniques, has enabled us to quantify and characterize the spatial-temporal pattern of LULC changes, notably those associated with forest cover changes in mountainous regions. The findings revealed that the terrain of the Chittagong Hill Tracts underwent LULC shifts between 1998 and 2018. Between those years, forest cover and water body area expanded while grassland and bare land shrank. The spatial pattern shift demonstrated that between 1998 and 2008 and between 2008 and 2018, the forest cover gained and lost area with varying annual intensities and dynamics.
The annual increase was significant during both periods, while deforestation was prominent during the first period (1998-2008) but inactive during the second period (2008-2018). Between 1998 and 2018, there was a net change of 2708.36 Km 2 in forest cover and a total change of 8284.62 Km 2 . The LULC structure, predicted based on the CA-Markov model that in 2048, the forest area would increase drastically with a consistent decrease in grassland and bare land areas. The linear regression results for the climate change factors model variables indicated that grassland was most affected by rainfall, the forest was significantly affected by rainfall and temperature, bare land was significantly affected by rainfall and humidity, and water was significantly affected by rainfall, humidity, and temperature. Rainfall was the main factor driving LULC changes. Currently, the area faces different environmental issues that threaten these resources, such as climate change, LULC change, disturbance of species diversity, ecosystem fragmentation, and flooding. It was found that massive population growth in settlements is an important factor influencing LULC changes and the implications of sustainable landscape management. Furthermore, the findings of this study have some important suggestions that a scientific forest management system should have a strict process. Additionally, efficient programs should be run to educate local communities about sustainable land management. Data Availability Statement: Data will be available on demand.