{'fraction_non_alphanumeric': 0.027205568972865465, 'fraction_numerical': 0.020812615428327887, 'mean_word_length': 4.9310490099705095, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 3, 'https://': 2, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Key indicators describing the evolution of landslides in the Zhuoshui River Basin caused by the Chi-Chi earthquake in Taiwan

Abstract This study selected the landslide caused by the 1999 Chi-Chi earthquake, which occurred once in a century, and analysed the spatio-temporal distribution of landslide evolution in the basin through the concavity, slope, and sediment delivery ratio (SDR) index of the landslide site. From the change in the determination coefficient of the relationship between the concavity index and the collapse rate, it can be shown that the strong earthquake and subsequent rainfall affect the basin collapse rate for about 8 years. The spatial distribution of landslide location types can be classified by the characteristics of slope and SDR into three categories, indicating that the evolution of landslide location types is gradually transforming from scattered off-bank landslides to concentrated near-bank collapses. Type 1 refers to the initial slope failure caused by the crustal uplift and Type 2 is the unstable slope caused by the earthquake and the subsequent rainfall-driven collapse. Near-bank landslides caused by the accumulated debris from Types 1 and 2 landslides are Type 3. Based on the huge and complex big data of the landslide evolution map, the focus of this research is to extract key indicators and build models to analyse the evolution mechanism of large-scale landslides caused by the earthquake. Highlights Topographic factors can be used as a criterion for judging the stream concavity of watershed. Concavity, slope, and sediment delivery ratio are key indicators, which can be used to extract sediment-prone areas for slope utilization and disaster prevention management.


Introduction
Taiwan is located at the junction of the Eurasian Plate and the Philippines Sea Plate (Huang et al. 2015). It has strong orogenic movement and the average uplift rates of the Hengchun Peninsula (5.3 ± 0.2 mm/year), Tainan area (4.3 ± 0.4 mm/year), and Coastal Range (5.0 ± 0.4 mm/year) have little difference from each other in the last 9000 years (Peng et al. 1977). The actual average uplift rate (minimum uplift rate plus denudation rate) of coastal terraces may be about 5.5 mm/year. If the 5.5 mm/ year ascent rate is also applied to the Central Range, then the height of the Central Range will not change much today, because the average denudation rate of the Central Range is also about 5.5 mm/year (Li 1976). Suppe (1981) studied the orogenic movement of the Taiwan terrane based on the spatio-temporal equivalence relationship and found that during the orogenic movement, the height of the central mountain range continues to rise until a balance is reached between the structural uplift of the orogenic movement and the material loss caused by denudation. This type of terrain is called steady state.
Large-scale earthquakes triggered by plate movement (internal forces) will uplift the crust and cause landslides or large-scale landslides (Petley 2012;Kennedy et al. 2015). A catastrophic earthquake with a Richter magnitude of 7.3 occurred at Chi-Chi of Nantou County, where the Zhuoshui River basin is located, on 21 September 1999. There were heavy casualties and extensive damage to buildings and property losses. A large number of landslides also occurred in Central Taiwan. There were more than 20,000 sites with a total area of 15,977 ha of landslides identified as a result of this quake (Lin et al. 2005) and Chang (2000) found that most landslides occurred at the outer edge or inner side of the terraces. The inventory of landslides caused by earthquakes is extremely important to expand our understanding of the relationship between earthquakes and landslides (Tanyaş et al. 2017). Marc et al. (2015) studied the evolution of the landslide rate in the epicentral area of four moderate to large earthquakes (Mw, 6.6-7.6), and pointed out that the landslide is related to the cumulative precipitation in a given time interval. Yunus et al. (2020) reconstructed the restoration of landslides after earthquakes from 2000 to 2018 based on the MODIS NDVI time series, and used the simple vegetation restoration rate (VRR) to quantify the vegetation restoration in the Wenchuan earthquake-stricken area, and the attenuation trend of landslide activity after the earthquake was discussed. These natural disasters are normal energy releases, but each large-scale disaster may take several years to return to the pre-disaster stable state (Bontemps et al. 2020). As climate and environmental changes complicate the types of disasters, it is necessary to understand the natural fluctuations of sediment-related disasters in watersheds (Odum et al. 1995).
Since the rate of plate tectonic movement is quite slow, the long-term characteristics of tectonic movement can only be reflected by the study of topographic features. Therefore, geomorphic indices have been used as a role to investigate tectonic activities and landform development. Anand and Pradhan (2019) used ASTER Global digital elevation model (DEM), which are based on stream length gradient index, hypsometric integral, asymmetry factor, basin shape, valley floor width to valley height ratio, mountain front sinuosity index, to assess active tectonics. Wang and He (2020) studied responses of stream geomorphic indices to fault activity in the Daqingshan area of China.  explored the relationship between channel profile and abrupt increase mountain uplift rate. Their work shows that channel steepness and basin elevation and relief are important factor to tectonic forcing on regional topography, and knickpoints and retreat distance, basin area, river length are correlated. In the process of structural evolution, erosion is the main force. Previous studies have pointed out that changes in the longitudinal profile of the river can reflect the changes in geomorphology (Hack 1957(Hack , 1973Snyder et al. 2000;Giaconia et al. 2012;Ayaz and Dhali 2020). Hack (1957Hack ( , 1973 suggests presenting the stream longitudinal profile in semi-logarithm form and the plot shows a straight line. However, most stream longitudinal profiles in active mountains present a convex curve rather than a straight line in the semi-logarithmic plot (Lee and Tsai 2010). The ups and downs in the base profile of the river can alter the course of the river and changes in channel geometry. The river valleys may get constricted or widened, and there can be an increase in the rate of deposition or changes in the stream energy which may cause the formation of fluvial landforms like alluvial fans, terraces, and bars. The changes in the base profiles of rivers could be in response to tectonics or climate. The relief of stream source area could be altered due to tectonics. This influences the gradient of the fan and sediment production. Local base level changes can also modify this relationship, influencing the erosional and depositional regime in the distal fan zone (Ayaz and Dhali 2020). Therefore, the concept of steady-state geomorphic indicators is utilized to quantify the characteristics of geomorphic structures and to predict the geomorphic evolution.
Previous studies of bedrock channel gradients (Hack 1957(Hack , 1973Howard 1994;Howard et al. 1994;Whipple and Tucker 1999;Whipple et al. 2000) show that the local channel slope S and the contributing drainage area A has the relationship in the form: where the parameter k s and h are called steepness index and concavity index, respectively, when the rate of uplift is about the same as of undercut erosion of the channel bedrock. The plot of S to A will be a straight line on a logarithmic paper. Therefore, the relationship of S-A can be used to judge if the terrain has reached steady state. Generally, a river channel from upstream to downstream can be classified into colluvial channel, bedrock fluvial channel, and alluvial channel. The transition area from colluvial channel to bedrock fluvial channel is about 10 5 $10 6 m 2 (Tarboton et al. 1991;Duvall et al. 2004).
Considering that the ground uplift caused by the long-term plate orogenic movement is reflected in the topographic features, the advantages of automatically extracting the topographic index, and the establishment of the river incision model requires a lot of manual processing and analysis. There are two types of landslide risk assessment methods, namely physical-based modelling and data-driven modelling (Liu et al. 2021). Physically-based models are widely used in landslide prediction considering various hydrological processes, such as precipitation and geomorphic factors in the basin (Schilir o et al. 2016). The data-driven model has the ability to predict landslides by identifying the relationships among various data sets, so the in-depth understanding of the physical processes of the watershed is not high. As a data-driven model, machine learning can discover patterns directly from data without predefined rules. McCarthy (2006) proposed a plan used the term 'artificial intelligence' (AI) to describe AI as allowing machines to perform human behaviours, showing the same intelligence as human thinking, but due to the concept of excessive innovation, it only caused a sensation at the time. Its renaissance began in the 1980s, when it was replaced by a large number of statistical theories such as probability and statistics, allowing computers to learn a set of skills on their own through data, called 'machine learning'. Deep learning is a sub-discipline of machine learning, which consists of continuous operations. With the rapid development of machine learning, as an advanced data analysis tool, deep learning has received many applications and inspirations in the analysis of geological disasters (Ma and Mei 2021). If topographic factors can be used to estimate the concavity index, the relationship of terrain parameters and concavity index in the sample watersheds can be established through multivariate analysis. This relationship can later be used to calculate the spatial distribution of the concavity index of the basin at different watershed scales, which can greatly reduce the complicated process of extracting the concavity index. To effectively explain the spatial distribution of landslide evolution, key indicators that reflect the landslide in the watershed were used to explore its evolution trend, and the commonly used classification method in the deep learning model-unsupervised learning was applied to classify a large number of landslide data through topographic features.
Previous studies mostly focused solely on the exploration of phenomena such as crustal changes and orogenic movements, or only analysed changes in spatial distribution for specific disaster events. This study is based on the Zhuoshui River Basin, using the stream-power river incision model as the basis to analyse the steady-state concavity index of the tributary watersheds, and then combines the topographic factors extracted from DEM to find their correlation. In addition, these factors are tested by multivariate analysis to evaluate the weight of their influence on the concavity index. Finally, the concavity is used in conjunction with the landslide map data of the Chi-Chi earthquake and subsequent annual surveys to illustrate the evolution of landforms and landslides, and integrate the concavity and slope/location indicators that are prone to landslide to provide references for disaster prevention management.

Study area
The study area is the Zhuoshui River Basin, which spans Changhua County, Nantou County, Yunlin County, Chiayi County, and other areas. The catchment area is 3162 km 2 and the most upstream is Wushe River, which gathers the water from the west drainage divide of Hehuan Mountain and flows down the rift valley in the southwest direction, until it converges with the several streams which flow into Zhangyun Plain afterwards. The topography of the basin varies greatly. The highest elevation is 3941 m. The main slope is the sixth grade (55-100%), and the western slope direction occupies the most area. In this study, the terrain uplift in the short and medium term (1-10 years) has no substantial significance. It can be regarded as a uniform profile for this equilibrium state. There are a total of 12 main tributaries of the Zhuoshui River. They are Tarotwan, Wanda, Liqi, Kashe, Danda, Junda, Zhuogun, Chenyoulan, Shuili, Qingshuigou, Beishi, and Qingshui Stream, respectively ( Figure 1). In this study, the 1999 Chi-Chi earthquake, which occurred once in a century, was selected as the base period for landslide evolution ( Figure 2).

Research process
The Zhuoshui River and its main tributaries with accompanied watersheds are extracted from the DEM before Chi-Chi earthquake with a spatial resolution of 40 m Â 40 m produced by the Forestry Bureau Aerial Survey Office. Based on the results of extraction, the quantitative and steady-state geomorphic indices can be obtained from Equation (1) by regression analysis. In addition, the longitudinal profiles of each tributary are extracted from DEM as well. These profiles are interpreted by mathematical fitting functions. The topographic factors and collapse rate of the catchment area are analysed by multivariate statistical analysis using SPSS software to evaluate the correlation between factors and the geomorphic index h.

Procedures of establishing river incision models
Step 1: Stream system extractionThe stream system of the basin is delineated according to the algorithm of non-depression slope direction, and DEM can be used to calculate the non-depression slope information and simulate the drainage direction of the surface (Lin 2000). According to Taiwan's 1/25,000 topographic map, the surface water system is digitized and the topographic information of the main river is extracted for topographic equilibrium analysis.
Step 2: S (Slope)-A (drainage area) extractionBeginning from the source point of a stream as the start point, the next point is selected where its elevation is about 100 m lower than the previous point along the stream and so on. The slope (S) is calculated as the elevation difference divided by the distance along the river section between these two points. The corresponding drainage area (A) is dissected from upstream catchment area according to the outlet location of that river section.
Step 3: Evaluating empirical values of geomorphic indices from bedrock modelCalculate the h and k s of each watershed, which can be obtained by linear correlation in the double logarithmic relationship diagram of S-A.
Step 4: Best fitting function of river profileTake the elevation as the vertical axis and the distance from the headwater as the horizontal axis of the river profile. The coefficient of determination (R 2 ) between the fitting function and the actual river profile is used to determine the best model of the river profile. Under relatively stable conditions such as less intense orogeny and less severe climate change, if considering only on degree and shape of the depression to reflect the evolution of the river longitudinal profile, the evolution sequence should be: linear profile ¼ > exponential profile ¼ > Logarithmic profile (Chen 2004).

Topographic factors extraction
Seventeen topographic factors are selected to be analyzed to determine their relationships with geomorphologic indices. These factors with their formulae or description are listed in Table 1.

Watershed landslide analysis
The landslide maps for the years since the Chi-Chi earthquake provided by the Soil andWater Conservation Bureau (1999-2003) and Forestry Bureau (https://tgos.nat.gov.tw). The Forestry Bureau used FORMOSAT-2 high spatial and temporal resolution multispectral image (https://www.nspo.narl.org.tw/history_prog.php?c¼20021802&ln¼zh_TW) to interpret the landslides in Taiwan. FORMOSAT-2 began taking images in 2004, and its telemetry cameras have a ground resolution of 2 m for black and white images and 8 m for colour images. The satellite can take images over Taiwan every morning, and its daily real-time image acquisition function has become an important reference for Taiwan's disaster relief units when major disasters occur. The original satellite image has undergone various corrections to improve the image quality. The landslide map data are produced through the expert-aided delineation system for collapse and shadow areas, and the 25 cm aerial photo captured by the Forestry Bureau Aerial Survey Office is used as the basis for accuracy verification. The average overall accuracy of landslide interpretation can reach 98%, and a landslide catalogue of Taiwan's island-wide consistency standard has been established (Lin 2013). The landslide map from 1999 to 2003 was produced in response to the Chi-Chi 921 earthquake and major sediment disasters, and was drawn based on aerial photo data. On the whole, the landslide map data are circled or checked for accuracy by the real data of aerial photos, and it has quite high accuracy and application value.
Geographic information system (GIS) software ArcGIS was applied to analyse the set of maps for extracting the differences of landslide location in the Zhuoshui River Basin. The collapse rate in a watershed is calculated as follows: The cumulative landslide spatial distribution of each survey year after the Chi-Chi earthquake (1999-2017) was used to extract the watershed collapse rate to explore the relationship between the collapse rate and the geomorphic index h.
Frequency distributions of slope and SDR of the landslide grid in each survey year after the Chi-Chi earthquake are used to describe the location and process of sediment transport. The calculation of SDR in watersheds is based on a GIS combined with other techniques such as grid computing and automatic watershed delineation tools. The calculation methodology is an extension of the theoretical framework constructed by Lin (2002). It is hypothesized that sediments on slopes are mainly driven by surface water, transported to channels (perennial streams), and then washed away. On the basis of SDR definition with the aforementioned hypothesis, SDR can be regarded as the sediment contribution of a specified grid point to the closest downstream channel outlet grid point. If the sediment contribution is expressed by the upstream inflow area of the grid, the larger the upstream inflow area is, the greater the runoff generation is, and the more sediment can be transported to the river channel. Thus, SDR at any grid point of watershed slope can be regarded as the ratio of the upstream inflow area of the grid point to that of the grid point which is the nearest downstream channel outlet.

Multivariate analysis
Multivariate analysis is an effective tool. It can find out the regularity of the complex relationships among various factors in the environmental system, simplify complex phenomena, extract the main messages, and analyse and judge the results of the data. With the display ability of GIS in space, it can rationally and systematically sort out, judge, and predict complex problems or phenomena (Lin 2000).
Discriminant analysis was first proposed by Fisher (1936). Its concept is based on the linear combination of independent variables as the basis for grouping observations. That is, find a set of linear discriminant functions and group the observations into the best group. The general formula of the linear discriminant function is: where D is discriminant score, X 1 , X 2 , … X p are independent variables, and B 0 , B 1 , … B p are classification function coefficients. The coefficient of the discriminant function is similar to the multiple regression coefficients. Because the variables are related to each other, the importance of individual variables cannot be evaluated. The relative importance of the variables can only be compared with the size of the coefficient. Standardized canonical discriminant function coefficients and structure matrix can be used to judge the influence of variables on discriminant function.
The maximum height difference of the watershed divided by the horizontal distance between the two points Width

Establishment and analysis of bedrock incision model
Based on the stream-power river incision model, the S-A double logarithmic graph drawn directly from the river slope (S) and the upstream catchment area (A) can calculate the concavity index (h) and steepness index (k s ) of each tributary watershed (Table 2 and Figure 3)

Concavity index and topographic factors
Concavity index and topographic information of the analysed watershed are shown in Table 3.

Landslide analysis
Based on the landslides of the Chi-Chi earthquake, using the landslide area information of the pre-and post-periods chronologically, the spatial distribution of newly formed landslides and the cumulative landslide distribution in each survey year are calculated (Figures 4 and 5). Figure 5 is then used to analyse the relationship between the collapse rate of the watershed and the concavity index in each stage. The impact of terrain uplift on the landslide is also explored by using the determination coefficient (R 2 ) and slope value (gradient index) of the correlation analysis ( Figure 6). The closer the waterfront area is to the grid point of the channel, the easier it is for the sediment on the slope surface to enter the channel and contribute sediment to the river channel, and the higher the SDR is. The value of SDR ranges from 0 to 1, close to 0 is ridge landslide, and close to 1 is near-bank landslide. The spatial distribution of SDR in the basin (Figure 7) overlaps with the landslide map, which can be used to indicate whether the landslide location is off-bank or near-bank.  Figure 3. Relationship of slope-area data for each tributary watershed.

Multivariate analysis
There are three main factors that affect the overall shape of the river's longitudinal profile. They are flow changes, river bed grain size, and sediment transport volume, respectively. Factors affecting local changes include lithology, the influx of tributaries, and neotectonic movement (Chen 2004). The flow rate can be reflected in the size of the catchment area and the slope of the river course. The particle size of the river bed is related to the erosion resistance of the river bed. However, in the short term,   the above-mentioned influence factors have very small changes in the landform and cannot be explained by sufficient quantitative data. Therefore, this study attempts to  use the topographic information reflected in the watershed to describe the changes in the h value.

Stepwise regression
Since the h value is highly correlated with k s , only the h value is discussed to explain the changes in geomorphic indices. Based on the data in Table 3, the h value of the watershed is used as the dependent variable (Y 1 ), and the topographic factors of the watershed are independent variables (X 1 $X 17 ). The results of stepwise regression analysis show that at the 5% level of significance, any factor in the watershed is not significant (p> .05), indicating that the individual topographic factor of the watershed could not effectively describe the concavity index. However, if the significance level of the individual factor is increased (p¼ .082), the interpretation rate of the model constructed by both form factor and stream order of the watershed can reach 41.3% (Table 4).

Ã The maximum absolute correlation between each variable and any distinguishing function. 3.4.2. Cluster analysis
The concavity index of each watershed of Zhuoshui River Basin has no regularity from upstream to downstream, so cluster analysis is used to classify and show the difference in the spatial distribution of the concavity index. The purpose of cluster analysis is to classify variables or observations, and classify the concavity indicators of analysed watersheds into three levels (high, medium, and low). After K-Means cluster analysis, cluster centres of various levels can be obtained (Table 5). The classification results of concavity indicators for each watershed are shown in Figure 8.

Discriminant analysis
The concavity index represents the degree of concavity of the balanced channel profile, and the correlation between the topographic factors of the watershed is used for multivariate analysis to evaluate the concave curvature of the channel profile. The area, perimeter, length, average elevation, relief, average slope, relief ratio, width, compactness, circularity ratio, elongation ratio, drainage density, stream frequency, and stream order of watershed's topographic factors are independent variables. Taking the level of the concavity index as the grouping variable, the eigenvalues and the structure matrix obtained after the canonical discriminant analysis are shown in Tables 6 and 7. The variables explained by the first function axis are relief ratio (X 9 ) and length (X 3 ). The second function axis mainly explains drainage density (X 15 ), average elevation (X 6 ), number of rivers (X 5 ), average slope (X 8 ), area (X 1 ), total length of river (X 4 ), and perimeter (X 2 ) of watershed. The larger the absolute value of the correlation coefficient in the structure matrix, the greater the influence of this variable on this function. Table 8 shows the classification function coefficients of canonical discriminant function. Observations can be directly classified, that is, each observation value is substituted into each classification function. The greater the discriminant score (F value) obtained, the higher the group is. The classification function is expressed as follows:      F 2 medium h value ð Þ ¼ 1:061X 1 þ 8:270X 2 À 10:445X 3 À 0:116X 4 À 10:602X 5 À 0:547X 6 þ 38:005X 8 À 63:834X 9 þ 2840:937X 15 À 1941:985 F 3 low h value ð Þ ¼ 23:12X 1 þ 8:174X 2 À 7:284X 3 À 1:118X 4 À 12:085X 5 À 0:501X 6 þ 33:194X 8 þ 732:227X 9 þ 2827:130X 15 À 1978:019 Table 9 shows the standardized canonical discriminant function coefficients. The discriminant function is expressed as follows: D 1 ¼ 16:865X 1 þ 2:161X 2 þ 2:103X 3 À 6:058X 4 À 11:555X 5 À 0:211X 6 À 1:930X 8 þ 4:875X 9 þ 1:287X 15 D 2 ¼ À8:060X 1 þ 4:962X 2 À 2:017X 3 þ 12:884X 4 À 9:843X 5 À 6:245X 6 þ 9:989X 8 À 1:912X 9 þ 2:723X 15 The discriminant score (D value) can classify the observation value. Those with a high score in the first discriminant function can be classified as a low h value; those with a lower score both in the first discriminant function and in the second discriminant function can be classified as high h value; and those with higher scores in the second discriminant function can be classified as a medium h value. Therefore, the first discriminant function can effectively distinguish the low h value group from other groups, and the second discriminant function can effectively distinguish the medium h value group from other groups (Figure 9). The discriminant analysis results are shown in Table 10. The application of topographic information can effectively classify the level of the h value of the watershed, and the overall accuracy of the interpretation can reach 100%. Therefore, the analysis of topographic factors can quickly understand the index of mainstream concavity, which can be used as a basis for judging the magnitude of the uplift and erosion of the watershed.

River incision models
Taiwan's mountainous orogenic activities are intense, and the annual rainfall is abundant. After a long period of geomorphic evolution of Zhuoshui River, the longitudinal The jumping data in the S-A double logarithmic map are affected by the influx of tributaries, changes in rock formations, local rainfall, artificial dams or reservoirs, and the evolution of knickpoints. It can be seen from the relationship between the h value and k s is shown in Figure 10, showing a high degree of positive correlation. Therefore, the greater the concavity index, the greater the steepness index of the watershed. Since the correlation between the two shows a very significant positive correlation, the subsequent analysis only uses the concavity index as a quantitative indicator of terrain uplift caused by plate movement. After a long period of topographic evolution, the longer the river develops, the greater the degree of concavity in the upper reaches of the river, and the shape of the river's longitudinal profile will also be closer to the profile of equilibrium (Davis 1909). Whipple (2004) pointed out that moderate concavities (0.4-0.7) are associated with actively uplifting bedrock channels in homogeneous substrates experiencing uniform (or close to uniform) rock uplift, and low concavities (<0.4) are associated either with short, steep drainages heavily influenced by debris flows or with downstream increases in either incision rate or rock strength, commonly associated with knickpoints. Table 2 shows that the average value of concavity index of the Zhuoshui River watershed is 0.5, a moderate concavity. The values of tributary watershed concavity are between 0.20 and 0.71 except the Kashe, Junda, Zhuogun, and Shuili stream watershed, which are less than 0.4. There are faults pass through downstrem of the Kashe, Junda, Zhuogun, and Shuili stream (Figure 11). These results are consistent with Whipple's description of concavities.
According to the longitudinal section of the 12 tributaries of Zhuoshui River, the mathematical model fitting results are mainly logarithmic/exponential profile types ( Table 2). The profile fitting function of rivers with larger concavity is mainly logarithmic, such as Tarotwan Stream, Wanda Stream, Liqi Stream, Chenyoulan Stream, Qingshuigou Stream, and Qingshui Stream, and the profile fitting function of streams with smaller concavity is linear, such as Kashe Stream. A study by Snyder et al. (2000) in the Mendocino triple junction area of California pointed out that the change of log k s value is related to the uplift rate of the ground body where the watershed is located. In addition, the k s is positively correlated with the concavity index. Therefore, in a watershed with relatively high h value, it can be seen that the uplift rate of the strata has a decreasing trend from east to west and from south to north (Figure 8), this is because the distribution of the orogenic vector field in Taiwan is southeast to northwest (Dadson et al. 2003;Huang et al. 2015). Figure 11 shows that the Lishan Fault passes through the upper reaches of the Zhuoshui River Basin (Lee and Tsai 2010), causing the Junda channel and the outlets of Danda and Kashe streams to absorb Plate compression energy due to the fracture. The fault absorbs the compression energy of the Plates, thereby reducing the uplift.
The h value is an indicator of the degree of concavity of the balanced channel profile. When the terrain is in equilibrium, h ¼ 0 indicates that the elevation of the balanced channel profile has a linear relation with the flow length, and h ¼ 1 indicates that the elevation of the balanced channel profile has a power-function relation with the flow length (Yang, 2007), the larger the h value, the steeper the topography.

Relationship of concavity index and topographic factors
The results of stepwise regression analysis showed that the topographic factors of the watershed did not reach the significant level of 5% (p < .05), indicating that the individual topographic factor of the watershed could not effectively describe the concavity index. However, if the significance level is changed to p ¼ .082, the interpretation rate of the model constructed by both the stream order and the form factor of the watershed can reach 41.3% (p ¼ .037) ( Table 4), showing that the concavity index of the watershed can be roughly described by the two factors. From the regression formula Y 1 ¼ 2.709X 11 þ 0.121X 17 À 0.350, the slopes of the form factor (X 11 ) and the stream order (X 17 ) are all positive values, and the weight of the form factor is larger. The form factor is the width of the watershed per unit length of the main flow. When it is equal to 1, the watershed is round. On the other hand, a form factor value of less than 1 indicates that the watershed is narrower and longer, and a value greater than 1 indicates that the watershed is flat. A watershed is a water collection unit, and its form factor affects the shape of the flow history, which means it affects the transfer of energy. The form factors of the tributary watersheds of Zhuoshui River are all less than 1, which indicate a long and narrow watershed. However, the change of the  form factor can reflect the magnitude of energy transport in the watershed. Assuming that the boundary conditions of the river are the same, when the length of the main flow is fixed, the larger the form factor, the larger the area of the watershed, and the greater the concavity caused by the undercutting of the river channel by the collected water. When the area of the watershed is fixed, the smaller the form factor, the longer the main flow length, the slower the energy transmission of the river channel, and the smaller the concavity caused by the undercutting of the river channel.
If we want to use the terrain parameters of the watershed to effectively describe the concavity index, cluster analysis can be used to display it. The concavity indices of the analysed watersheds are divided into three categories (high, medium, and low), and then topographic factors can be used to effectively extract the concavity classification. Discriminant analysis shows that 100% of the original group can be classified correctly (Table 10), and the classification function coefficients of the canonical discriminant function variables (Table 8) can be used to judge the classification of the concavity index of any watershed with balanced channel. Figure 12 shows that the greater the number of classifications, the lower the classification accuracy. If the analysed watersheds are divided into nine categories, the classification accuracy is reduced to about 82%, which can still meet the management requirements. Rivers have directional headwater erosion. The older the river, the longer its water system will develop. In fact, there is no significant correlation between the concavity index and the mainstream length. Therefore, the evolution time is not a factor that affects the size of the concavity index. Rãdoane et al. (2003) studied the rivers of Carpathians and pointed out that the profile of the river is not affected by time, but is instead influenced by geological uplift time.  Figure 14. Frequency distribution of slope (left) and SDR (right) of the landslide grid in each survey year after the Chi-Chi earthquake.

Value-added application
After the Chi-Chi earthquake, every typhoon torrential rain resulted in a large amount of sediment output, which is an important process of geomorphic evolution. The correlation of collapse rate and the concavity index showed that the concavity index was positively correlated with the collapse rate ( Figure 6). The greater the concavity of the river profile, the higher the collapse rate, and the two are mutually causal. However, the coefficient of determination (R 2 ) of correlation analysis in different periods shows that the interpretation rate of the influence of the river concavity index on the collapse rate of the analysed watershed increased rapidly from 0.15 in the earthquake year to 0.47 in 2001, and then gradually decreased. The interpretation rate of the concavity index for the collapse rate can be reduced to its original value about 8 years after the earthquake (Figure 13left). The slope change of the correlation equation also shows the same trend ( Figure 13right), which indicates that the ground can collapse quickly after the earthquake to adjust the slope instability caused by the uplift. The correlation between R square and slope shows that the higher the determination coefficient, the larger the gradient index, and the higher the determination coefficient, the closer the linear relationship of the gradient index.  It can be seen that the concavity index can effectively reflect the evolution trend of watershed landslide over time. The slope (steepness), the SDR (the degree of closeness to the river system), and the stream order (upper, middle, and lower reaches of the basin) of the landslide in different periods were used for cluster analysis. The slope, SDR, and stream order are combined with their respective mean value (l) and coefficient of variation (cv) for a total of six features, and a 5-cluster variance analysis is performed. Because the stream order did not reach a significant level of 5% (Table  11), it was discarded and only slope and SDR were used for subsequent analysis.  Frequency distribution of slope and SDR of the landslide grid in each survey year after the Chi-Chi earthquake ( Figure 14) shows that the initial landslide-slope distribution curve after the earthquake is on the left, mainly due to the flat ridge landslide. The subsequent heavy rain landslides tend to move to the right side (more steeply), and recently changed back to the middle. The landslide area tends to gradually decrease over time. The slopes prone to landslides are distributed between 20 and 50 . This is basically consistent with Liu's study of the northern Chenyulan basin in 2001. The SDR distribution shows that from the off-bank landslide at the beginning of the earthquake to the near-bank landslide on the right side, there is a recent trend back to the middle, which is similar to the distribution trend of the slope. The prone areas are mainly distributed in SDR < 0.5, but concentrated between 0.1 and 0.4.
From the feature of l and cv for annual landslide's slope and SDR, the characteristics of cluster 1 with low slope-l, high slope-cv, low SDR-l, and high SDR-cv are mainly off-bank landslide; cluster 5 has the opposite characteristics of cluster 1, which mainly occurred near-bank landslide, and the characteristics of the middle three clusters show a transitional value (Table 12). The classification results of landslide in each survey year are shown in Table 13. The clusters of landslides observed within 7 years after the 1999 earthquake vary significantly, but the clusters of landslides more than 7 years after the earthquake have little difference. Since the temporal distribution of landslide clusters has the same trend as the concavity index ( Figure 13left), and the spatial distribution of landslide clusters (Figure 15) shows that the order of the clusters is the same as the order of evolution, this proves that the concavity, slope, and SDR index of the basin are key indicators, which can effectively explain the spatial distribution of landslide evolution in the basin. Based on cluster category in Table 12, the spatial distribution of landslide location types can be classified into three categories. Type 1 (cluster 1) refers to the initial slope failure caused by the crustal uplift only, and Type 2 (cluster 2-4) is the unstable slope caused by the earthquake and the subsequent rainfall-driven collapse. Type 3 (cluster 5) refers to near-bank landslides caused by insufficient river crosssection due to debris accumulated from Types 1 and 2 landslides. Table 13 shows the gradual transformation from scattered off-bank landslides to concentrated nearbank collapses.
The trajectory of the changes in the landslide from the Chi-Chi earthquake to the present can provide an important reference for the management of watershed landslide. Plate movement (internal force) will cause an earthquake to lift the crust and cause the slope to become unstable. The ground surface will collapse due to the stability of the slope by external forces such as gravity, weathering, and precipitation. In the process of orogenic movement, Taiwan faced different types of sediment disasters such as landslides, dammed lakes, and debris flows ( Figure 16). Watershed landslide caused by earthquakes and rainfall will affect protected targets. It is necessary to understand the occurrence and evolution mechanism of watershed landslide, screen relevant environmental indicators, establish risk models for estimating disaster hotspots, and explore the possibility and scale of hotspot disasters. The concavity, slope, and SDR proposed by this research can be used to extract areas prone to sediment disasters, and can be used as quantitative indicators for slope land use and disaster prevention management. Based on this, environmentally sensitive areas can be delineated for reference in policy governance.

Conclusion
The natural forces of Taiwan's landslide evolution are mainly affected by the interaction of crustal orogeny, precipitation, and runoff, and cannot be discussed by a single factor. In this study, the Chi-Chi earthquake occurred once in 100 years in 1999 as the base period, and the landslide data of 17 survey periods (1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017) in the Zhuoshui River Basin were used to explore the spatial distribution of landslide evolution based on its environmental characteristics. On the whole, the Zhuoshui River Basin has reached topographical balance. The results of this study show that under the condition of topographic balance, h and ks are highly positively correlated, indicating that the greater the concavity of the longitudinal profile of the watershed, the greater the topographic uplift, and the stronger the water erosion. The size of the concavity and the collapse rate is causal to each other. The determination coefficient and slope (gradient index) of the correlation analysis show that the ground can quickly collapse after the earthquake to adjust the slope instability caused by the uplift, and the current topography is the result of long-term geomorphological evolution. This research uses the concept of geomorphic balance to obtain quantitative indicators and proposes the correlation with topographic information. It provides a simple and fast method for distinguishing the spatio-temporal distribution of unstable slopes in watersheds.