Coastal Vulnerability Assessment of Bali Province, Indonesia Using Remote Sensing and GIS Approaches

Abstract

Coastal zones are considered to be highly vulnerable to the effects of climate change, such as erosion, flooding, and storms, including sea level rise (SLR). The effects of rising sea levels endanger several nations, including Indonesia, and it potentially affects the coastal population and natural environment. Quantification is needed to determine the degree of vulnerability experienced by a coast since measuring vulnerability is a fundamental phase towards effective risk reduction. Therefore, the main objective of this research is to identify how vulnerable the coastal zone of Bali Province by develop a Coastal Vulnerability Index (CVI) of areas exposed to the sea-level rise on regional scales using remote sensing and Geographic Information System (GIS) approaches. This study was conducted in Bali Province, Indonesia, which has a beach length of ~640 km, and six parameters were considered in the creation to measure the degree of coastal vulnerability by CVI: geomorphology, shoreline change rate, coastal elevation, sea-level change rate, tidal range, and significant wave height. The different vulnerability parameters were assigned ranks ranging from 1 to 5, with 1 indicating the lowest and 5 indicating the highest vulnerabilities. The study revealed that about 138 km (22%) of the mapped shoreline is classified as being at very high vulnerability and 164 km (26%) of shoreline is at high vulnerability. Of remaining shoreline, 168 km (26%) and 169 km (26%) are at moderate and low risk of coastal vulnerability, respectively. This study outcomes can provide an updated vulnerability map and valuable information for the Bali Province coast, aimed at increasing awareness among decision-makers and related stakeholders for development in mitigation and adaptation strategies. Additionally, the result may be utilized as basic data to build and implement appropriate coastal zone management.

1. Introduction

The coastal area is considered to be the most vulnerable area affected by climate change [1,2,3,4]. Climate change, with its associated rise in sea level and possible increases in the frequency and/or the intensity of storms and changes in wave climate, can be expected to significantly increase the risk of coastal erosion and flooding in most coastal locations [5,6,7,8,9,10,11,12,13]. According to Nicholls et al. [14], more than 200 million people worldwide are vulnerable in the case of flooding by extreme sea levels. As sea levels rise, coastal populations and people living in low-elevation coastal zones and small islands can expect more frequent and more severe high tides, flooding, and storm surges [15,16,17,18]. In addition, SLR seriously impacts densely populated coastal zones with a great deal of resources [19,20]. The loss of coastal ecosystems will have negative effects on tourism, infrastructure, freshwater supplies, biodiversity, and fisheries [8,21,22,23,24,25].

Since the early 1990s, high-precision altimeter satellites have routinely been used to measure sea level, showing that the global mean sea level is rising by 3.4 mm/year [26], and has seen a significant rise by 21–24 cm since 1880 [27]. The rising sea level is mostly due to the combination of ocean thermal expansion and the melting of sea ice. The change in sea level at regional and global scale may differ from local levels on a particular coast. Therefore, producing locally based studies is essential, particularly when their findings can provide more useful information for decision makers.

The effects of rising sea levels endanger several nations, and it potentially affects coastal populations and natural environments [28,29]. In this context, identifying the vulnerability of various coastal sectors to the impact of rising sea levels is indispensable for coastal zone management. Various methods have been developed to identify coastal vulnerability over the past two decades. The most common and simple method for assessing coastal vulnerability concerning sea-level rise is that based on the calculation of an index that aggregates a set of parameters representing several spatial entities (geographic data) that influence coastal vulnerability. The Coastal Vulnerability Index (CVI) was introduced by Gornitz [30], and has been adopted and/or applied to assess coastal vulnerability globally [4,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53]. It uses a ranking system; therefore, it is easy to identify the degree of vulnerability. The model developed by Gornitz [30] is composed of seven variables to determine physical vulnerability in the USA to SLR impacts, consisting of relief (elevation), rock type (geology), landform (geomorphology), vertical movement (relative sea-level change), shoreline displacement, tidal range, and wave height.

A major focus of CVI studies is addressing geophysical vulnerabilities on the basis of remotely sensed data processed by means of the GIS methodology. The combination of remote sensing and GIS technology has been proven to provide valuable data for analysis of a given scenario [54], and can be used to help identify disaster-prone regions [55]. This method has been developed and used in several countries to assess coastal vulnerability, e.g., Canada [31], USA [56], UK [57,58], Spain [6,59,60,61], Portugal [62,63,64], Italy [65], Croatia [66], Greece [67], Turkey [45], Mediterranean coast [22,50], Lebanon [68,69], Nigeria [70], Ghana [11,71], South Africa [72], India [47,54,73,74,75,76,77,78], Bangladesh [79,80,81], South Korea [82], China [52,83], Malaysia [84,85], and Australia [46].

