Risk Assessment of Land Subsidence in Kathmandu Valley, Nepal, Using Remote Sensing and GIS ()
1. Introduction
1.1. Background
Land subsidence is just a geological phenomenon either triggered by natural or anthropogenic activities but when this phenomenon has the probability of resulting harmful consequences or the expected loss (of lives, property, livelihoods, economic activities or environment) then it is considered as risk [1] . Risk factors are compounded by rapid increase in urban population and economic development [2] . The physical damage caused by land subsidence can be mainly categorized into two forms: damage on artificial (manmade) infrastructures and damage on natural systems. Significant damage is seen in areas corresponding to land subsidence occurrence.
The main damage on manmade infrastructure reported worldwide is mostly related to water transport structures [3] [4] . Since land subsidence affects the elevation of the ground, and because water transportation infrastructures are very sensitive to minor gradient changes, subsidence can hugely affect such structures. Other reported problems include damage to buildings and transportation facility (i.e. roads, bridges, railways).
The damage to manmade infrastructures is more emphasized and noticeable unlike to the damage to natural systems, which is invisible and generally more threatening. The main reason is that artificial infrastructure damages can generally be repaired opposed to natural system damage which is generally permanent. Some of the examples of natural structure damage are permanent compaction of aquifer system, change in topography which ultimately affects the river patterns and low lying areas.
The other main factor affected by land subsidence is damage to the social environment which includes the human society and the economic development level. The physical damage caused by land subsidence will eventually affect the social environment directly or indirectly but the intensity is determined by the recoverable capability of life, property and various economic activities in the disaster affected areas. Remarkable economic losses have been caused by land subsidence throughout the world [5] [6] [7] .
Kathmandu valley, the capital and the urban core of a developing country Nepal is lagging in terms of data documentation and research work regarding land subsidence and its risk assessment. The factors that make a location prone to land subsidence risk (i.e. geology and groundwater extraction characteristics) are in favor of the valley, yet research is not being conducted. Also, the valley is experiencing rapid increase in population and economic development in the past few decades that will ultimately contribute to increase in risk of damage induced by land subsidence if no counter measures are considered.
Therefore, it is necessary to assess land subsidence risk for decision and policy makers to prevent a huge potential disaster. Risk assessment is simply an application of a methodology for evaluating risk, where risk is defined as the probability and frequency of occurrence of a hazardous event, exposure of people and property to the hazard and consequences of that exposure [1] .
Most frequently deployed approach for land subsidence risk assessment are by the means of Geographic Information System (GIS) techniques and Disaster Risk Index Method and Analytic Hierarchy Process (AHP) [1] [7] [8] [9] . GIS provides robust tools for inclusive spatial modeling and analysis. Disaster Risk Index method is an approach where the hazard, the vulnerability and the capability of disaster prevention and reduction are considered for the quantitative evaluation of a risk. AHP is a multi-criteria mathematical evaluation method used for decision making where hierarchical structures are used to quantify relative priorities for a given set of elements on a ratio scale set by the user [1] .
The main objective of this study is to assess land subsidence risk in Kathmandu valley, Nepal, by using Geographic Information System (GIS) techniques, Disaster Risk Index Method and Analytic Hierarchy Process (AHP). Land subsidence map of Kathmandu valley for 2007 to 2010 was generated by applying the Differential Synthetic Aperture Radar Interferometry (D-InSAR) technique before conducting this study. Peer-reviewed paper related to this study can be found at http://www.mdpi.com/2073-445X/6/2/39.
1.2. Study Area
Kathmandu valley is the capital city of Nepal and hence has been a center of ever growing economic activities from a very long time. Kathmandu valley is an urban agglomerate with a core urban center surrounded by extended urban economic zones. The lack of decentralization of developmental activities has propagated Kathmandu valley to be one of the most desired city to live in the country consequently, increasing the internal migration rates. The population density of Kathmandu valley is 2793 people per square kilometer as per the 2011 census [10] . The increasing population and failure in implementing strict regulation has resulted in haphazard development of the valley both in terms of infrastructure and economy. The Landsat image of Kathmandu valley is shown in Figure 1.
2. Methodology
2.1. Data Used
Land subsidence volume and velocity obtained from D-InSAR processing (a remote sensing technique) of ALOS PALSAR data from 2007 to 2010 was used in this research. These results were obtained by the authors in a previous study considered as the first part of this study (refer http://www.mdpi.com/2073-445X/6/2/39).
Groundwater exploitation intensity data was used for hazard mapping of land subsidence and was provided by the Kathmandu Valley Water Supply Management Board (KVWSMB), Ministry of Water Supply and Sanitation, Government of Nepal.
Population density, Gross Domestic Product (GDP) and Construction Land Proportion data was used for Vulnerability mapping. Population density data was obtained from the Kathmandu Valley Development Authority (KVDA), Government of Nepal. GDP data for the study area was obtained from National Accounts Section, Central Bureau of Statistics, Government of Nepal. Construction land proportion data was obtained from the Ministry of Land Reform and Management, Nepal.
2.2. Methodology
2.2.1. Methodology for Subsidence Volume Estimation
Through various literature review it was found that in most of the cases the shape of subsidence is very much like a cone shape. Therefore, the subsidence volume of each subsidence zone can be estimated by an assumption that the border extremities of a subsidence zone are linearly moving at a constant rate. A simple cone model designed by [11] is used to estimate the land subsidence volume (Figure 2).
Subsidence volume represented by the shaded portion in Figure 2 can be
![]()
Figure 1. Landsat image of Kathmandu valley.
![]()
Figure 2. Cone model for volume estimation of land subsidence.
estimated by applying the formula of a volume of a cone which is shown in Equation (1) [12] .
(1)
where, V is the subsidence volume to be estimated, A1 is the upper base area, A2 is the lower base area and h is the height or the perpendicular distance between the surface A1 and A2. Area A1 and A2 was calculated in ArcGIS (Version 10.4.1) by converting each subsidence zone obtained by DInSAR processing into a shapefile. The subsidence depth (h) and the velocity of land subsidence was also obtained from the DInSAR processing result. The methodology for D-InSAR processing to map land subsidence in Kathmandu valley from 2007 to 2010 has been explained in detail in the first part of this study which can be found at http://www.mdpi.com/2073-445X/6/2/39.
2.2.2. Methodology for Risk Assessment
The main objective for assessing risk of land subsidence is to link the subsidence phenomenon with the damage it causes to the physical as well as social environment. The main factors that help to determine risk are hazard and vulnerability. Therefore, a hazard map and vulnerability map were generated and then combined to obtain the final risk map. The detailed methodology is explained as follows.
Disaster Risk Index Method
Risk assessment is an approach for evaluating risk, where risk is defined as the probability and frequency of occurrence of subsidence, exposure of people and property to the subsidence and consequence of that exposure [1] . The degree of risk of land subsidence significantly depends on two factors: hazard and vulnerability [7] . As per the disaster risk index method, quantitative risk can be estimated by using Equation (2) [13] .
(2)
where, DR is the disaster risk, H is the hazard, V is the vulnerability and C is the capability of disaster prevention and reduction.
Hazard, in general can be defined as a phenomenon that has the potential to disrupt and damage people, property and their immediate environment. The hazard of land subsidence refers to the intensity and the probability that land subsidence will occur in a certain area in a certain period. As defined by the Asian Disaster Preparedness Center, Land subsidence hazard evaluation simply is the process of determining the degree of severity and the extent of the impact area. In this research, accumulative subsidence volume, land subsidence velocity and groundwater exploitation intensity was used as the indicators to evaluate hazard in the study area [14] . The former two indicators were obtained from the DInSAR processing results and the groundwater exploitation intensity data was obtained from the Kathmandu Valley Water Supply Management Board (KVWSMB), Ministry of Water Supply and Sanitation, Government of Nepal.
Vulnerability, in general can be defined as a concept that describes the factors (including economic, social and physical) aiding to reduce the ability to cope with the potential hazards impacts. The vulnerability of land subsidence refers to the measure of susceptibility to physical harm or damage caused due to land subsidence. Vulnerability includes the ability of the human society and the economic development level of the society to cope with the disaster caused by land subsidence. Land subsidence vulnerability evaluation is the process of assessing the sensitivity of the economy, population and physical infrastructure to the land subsidence phenomenon. In this research, population density, gross domestic product (GDP) and construction land proportion data was used as the indicators to evaluate vulnerability of the study area [15] . Population density refers to the number of people per unit area and this data for the study area was obtained from the Kathmandu Valley Development Authority (KVDA), Government of Nepal. GDP is one of the primary indicators used to evaluate the economic condition of a country and this data for the study area was obtained from National Accounts Section, Central Bureau of Statistics, Government of Nepal. Construction land proportion data gives the information of the proportion of built-up space and open space. This data was obtained from the Ministry of Land Reform and Management, Nepal.
The capability of disaster prevention and reduction refers to ability of the country to prevent or reduce the effect of potential land subsidence on life, property and economy. However, the land subsidence monitoring of the study area being very poor it was assumed that the country has no ability to control land subsidence at present.
Analytic Hierarchy Process (AHP)
Analytic hierarchy process is a multi-criteria mathematical decision making process developed by Professor Thomas Saaty in (1977). This process uses hierarchical structures to derive relative priorities for criteria (indicators) employing pair wise comparisons. In this research, this process was used to give weights to the indicators identified for evaluating hazard and vulnerability in MSExcel. The basic procedure includes the following steps:
Pair-wise comparison: Pair-wise comparison matrix for hazard and vulnerability was developed separately to establish priorities among the indicators. The result of comparison is derived in terms of integer.
Normalization: The integers obtained from the above step is normalized to compute the priority vector which gives the relative weights among the indicators and ultimately helps to decide which indicator is relatively more important in determining land subsidence risk. Normalization generally means to average the values in each row to compute the corresponding weight.
Consistency Analysis: The main objective of this step is to check if the preference ratings made in the pair-wise comparison are consistent. This is measured in terms of Consistency Ratio (CR), which can be calculated using Equation (3) [16] .
(3)
where, CI is the consistency index and RI is the random inconsistency indices. RI is provided for each order of matrix by [16] . In this research, since the order of matrix was 3, the corresponding RI value 0.58 was used. The CI value can be calculated using the Equation (4).
(4)
where,
is the value obtained from the summation of product of each normalized weight and sum of columns of the reciprocal matrix, n is the number of indicators used. Saaty, (1980) suggests that the CR value equal to 0.1 or below shows that the comparison is consistent and hence acceptable. (For detail description of the methodology refer “The analytical hierarchy process” Saaty, 1980).
After obtaining the weights for each indicators Hazard map and Vulnerability map was generated in ArcGIS (Version 10.4.1). These two maps were then utilized to obtain the final Risk map by using Equation (2). The methodology flowchart for risk assessment is shown in Figure 3.
3. Result and Discussion
3.1. Volume Estimation of Land Subsidence
Figure 4 shows the interferogram image for Kathmandu valley generated by utilizing DInSAR processing to ALOS PALSAR data acquired on 2007/11/02 and 2010/02/07. The locations of subsidence are indicated by L1, L2, L3, L4, L5, L6, L7 and L8. Location L1 represents the central Kathmandu area, L2 is Chauni area, L3 is Lalitpur area, L4 is Imadol area, L5 is Thimi area, L6 is Madhyapur Thimi, L7 is New Baneshwor and Koteshwor area and L8 is Gothatar area. The details of this result have been described in “Detection of Land Subsidence in Kathmandu Valley, Nepal, using DInSAR Technique (Bhattarai et al., 2017)” (refer http://www.mdpi.com/2073-445X/6/2/39 ) considered as the first part of this study.
The subsidence volume for each location was estimated using Equation (1). The subsidence volume estimation for each location is shown in Table 1. These values were further used for hazard mapping of land subsidence in this research.
![]()
Figure 3. Methodology flowchart for risk assessment.
![]()
Figure 4. Differential interferogram for Kathmandu valley from 2007/11/02 to 2010/02/07 showing the subsidence locations L1 - L8.
![]()
Table 1. Volume estimation of land subsidence locations.
3.2. Risk Assessment of Land Subsidence
An indicator framework based on the Disaster Risk Index method was prepared with two significant factors for risk assessment along with three indicators for each factor. The indicators were weighted by its significance to govern land subsidence risk using the AHP multi-criteria decision making process. The values of each indicator were classified into three classes Low, Medium and High where low indicates the lowest hazard/vulnerability level and high indicates the highest hazard/vulnerability level. Table 2 shows the weighted values and grade and values for the respective indicators.
3.2.1. Land Subsidence Hazard Evaluation
Kathmandu valley was divided into three zones corresponding to the major cities Kathmandu, Bhaktapur and Lalitpur along with the subsidence locations obtained from previous processing for hazard evaluation. The data values of indicator for each zone and location is shown in Table 3. The volume and maximum subsidence velocity for the three zones were considered to be nil as no subsidence was detected in these areas.
![]()
Table 2. Weighted values and grade and values for indicators of land subsidence risk assessment.
![]()
Table 3. Data values of indicators used for land subsidence hazard evaluation of Kathmandu valley.
The land subsidence hazard map of Kathmandu valley is shown in Figure 5. Land subsidence in Kathmandu valley mainly situates in low hazard area (indicated by yellow colour in Figure 5). High hazard areas can be seen scattered in the subsidence locations L1, L3 and L4 (indicated by red colour in Figure 5) whereas medium hazard areas can be seen scattered in the subsidence locations L2, L5, L6, L7 and L8 (indicated by green colour in Figure 5).
This indicates that there is a high probability that land subsidence occurrence is bound to intensify in locations L1, L3 and L4 considering there is no reduction in the average groundwater discharge.
3.2.2. Land Subsidence Vulnerability Evaluation
In the same manner to hazard evaluation, Kathmandu valley was also divided into three zones corresponding to the major cities Kathmandu, Bhaktapur and Lalitpur along with the subsidence locations obtained from previous processing. The data values of indicator for each zone and location is shown in Table 4. The individual subsidence locations were categorized into their respective municipality since it is difficult to acquire social and economic data (i.e. population density and GDP) for such small zones.
The land subsidence vulnerability map for Kathmandu valley is shown in Figure 6. The evaluation results show that subsidence location L1, L2, L3 L7 and L8 (indicated by red colour in Figure 6) are highly vulnerable areas. These locations are situated in the main urban core of the Kathmandu valley where the population density is the highest and the economy is the most developed. Kathmandu zone, Lalitpur zone and subsidence locations L4, L5 and L6 are found to be in medium vulnerable zone (indicated by green colour in Figure 6). Lowest
![]()
Figure 5. Land subsidence hazard map of Kathmandu valley generated through GIS processing.
![]()
Table 4. Data values of indicators used for land subsidence vulnerability evaluation of Kathmandu valley.
![]()
Figure 6. Land subsidence vulnerability map of Kathmandu valley generated through GIS processing.
vulnerability is seen in Bhaktapur zone as the population density and economic activity is lowest in this area. The result indicates that Location L1, L2, L3 L7 and L8 are most sensitive to damage caused by land subsidence.
3.2.3. Land Subsidence Risk Assessment
Land subsidence risk map of Kathmandu valley was generated based on the hazard and vulnerability evaluation in ArcGIS (Version 10.4.1) using Equation (2) (Figure 7). The areas where land subsidence was detected through DInSAR processing were found to be in high (Location L1 and L3; indicated by red colour in
![]()
Figure 7. Land subsidence risk map of Kathmandu valley generated through GIS processing.
Figure 7) and medium (Location L2, L4, L5, L6, L7 and L8; indicated by green colour in Figure 7) risk areas. However, rest of the Kathmandu valley was found to be at low risk of land subsidence (indicated by yellow colour in Figure 7).
Shrestha, et al., (2017) [17] also predicted that Kathmandu valley is at low risk through a model-based estimation of land subsidence. Even though the prediction has not been validated, it resembles the result of this study.
Due to lack of proper scientific data and research, evidence of damage caused by land subsidence has not been reported in Kathmandu valley till date. Case studies from around the world can be referred to utilize the knowledge and experience in planning and policy making to reduce if not prevent the disastrous effect of land subsidence. The location of Kathmandu valley has a very close resemblance with the location of Mexico City where land subsidence has been documented very well. Hence, case study of land subsidence in Mexico City is discussed here.
Mexico City is the capital city of Mexico located in the valley of Mexico Basin. It is surrounded by volcanic chain and mountains with elevations reaching up to 5000 m [18] . Like the origin of Kathmandu Valley, Mexico City emerged where there was once Lake Texcoco, therefore the geology consists of highly saturated clay [19] . The geology is classified into three zones namely Foothill zone, Transition zone and Lake zone. The Foothill zone comprises heterogeneous volcanic deposits and lava. The Transition zone mainly comprises of sand and gravel alluvial deposits along with volcanic materials. The lake zone comprises of highly compressible lacustrine clays [20] .
The city has been suffering from groundwater extraction related subsidence since decades [21] . The city has been reported with subsidence reaching up to 38 cm/yr. [20] . The main cause has been identified to be the drying and compaction of soft clay layers which has less permeable capacity and is mainly triggered by excessive groundwater exploitation from the aquifers [22] .
The main problems reported relating to land subsidence in Mexico City are
・ Decrease in runoff and wastewater drainage ability, that ultimately results in flooding during rainy season [19] .
・ Disruption to the water transportation structures (i.e. pipelines, canals) resulting in interference with water supply [20] .
・ Damage to the stability of manmade infrastructure (Buildings, transportation facilities like highways, roads and bridges) due to changes in surface gradients [20] .
Since, there is a close resemblance in Kathmandu Valley and Mexico City, it can be expected that these problems may be encountered in Kathmandu Valley in near future if current situations prevail. The high-risk areas indicated by red color in Figure 7, Central Kathmandu and Lalitpur which are the main urban core of the valley has the highest probability to suffer from such damages. Other locations that are at medium risk to such damages are Chauni, Imadol, Thimi, Madhyapur Thimi, New Baneshwor and Koteshwor and Gothatar. Since, these damages are not limited to a point location, the periphery areas are also bound to suffer.
4. Conclusion
Land subsidence risk assessment based on current data revealed that the identified subsidence areas are at high and medium risk of suffering from subsidence induced damages like decrease in runoff and wastewater drainage ability, disruption to water transportation structures and damage artificial infrastructure stability whereas the rest of the Kathmandu valley is at low risk under the circumstance that similar conditions prevail. The outcomes of this research even though it is not directly validated can be used as a base for further detailed study. The results could also serve beneficial for developing disaster prevention policies. In addition, a more comprehensive risk assessment of land subsidence can be done by considering other indicators like the geological characteristics and land use type of the location. This would give a more detailed outlook on factors that can be controlled to reduce if not prevent a huge disaster. Also, in this study, it was assumed that the study area had no capacity for disaster prevention and reduction for risk assessment. However, different case scenarios like with government action to reduce groundwater exploitation and with government action to reduce or prevent construction land proportion can be employed to judge the situation in a different perspective.
Acknowledgements
The authors would like to thank Kathmandu Valley Water Supply Management Board (KVWSMB), Ministry of Water Supply and Sanitation, Kathmandu Valley Development Authority (KVDA) and National Accounts Section, Central Bureau of Statistics, Government of Nepal for providing necessary data to conduct this study. The authors would also like to express a hearty gratitude to Mr. Suyesh Bhattarai for his help in obtaining the required data.