In Indonesia, coastal vulnerability has been studied at national and regional scales, whether conducted as individual research or in the context of group research. Mostly, assessments have been conducted using an index-based method. The national assessment of coastal vulnerability was performed by the Marine Research Center—Ministry of Marine Affairs and Fisheries (MMAF) (2009) and Badan Perencanaan dan Pembangunan Nasional (Bappenas, National Development Planning Agency) (2018) using CVI, even though there were differences in their determination of the index. Husnayaen et al. [86] conducted the physical assessment of coastal vulnerability under enhanced land subsidence, and the study proved that land subsidence had a significant influence on coastal vulnerability in Semarang. Imran et al. [87] revealed the coastal vulnerability index in the aftermath of the tsunami in Palu Bay, and concluded that, generally, the level of vulnerability corresponded to moderate vulnerability, both before and after the tsunami. In 2021, Irham et al. [88] determined the vulnerability of the west coast of Aceh Besar, Aceh Province, Indonesia using only four geological parameters, and concluded that the area possessed very high vulnerability, being generally formed of sandy beaches with a very gentle slope, while the areas with very low vulnerability had a high elevation and cliff beaches.

Bali Province is a well-known tourist destination that is dependent on sun-and-beach recreation activities. However, about 86 km or 20% of the length of the existing beaches has been eroded [89], and environmental degradation [90] has occurred due to natural factors as well as human activities. As a natural coastal defense system, beaches play an important role in reducing the risks of coastal erosion. Thus, their retreat and eventual disappearance increases their vulnerability to hazard events. In addition, beach narrowing threatens beach environmental services that are critical to the economy of tourist destinations, since recreational activities are dependent on the beach backshore [91,92,93]. Considering the threat of sea-level rise in coastal areas and on small islands, it is necessary to conduct a study to determine the degree of vulnerability experienced by a coast, since measuring vulnerability is a fundamental phase in achieving effective risk reduction [94]. Therefore, this study aims to quantitatively assess the vulnerability of the coastal zones of Bali Province to the actual rate of SLR by developing a CVI, considering the geological and physical characteristics of coastal processes.

2. Study Area

The study area (Figure 1) is located in the coastal zone of Bali Province, Indonesia, astronomically settled at 8°03′40″S–8°50′48 S and 114°25′53″E–115°42′40″E. Bali Province has a beach length of ~640 km, and about 18% are coral beaches with white sand. This small island makes up only 0.3% of Indonesia’s landmass and is home to 1.4% of Indonesia’s population [95]. Bali Island, as the mainland of Bali Province, has a beach length of ~593 km, and about 18% are coral beaches with white sand. Most beaches in Bali Province are characterized by sandy beaches, whereas the south of Bali Island and Nusa Penida are composed of cliff coast. Some particular areas, such as Denpasar Regency, Gilimanuk Bay, and the northern part of Nusa Lembongan Island, are characterized by vegetated beaches with mangroves.

Bali Island is dominated by a slope greater than 15%, which is predominantly found in the central section. This covers the mountain area that stretches from west to east, through Jembrana, Tabanan, Klungkung, Bangli, and Karangasem regencies, while a terrain slope of less than 15% is found in the Denpasar, Gianyar, and Badung regencies. There are five main types of soil, i.e., alluvial, regosol, grey-brown andosol, latosol, and Mediterranean. Bali’s climate is tropical. The wet season occurs between November and April, leaving May through October typically dry. The tides are of the mixed semidiurnal type, and are influenced by the Indian Ocean.

3. Methodology

The present study adopts an index-based methodology applying the CVI. CVI represents a combined result of the parameters influencing coastal vulnerability in the coastal area. The data were compiled from various sources to initially build a database of the parameters. The database is based on that used by Thieler and Hammar-Klose [32,33,34,35], Rao et al. [44], and Pendleton et al. [96], and loosely follows an earlier database developed by Gornitz [30].

Table 1. Sources and period of the different usage parameters in the compiling of CVI.

Parameter	Data Source	Resolution	Time Period
Geomorphology	Land use data and geology map by BIG https://portal.ina-sdi.or.id/downloadaoi/ (accessed on 20 December 2020)	Scale 1:25 k	2005
Shoreline change rate (m/year)	Sentinel—2A imagery https://earthexplorer.usgs.gov/ (accessed on 15 July 2020)	10 m	2015 and 2019
Elevation (m)	DEM imagery by BIG tides.big.go.id/DEMNAS/ (accessed 23 July 2020)	0.27 arcsecond ~8.1 m	-
Sea level change rate (mm/year)	Tide gauge and satellite data https://ccar.colorado.edu/altimetry/index.html (accessed on 13 June 2020)	1/6th deg	1992–2019
Tidal range (m)	Tide gauge data https://ccar.colorado.edu/altimetry/index.html (accessed on 2 August 2019)	-	1998–2019
Significant wave height (m)	Marine Copernicus data https://marine.copernicus.eu (accessed on 29 September 2020)	0.2 deg ~22.2 km	1993–2019

provides information on the data source and period for each parameter used in this study.

The data values of parameters were assigned a vulnerability ranking based on values ranges contributing to coastal vulnerability, while the non-numerical geomorphology parameter was ranked qualitatively according to the relative resistance of a given landform to erosion. Each parameter input was assigned an appropriate risk level (

Table 2. Risk rating assigned for different parameters.

Parameters	Very Low 1	Low 2	Moderate 3	High 4	Very High 5
Geomorphology [44,96]	Rocky, Cliff coast, Fjords	Medium cliffs, Intended coasts	Low cliffs, Glacial drift, Alluvial plains	Cobble beaches, Estuary, Lagoon	Barrier beaches, Sand beaches, Saltmarsh, Mudflats, Deltas, Mangroves, Coral reefs
Shoreline change rate (m/year) [30]	≥2.1 Accretion	1.0 to 2.0 Accretion	−1.0 to 1.0 Stable	−1.1 to −2.0 Erosion	≤−2.0 Erosion
Elevation (m) [30]	≥30.1	20.1–30.0	10.1–20.0	5.1–10.0	0.0–5.0
Sea level change rate (mm/year) [30]	≤−1.1 Land rising	−1.0 to 0.99 Land rising	1.0 to 2.0 within range of eustatic rise	2.1 to 4.0 Land sinking	≥4.1 Land sinking
Tidal range (m) [30,44]	≤0.99	1.0–1.9	2.0–4.0	4.1–6.0	≥6.1
Significant wave height (m) [44,96]	<0.55	0.55–0.85	0.85–1.05	1.05–1.25	>1.25

) based on its ability to cause very low, low, moderate, high, and very high damage, respectively, for a particular area of the coastline) [30,44,96]. Later, the key parameters were integrated into a single index and categorized based on the relative intensity of risk to the coast. A flow chart of the methodology to obtain the CVI map of Bali Province is provided in Figure 2.

The parameters used are dynamic, and have varying resolutions from different sources. Hence, each parameter was measured in a grid cell known as a shoreline grid. This shoreline grid had a size of 1 km × 1 km (Figure 1). Thus, 681 shoreline grids with an identification number were created to assign the vulnerability scale for each specific parameter. We assumed that this size would be sufficient for uniform data scaling. As a result, certain data had to be upscaled or downscaled. This grid cell is widely used to conduct the coastal vulnerability assessment [49,61,79,86,97].

In this study, six parameters were considered for developing the CVI: geomorphology, shoreline change rate, elevation, sea level change rate, tidal range, and significant wave height. The importance of each parameter and the procedure for generating the values in assessment CVI are given in the following section.

3.1. Geomorphology

The geomorphological mapping was derived from the land use data and land cover (LULC), since geomorphology plays an important role in the land surface. Hence, the association of land use, land cover and geomorphology were used for the visual interpretation and integrated with geology maps. In some cases, geology and geomorphology were included as one variable, because geomorphology often includes both the landform and the landform rock type [83]. The LULC and geology of Bali Province are shown in Figure 3.

Geology, in other words, rock-type variables, is a parameter expressing the relative erodibility of different landform types [40]. Bedrock lithology, shore materials, and coastal landforms vary substantially in terms of their resistance to erosion [61]. A generalized scale of lithologic and geomorphologic resistance to erosion was discussed by Gornitz and Kanciruk [98]. Geomorphology is the study of landforms and geologic processes that lead to these landforms, and identifies similarities among landforms [74]. Of particular importance to coastal vulnerability is the fact that geomorphology defines the erodibility of different landform types [96]. Erosion and rising sea levels will redistribute coastal landforms such as intertidal flats, marshes, dunes, etc. because of increased wave erosion and storm surges [74].

Coastal geomorphology refers to the surface type of a coastal area. Various surface types in coastal areas respond differently to coastal erosion. For example, sandy coasts are more vulnerable to erosion than cliff coasts. The geomorphology ranking is based upon the classification of Gornitz [30], shown in Table 2, according to which the following five coastal types were recognized along the coast.

3.2. Shoreline Change Rate

Coastal shorelines are naturally subjected to change due to certain coastal processes, including wave characteristics, near-shore circulation, sediment characteristics, or beach form [74]. This variable relates shoreline change to the rate of erosion or accretion to coastal vulnerability. Coastlines subjected to erosion are considered more vulnerable because of the loss of private land, infrastructure, and coastal habitats (e.g., beaches, dunes, and wetlands) [74]. Erosion processes also reduce the distance between coastal populations and the ocean, thus increasing the exposure of coastal populations to coastal vulnerabilities [74,99]. Coastlines subjected to accretionary processes are less vulnerable as additional land is created (e.g., accreting deltas, and sedimentations) [44].

A change in the location of the shoreline is an indication of the sensitivity of the coast to erosion. Coastal erosion is considered a risk not only because it threatens buildings and infrastructure, but also because it degrades and diminishes the extent of the beach, potentially impacting negatively on tourism [100]. This is unless the shoreline and beach are not obstructed from moving landward, which is often not the case in regions with a developed coast such as tourist destinations. Shoreline change rate can be defined as the rate at which the shore becomes eroded or accretes due to wave-action, sea-level rise or other hazards and processes that affect the land. The shoreline change rate was obtained using data for the years 2015 (as an initial condition) and 2019 (to describe the current condition) from the Sentinel-2 satellite imagery. At each grid cell along the shoreline, the shoreline in 2015 and 2019 were compared and estimated the rate of shoreline change. The rates of shoreline change correspond to the very low-, low-, moderate-, high-, and very high-risk rating categories (Table 2). According to this classification, grid cells (areas) with a shoreline change rate less than −2.0 m/year (erosion) are given a very high-risk rating, and vice versa.

3.3. Elevation

The elevation is the average elevation above the mean sea level (MSL) of a particular coastal area. The coastal elevation is important to include in vulnerability assessments because elevation can be utilized to (1) identify and estimate the potential extent of land threatened by inundation from sea-level rise, (2) estimate potential available land for wetland migration, and (3) identify sea-level rise impacts on human populations and infrastructure [74]. Coastal areas with higher elevations are considered less vulnerable, because higher elevations are more resistant to inundation from sea level rise or storm surges [74]. Therefore, a lower-elevation coastline is highly vulnerable to inundation and erosion.

Digital Elevation Model (DEM) is a digital cartographic dataset model of the earth’s surface in three (XYZ) coordinates and has been derived from contour lines or photogrammetric methods [101]. DEMs of Bali Province with 0.27 arcsecond (~8.1 m) resolution are distributed and packaged in small tiles to keep file sizes and processing times manageable [102]. The tile of DEM data has a number ID to download and is available with *.tif format. DEM data was processed to obtain the coastal elevation, where the elevation value is classified by the determination of the vulnerability index (Table 2).

Coastal areas with higher elevations are considered less vulnerable because higher elevations are more resistant to inundation from the sea-level rise or storm surge. Therefore, lower-elevation coastlines are highly vulnerable to inundation and erosion. The coastal elevation is an indicator of not only the relative risk of inundation, but also the potential speed of shoreline retreat. The coastal elevation is the main aspect to be measured along with the coastal morphology in approximating the coastal vulnerability [44]. The coastal elevation was used to replace the coastal slope and used by Gornitz [30], and has also been included in several CVI studies [31,74,103,104] to characterize the vulnerability of coastal areas due to inundation.

3.4. Sea Level Change Rate

Mean sea level is the average height of the sea surface over a longer period (usually a month or year), with the short-term variations associated with tides and storm surges averaged out [105]. Sea level rise shifts the wave action zone to higher elevations due to recession of coast to inland. The major physical impacts of sea level rise include significant shoreline retreat, inundation of deltas as well as flooding and loss of many salt marshes and wetlands [8,106,107]. Relative sea-level rise is a combination of eustatic (global), regional, and local changes in sea level [83,108]. Coastlines subjected to high relative sea-level rise rates are considered highly vulnerable areas due to the potential inundation of coastal land [74]. Coastlines with low relative sea-level rise rates are less vulnerable to inundation, and vice versa. Sea-level rise impacts shoreline change, geomorphology, land use, land cover, and groundwater resources, and also increases the inundation of coastal areas and the risk of flooding by storm surges [109].

Relative sea-level change data are a historical record. The mean sea level variation in Bali Province for 27 years was retrieved from satellite altimetry. Sea level data can be obtained not only from satellite altimetry, but also from the tide gauge. Satellites provide another independent source of information on sea level, in addition to tide gauge measurements. However, reliable satellite measurements of global sea levels did not become available until the early 1990s. Like the tide gauge data, the global mean sea level (GMSL) from satellite observation shows a clear upward trend in sea level, with year-to-year variability. Additionally, compared with the tide gauge data, the correlation of satellite data is high enough and still reliable, which is 0.86 or 86%. The trend from the satellite observations over this period is 3.3 mm/year.

The mean sea-level rise data were obtained from the Sea Level Explorer (https://ccar.colorado.edu/altimetry/index.html accessed on 13 June 2020) in *.csv format. The monthly mean sea level data from 1992 to 2019 were downloaded for Benoa Station and were computed to calculate the sea level change rate, and risk ratings were assigned.

3.5. Tidal Range

The tidal range of a coastal area is the vertical distance between the highest maximum high tide and the minimum low tide and is related to tidal flooding [42,43,74,83]. The tidal range varies temporally and spatially. The tidal range is linked to both permanent and episodic inundation hazards. The tendency in coastal vulnerability assessments is to designate coastlines with high tide ranges as highly vulnerable [74,97], which can cause the erosion and transportation of sediment in such a way. According to previous studies, the mean tidal range is a determinant of coastal vulnerability [11,42,43,110,111].

Conventionally, tides are observed by using a tide gauge. A tide gauge can provide accurate coastal tide data due to its spatial location [112]. Tidal data were obtained from the Geospatial Information Agency (BIG) website, and were generated from permanent tide stations. There are four tide stations located in Bali Province: Port of Benoa, Denpasar; Nusa Penida, Klungkung; Pengambengan, Jembrana; and Celukan Bawang, Singaraja. Hourly data from January 1998 to December 2019 from four tide stations were downloaded and annually averaged to produce the mean tidal range.

3.6. Significant Wave Height

Significant wave height (SWH) is similar to mean wave height. However, significant wave heights are the average height (trough to crest) of the highest one-third of waves within 12 h [74]. Mean significant wave height is used as a proxy for wave energy, which drives the coastal sediment budget [4,36,37,38,39]. As wave energy increases, the mobilization and transport of coastal material also increase. Coastal areas with higher significant wave heights are considered more vulnerable than coastal areas with lower significant wave heights [74]. Wave energy is directly related to the square of wave height:

$$\mathrm{E}=\frac{1}{8}{\mathsf{\rho }\mathrm{gH}}^2$$

(1)

where E is energy density, H is wave height, ⍴ is water density and g is the acceleration due to gravity. Thus, the ability to mobilize and transport coastal sediment is a function of wave height square. The significant wave height data was collected from Marine Copernicus. Down-scaling was applied to infer high-resolution information from 0.2° (~22.2 km) to 9 km.

3.7. Coastal Vulnerability Index (CVI) Analysis

The CVI analysis was carried out using the software package QGIS version 3.10, a free and open-source GIS application that makes it possible to create, visualize, process, and analysis geospatial data. Once each section of the coastline has been assigned a risk value for each parameter, the key parameters are integrated into a single index, through a mathematical formula by CVI. The CVI is calculated as the square root of the product of the ranked parameters divided by the total number of parameters [30] and represented as shown in Equation (2).

$$\mathrm{CVI}=\sqrt{\frac{\left(\mathrm{a}\ast \mathrm{b}\ast \mathrm{c}\ast \mathrm{d}\ast \mathrm{e}\ast \mathrm{f}\right)}{6}}.$$

(2)

where a = risk rating assigned to coastal geomorphology; b = risk rating assigned to shoreline change rate; c = risk rating assigned to coastal elevation; d = risk rating assigned to sea-level change rate; e = risk rating assigned to tidal range; f = risk rating assigned to significant wave height.

The range value of CVI can be classified into four equal parts [43,113], as quartiles [36,46] or percentiles ranges [31,77,114]. In this study, quantile classification was used to classify the CVI range into four risk level. It divides an equal-sized group of CVI range values into quarters. Therefore, the CVI scores in the lowest range were assigned low vulnerability, followed by moderate vulnerability, high vulnerability, and very high vulnerability for the highest range of CVI values. As a relative quantitative method, this method cannot directly express specific physical processes and effects, but it can be used as a diagnostic tool to determine which regions are most vulnerable to sea-level rise.

3.8. Model Validation

The fieldwork was designed as a model validation to integrate and apply evidence-based information from the assessment to identify and define the current condition of the study area. The fieldwork was conducted from 23 to 27 November 2020, and the sampling locations are shown in Figure 4. Before fieldwork is performed, sampling locations should be prepared in order to be able to achieve effective sampling and inspection. The Fitzpatrick–Lins equation [115] was applied to determine how many site locations to inspect.

$$\mathrm{N}=\frac{{\mathrm{Z}}^2\mathrm{pq}}{{\mathrm{E}}^2}$$

(3)

where: N = number of samples; Z = level of confidence according to the standard normal distribution (level confidence of 80%, Z = 2; level confidence of 95%, Z = 1.96; level confidence of 99%, Z = 2.575); p = estimated proportion of the population that presents the characteristic (when unknown p = 50); q = 100 − p; E = allowable error.

The survey design that was used to take a sample of beach information was purposive sampling based on the existing vulnerability class map. The purposive sampling technique is a type of non-probabilistic sampling method that is most effective when studying a certain characteristic with knowledgeable export within.

To determine the vulnerability class in fieldwork, data and information were collected by visual observation in consideration of the accessibility and predominant coastal conditions based on the aerial photographs (Google Earth), which can represent an area within a 1 km grid. Each sampling location had a database containing all information related to physical aspects (land use and/or geomorphology, substrate, elevation, and coastline condition). Simple tabulation was used to display the numerical data by logically using rows and columns to provide visual relationships and connections and transform raw data into information. Physical data (i.e., location, elevation, photos) were gathered using a GPS map.

To explain the suitability for qualitative items, the kappa coefficient of the agreement was introduced to remote sensing in the early 1980s as an index to express the accuracy of an image classification [116,117]. Similarly, Lillesand et al. [118] note that the Kappa coefficient is widely used in remote sensing as an accuracy assessment method to measure the agreement between interpretation and real conditions in the field. The kappa coefficient of agreement ($\hat{\mathrm{K}}$) was estimated, following Cohen [119], using the equation below.

$$\hat{\mathrm{K}}=\frac{{\mathrm{p}}_{\mathrm{o}}-{\mathrm{p}}_{\mathrm{c}}}{1-{\mathrm{p}}_{\mathrm{c}}}$$

(4)

where ${\mathrm{p}}_{\mathrm{o}}$ is the proportion of cases correctly classified (observed agreement) and ${\mathrm{p}}_{\mathrm{c}}$ is the expected proportion of cases correctly classified by chance (chance agreement).

The magnitude of $\hat{\mathrm{K}}$ lies on a scale from −1 to +1, but interest is typically focused only on positive values, because negative values are generally interpreted as ‘no agreement’ [120].

4. Results and Discussion

As a result of the CVI determination, it is possible to classify the physical vulnerability of the Bali Province coastal zone on the basis of six parameters. Figure 5 presents the results of the vulnerability ranking from each parameter.

A more detailed description of the coastal vulnerability rank of each parameter is given below.

4.1. Geomorphology

The geomorphology of Bali Province varies from very low vulnerability to very high vulnerability (Figure 5a). The classification on the basis of geomorphological characteristics indicates that 152.34 km (24%) of the coastline is very highly vulnerable, another 24% (155.28 km) is highly vulnerable, 228.09 km (36%) is moderately vulnerable, 25.68 km (4%) is at a low level of vulnerability, and 77.87 km (12%) has a very low vulnerability rank. The most dominant geomorphology of Bali Province is alluvial plains and sandy beaches.

The coastline of Bali Province is segregated into three distinct types: a rocky and cliff coast in the southern part of Bali Island and Nusa Penida Island; a vegetation coast in the center and west of Bali Island and northern Lembongan Island; and sandy coast for the remainder, which is located in the most coastal area of Bali Province. Coastal areas with rocky and cliff shorelines are highly resistant to erosion [44], as is the case on the southern coast, indicated with blue color, and are therefore less vulnerable, and classified with a very low rank. Cliff coastlines typically occur along mountainous and hilly coasts, where the offshore slope is steep, and little sediment is available for coastal progradation [121].

Sandy beaches, mangroves, and wetlands are found in areas with high and very high vulnerability. Coastal vegetated wetlands and mangrove communities are sensitive to climate change and long-term sea-level change, as their location is intimately linked to sea level. These areas are the least resistant, and are therefore extremely vulnerable to the effects of erosion and sea-level rise [44,74].

4.2. Shoreline Change

According to the results of this analysis, each regency has experienced both erosion and accretion (Figure 5b). Over 630 km of coastline was evaluated, with 163.56 km (26%) of coastline having a very high risk level of vulnerability, while 136.74 km (21%) is highly vulnerable, 56.62 km (9%) is moderately vulnerable, 100.82 km (16%) is at the low level of vulnerability, and 181.50 km (28%) has very low vulnerability. The greatest erosion changes were observed in the Klungkung, Badung, Jembrana, Buleleng and Karangasem Regencies. The study revealed that about 300 km of coastline in Bali Province experienced erosion at various rates due to high waves, coastal sediment mining, construction of beach arching infrastructure, and construction of inappropriate beach structures [122]. Coastlines subjected to erosion are considered to be more vulnerable because of the loss of coastal habitats, infrastructure, and land.

The shoreline changes between 2015 and 2019 are significant due to sediment transport, erosion, and accretion, as well as the impact of reclamation. The impact of reclamation can be noted in Benoa Bay. Many land-reclamation projects are also found in coastal cities that are short of space for expansion, particularly for port activities, such as the expansion of I Gusti Ngurah Rai International Airport. The land reclaimed will be used for the extension of the runaway, as well as additional aprons and terminals. Land expansion also took place at Benoa Harbour and Serangan Island. Land reclamation causes significant negative impacts on coastal habitats and the ecosystem services they provide [123,124]. Changes in wave regimes may affect the stability of the sandy shoreline. Large-scale reclamation can affect the wave regime (wave reflection and diffraction), and changes in current patterns, causing coastal erosion in the coastal area and its surroundings, as happened in Serangan Island [89].

The most significant factor causing the coastal erosion along most of the coastline in Bali Province is large waves and storms, which result in the occurrence of high sand piles as well as damage to coastal infrastructure. In the eastern part of Bali, strong currents and waves on the northeast coast of Bali cause significant coastal erosion processes. The waves are influenced by the level of distance and wind in the fetch. The fetch is the area in which ocean waves are generated by the wind. The coastal fetch in Buleleng Regency predominantly comes from the northeast and the west. The waves occurring from the east and west, depending on the wind season, occurring on Buleleng’s coast cause the movement of alongshore current processes. These processes allow the water to transport sediment along the coast.

Shoreline changes are greatly influenced by processes occurring in the surrounding beach area (nearshore processes), where the beach always adapts to various conditions that occur [125]. This complex process is influenced by three factors, namely the combination of waves and currents, sediment transport, and beach configuration, which mutually affect each other. The changes in these factors vary spatially and temporally, lasting a long time.

4.3. Elevation

Coastal elevation vulnerability rank based on the spatial distribution of coastal elevation is shown in Figure 5c, indicating that a length of 374.14 km (59%) of the 640 km of coast of Bali Province has very high vulnerability, indicated with red color, with an elevation of less than 5 m. Meanwhile, 136.74 km (21%) is highly vulnerable, with a range of elevation between 5 and 10 m. The elevation within 1 m of the shoreline faces the highest probability of permanent inundation, while the coastal strip within 5 m above the normal tides of the present shoreline is also at high risk from severe storm surges.

In accordance with the results of the analysis in Figure 5c, the southern part of the Bali Province area, which is known to possess high cliffs, has the highest coastal elevation, at more than 30 m, and is classified as having very low vulnerability. On the other hand, most of the coastline in Bali Province has the lowest elevation, and is classified as having very high vulnerability.

Using geomorphology and elevation as indicators of sea level vulnerability, the impacts of shoreline change and flood hazards are likely to have a great effect [126]. Gornitz [30] stated that the hazard decreases progressively with higher average elevation. Of the coastline, 85.74 km (13% of coastline) comprises a very low vulnerability zone, indicated with blue color, corresponding to the cliff coast, which has an elevation of more than 30 m and is considered to be resistant to inundation from the sea-level rise or storm surge. The area with low vulnerability, indicated with green color, is also a coastal zone with moderate cliff coast. Of the coastline, 19.62 km (3%) has low vulnerability and 56.62 km (9%) has moderate vulnerability.

4.4. Sea Level

Based on the range value of the sea-level change rate of Bali Province and the scoring table, Bali’s coastline is subjected to very high relative sea-level rise rates and can be considered a highly vulnerable area due to the potential inundation of coastal land (Figure 5d). The range value of mean sea-level variation in Bali Province for the period May 1992–June 2019 is 4.5 to 5.2 mm/year, which is higher than the global average.

Sea-level rise is a relatively weak-amplitude signal, with changes in the range of mm to cm over decades. For this reason, tide gauge datum stability is one of the key observation considerations. The rate of sea-level rise at any given location also varies. According to [127], sea-level changes can be due to two phenomena—one is global, while the other is a regional contribution from land and ocean processes. The variable of relative sea level change naturally includes land subsidence and tectonic motion [128]. Higher sea levels mean deadly and destructive storm surges, which also means more frequent nuisance flooding and also accelerated coastal erosion as sea level rises [13,129,130].

4.5. Tidal Range

Based on the tide data obtained from the four permanent tide stations, the mean tidal range in Bali Province is between 2.133 and 3.265 m. The mean tidal range in Celukan Bawang, Singaraja (northern part) is 2.133 m while in Port of Benoa, Denpasar (southern part) is 3.189 m. The highest mean tidal range is located in Pengambengan, Jembrana (western part), at 3.265 m. The tidal range on Bali’s coast is classified as mesotidal. Remarkably, it is claimed that vulnerability evaluations for the tidal range rank at a moderate vulnerability level, which is indicated with yellow color (Figure 5e). The tidal range, which is incorporated in the CVI as an indicator of the zone impact, with microtidal coastlines ranked as less vulnerable than macrotidal coastlines.

On the basis of the measured mean tide level, it was found that the tide cycle on Bali’s coast is mixed semi-diurnal. A mixed semidiurnal tide cycle is a cycle with two high and low tides of different sizes every lunar day, with a maximum range of approximately 2.24 m. As an archipelagic country, Indonesia can be regarded as an example of tidal complexity. Wyrtki [131] constructed an overview of the tidal pattern occurring in Indonesian seas (Southeast Asian water), which is generally a mixed type that tends to be semidiurnal, being influenced by the Indian Ocean. In most of the Java Sea, the tidal cycle type that occurs is a mixed type (predominantly diurnal). Meanwhile, the semidiurnal type generally occurs in the Malacca Strait and the diurnal type occurs only in a few places. The tides throughout the eastern archipelago, as well as in the adjoining Pacific, are mixed and predominantly semi-diurnal.

4.6. Significant Wave Height

The present study revealed that the mean SWH ranges between 0.55 and 1.86 m. Most of the northern coast is in the very low vulnerability rank, while the southern coast is in the very high vulnerability rank (Figure 5f). This is due to the northern coast being a closed ocean (Bali Sea), in contrast with the southern coast, which directly faces the open ocean. The largest waves occur where there are large expanses of open ocean (Indian Ocean) that can be affected by wind. Moreover, the waters of southern Bali are close to existing passages of the Indonesian Throughflow (ITF) into the Indian Ocean. The waves are primarily caused by the wind blowing over the surface of the sea and transferring energy from the air to the water, forming waves [132,133]. Wave heights depend on the characteristics of the wind that produces them.

The significant wave height can be used instead of the wave energy, which indicates the influence of waves on coastal erosion. Stronger wave energy corresponds to more intense and accelerated coastal erosion [52,134]. Hence, the distribution of wave energy along coastlines is an important control of spatial variability of coastal erosion [135,136].

4.7. Coastal Vulnerability Index of Bali Province

Over 640 km of shoreline is evaluated and ranked along the Bali Province. The classes of CVI values are divided into “low vulnerability” (green), “moderate vulnerability” (yellow), “high vulnerability” (orange) and “very high vulnerability” (red) categories, based on the quartile ranges and visual inspection of the data. The calculated CVI values range from 3.54 to 39.53. CVI values below 8.9 are assigned to low vulnerability. Values from 8.94 to 12.25 are considered moderate vulnerability. High vulnerability corresponds to values between 12.25 and 17.68. CVI values above 17.68 are classified as very high vulnerability. About 138 km (22%) of the mapped shoreline is classified as being at very high vulnerability, and 164 km (26%) of shoreline is in the highly vulnerable category. The remaining shoreline, 168 km (26%) and 169 km (26%) are in the moderate and low vulnerability categories, respectively. Figure 6 shows a histogram of the percentage of Bali Province shoreline in each vulnerability category based on CVI analysis.

Figure 7 shows the overall ranking category distribution for the classes’ CVI values, and also the CVI map of Bali Province. Regions with low vulnerability are formed by variations in geomorphological parameters, with cliff coast and elevations more than 30 m. The rate of shoreline change between accretion and erosion is in the range of −1 to 1 m/year, and is categorized as stable due to the balanced condition of accretion and erosion. Regions with low vulnerability spread and dominate along the coastal regions in the southern of Bali Island (Badung Regency, 22.35 km), West Bali National Park (Buleleng Regency, 50.31 km), the southern part of Nusa Penida (Klungkung Regency, 65.25 km), Jembrana (9.13 km), and Karangasem (21.60 km). These regions are relatively safe from sea-level rise threats, despite their mean significant wave heights indicating very high vulnerability.

Regions with moderately vulnerable zones are formed by variations in the geomorphological parameters with alluvial plains and low cliffs, followed by an elevation of 10–20 m. The rate of shoreline change due to accretion and erosion is in the range of −1.1 to −2 m/year and is the most dominant factor. Regions with moderate vulnerability are located in the eastern part of Buleleng (74.49 km), Jembrana (38.21 km), Karangasem (13.50 km), Klungkung (11.57 km), Gianyar (11.33 km), Tabanan (7.05 km), Badung (5.94 km), and Denpasar City (1.11 km). Regions with moderate vulnerability can also be categorized as safe regions; however, if those regions are not managed appropriately, they could become high or very high vulnerability regions.

Regions with high vulnerability are formed from variations in geomorphological parameters including estuary, sandy beaches, and mangroves, with elevation heights of less than 5 m. The rate of shoreline change is in the range of −1.1 to −2 m/year in some locations at Buleleng Regency and Perancak estuary (Jembrana Regency), as well as more than −2 m/year in some locations at Jembrana Regency and Karangasem Regency. Regions with high vulnerability are scattered throughout moderate vulnerability regions, and are present in the Buleleng (40.56 km), Karangasem (36.57 km), Jembrana (27.75 km), Denpasar (12.96 km), Klungkung (12.77 km), and Tabanan Regency (5.82 km). Regions with high vulnerability require attention to protect them from potential damages or the degradation of coastal areas.

Regions with very high vulnerability zones are formed from variations in geomorphology parameters including sandy beaches, mangroves, and coastal buildings, and have elevations with heights of less than 5 m. The rate of shoreline change is more than −2 m/year, which means that such regions experience coastal erosion. Regions with very high vulnerability are found at Klungkung (32.07 km), Badung (30.49 km), Denpasar (30.07 km), Tabanan (19.99 km), Jembrana (15.29 km), Karangasem (11.76 km), and Gianyar Regency (1.95 km). Regions with very high vulnerability also require more attention, similar to high vulnerability regions, in order to avoid the damage and degradation of coastal areas.

The present study is an attempt to categorize the coast of Bali Province according to its coastal vulnerabilities. The most important parameters for the study area are geomorphology, shoreline change rates, elevation, and significant wave height, since the other variables are constant (classified into the same class). The CVI method is commonly used to assess the effects of sea-level rise on coastal areas. The CVI method is simple, as it uses a ranking system; therefore, it is easy to identify regions with high vulnerability. It is easy for policy and decision makers to make decisions regarding proper management programs for coastal regions with high vulnerability in order to prevent the effects of sea-level rise [24,97,137]. However, the CVI method has several disadvantages as well; for one, it is only based on geological and physical parameters. It does not consider the effects of social/human activities on ecological and physical changes, and a limited number of parameters are used as input for assessing vulnerability. Therefore, further studies in this area need to include the socio-economic aspect, local factors at a detailed scale, etc., as additional inputs to produce better and more accurate results with respect to coastal vulnerability assessment.

This study proves that the quality of the results is highly dependent on data quality. Access to reliable data was one of the most challenging factors throughout this study, since the available free data have low positional accuracy and detail due to the scale of the products that are made available to the public through the website. In these circumstances, it is necessary to apply the primary data from the local related institution as a comparison to make a high-accuracy vulnerability map.

Despite these limitations, the grid cells were used for detailed assessment and validated on 35 points. The fieldwork performed as part of the model validation was conducted along Bali Island beach. However, some areas were prohibited from being accessed due to being private areas or being subject to ongoing construction, and some areas were difficult to visit due to there being no access roads leading to the beach area. Nevertheless, representative locations were selected and visited in order to perform a ground check and validation. By following Equation (4), the agreement between fieldwork and CVI was calculated to be 0.80, indicating good agreement.

The methodological approach is simple, robust, and easy to implement for nation-wide mapping, since it is based on well-defined criteria through an index-normalized formulation, and additional parameters can be included with the determination of different weights. Despite the fact that further improvements in the methodology are required in order to be able to assess coastal risk for sea-level rise mitigation and adaptation measures, the present results are an important contribution to the identification of coastal vulnerability, constituting an additional instrument for decision-makers with responsibilities of management and planning of areas exposed to the sea-level rise.

Despite the coastal vulnerability assessment at the national level using CVI having already been calculated in Indonesia, the method for acquiring CVI was different from that used in the current study. The national assessment did not calculate the rate of the shoreline change based on satellite data or regional changes to the coastline. The risk level of coastal erosion/accretion (shoreline change) was determined on the basis of the results of proxies of coastal geomorphology, coastal slope, and significant wave height. Moreover, the scale values of the shoreline changes were principally determined by weighting the scale values of those parameter proxies. Coastal areas with very high levels of vulnerability are marked by a slight coastal slope, a high level of vulnerability in terms of geomorphology, and high ocean waves. The current study used updated and high-resolution satellite imagery data to obtain the shoreline change, and conducted an analysis based on grid cells to acquire the precise index score. It is important to highlight that for index-based methodologies such as CVI, the availability of reliable and up-to-date databases is crucial [65].

5. Conclusions

The assessment of coastal vulnerability can be performed by CVI with remote sensing and GIS approaches to measure the risk level of coastal areas in Bali Province efficiently and appropriately, in consideration of dynamic processes and both geological and physical parameters. By using a 1 km grid cell and recent data, this study provides an up-to-date coastal vulnerability map of Bali Province that is more precise than in previous studies.

This study indicated that geomorphology, shoreline change rate, coastal elevation and significant wave height are the most contributing parameters determining coastal vulnerability, since the sea level change rate and tidal range were given the same risk level along the coast. The most contributing parameter could be further improved by means of weighted determination.

In addition, the proposed index shows the feasibility of the index for making assessments on coastal vulnerability in other coastal areas when dealing with climate change. The CVI map of Bali Province presented in this study can be used by decision makers and related stakeholders in disaster mitigation and management in order to mitigate the effects of climate change.