Sand dune risk assessment in Sabha region, Libya using Landsat 8, MODIS, and Google Earth Engine images

Abstract Globally, sand dunes are a major environmental problem that causes damage to urban areas, transportation, and population. The current study proposes a comprehensive investigation on sand dune risk modeling in Sabha located in the southwestern part of Libya. Data from various sources were collected and prepared in a GIS database. Data from 2016 were used to derive several controlling factors, such as altitude, rainfall, soil texture, wind direction and speed, land cover, and population density. Next, sand dune susceptibility, hazard and vulnerability assessments were performed. Finally, a risk map was produced. Results indicate that land use and soil are the most influential factors affecting the sand dunes in the study area, whereas rainfall is the least significant factor. Results indicate that, southern part has a higher chance of sand dune occurrence than the northern part, whereas the highest risk zone is located in the middle part, where the urban and agricultural lands are present. More than 200 km2 of the study area are under high and very high risk zones. Overall, this study provides an effective tool for assessing sand dune risk in Sabha, which can be useful for land management.


Introduction
Sand movement is one of the significant causes of environmental problems in dry and semi-dry areas worldwide (Iranmanesh et al. 2016). Nearly 77% (42,779 km 2 ) of the total drylands of the globe (55,558 km 2 ) experience a dry or semi-arid climate, and sand dunes generated by wind activities (Darraz 1995) affect a large proportion of that area. Sand dunes cover vast areas of our planet, which sometimes poses threat to vibrant cities, villages, road networks, agricultural lands, water resources, and human health (Walker and Cronon 2009). Despite significant efforts to reduce sand dune movement, recent research indicates that dune movements are continuing and its risk to the environment and humans is becoming increasingly serious (Lamqadem et al. 2018;Masoudi et al. 2018). The advancement of sand dunes in nearly all of the North African countries creates severe environmental problems, including desertification (Schuster et al. 2006). Sand dunes creep up agricultural lands and urban environments, as well as economic and social activities in the area (Darraz 1995). Particularly, Libya is largely affected by sand dunes because a large part of the southwest of Libya has been covered with sand, which is approximated to be over 250,000 ha (Koja 2015). Consequently, many cities in Libya, including its capital city Tripoli, are exposed to sand movements. In such cities, sand movements have not only affected people but also have threatened road networks and other transportation systems.
Two main factors generally contribute to the moving sand in African states: humans and the environment (Darkoh 1998). The first factor (i.e., humans) that affects the sand dunes in Libya is due to the exploitation of soil and removal of vegetation (Darraz 1995). The second factor includes the effects of climate change, reductions and fluctuations of precipitation amounts, and severe winds. These factors increase the proneness of dry and semi-arid areas to sand movements. The old methods of assessing sand dunes in African states include simple estimations of the areas that are covered by sand and the potential agricultural and urban areas under risk of sand dunes. The estimation in 1913 demonstrated that the sand dune movement within Libya ranged between 15 and 20 m/yr (Darraz 1995;Darkoh 1998). The Ministry of Planning suggested implementing necessary plans to afforest the areas with high sand flows. Given the rapid and increasing growths of population and urbanization in the dry and semi-arid areas, the potential of sand dunes has increased as well; reforestation may be inefficient at this point (Tesfay and Tafere 2004). Although reforestation could have been an effective solution in the past decades, current remote sensing and GIS solutions provide advanced tools for better planning and mitigation of sand dunes.
Remote sensing provides satellite or aerial imagery that can be used to generate numerous forms of spatial information, such as transportation networks, urban areas, land use types, and other information regarding utilities and infrastructures in a particular region (Althuwaynee et al. 2012(Althuwaynee et al. , 2014aBui et al. 2017;Chen et al. 2017;Nampak et al. 2018;Mezaal and Pradhan 2018). By contrast, GIS provides efficient and useful solutions to create, store, analyze, and visualize those data acquired using remote sensing (Pullar and Springer 2000). Remote sensing and GIS provide solutions to the various stages of dune movement and risk analysis, including susceptibility mapping, hazard assessment, and risk modeling (Ghadiry et al. 2012;Regmi et al. 2014).
First, sand dune susceptibility is theoretically the spatial probability of the occurrence of sand movement in an area at a given period (Hugenholtz et al. 2012). This analysis can be accomplished by providing certain controlling factors that affect the occurrence of sand movement, such as wind speed, wind direction, altitude, land use, and vegetation density (Darraz 1995). These factors are then overlaid considering their importance based on the weights given to each factor by experts or by analyzing historical sand dune events. The final output of this analysis is a map that shows the areas that are prone to sand dune with the degree of the proneness, which is normally categorized as follows: very low, low, moderate, high, and very high (Selvam et al. 2015). Second, sand dune hazard defines the possibility of sand dunes occurring throughout a given area. An ideal sand dune hazard map presents not only the chances that a sand movement may form at a particular place but also the chance that it may travel from one place to another. Sand dune hazard maps are often generated by multiplying the sand dune probability, which is estimated in susceptibility analysis, and a triggering factor or a combination of many triggering factors, such as the lack of precipitation and severe wind activity. Third, risk modelling of sand dunes aims to define the element at risk of sand movements. This stage often includes vulnerability assessment, which concerns the people and land use under risk of sand dunes. The risk is defined as the vulnerability multiplied by the hazard. The final output of risk is also categorized into five classes. This process ensures the effectiveness of the information presented, which decision-makers can easily use in their planning.
Local administrations in Libya have scarcely paid attention toward reforestation. Therefore, the risk of sand dunes can be high. Subsequently, a need exists to identify areas at risk of sand dune movements and seek the potentially affected land cover types and population. Thus, this study presents a comprehensive risk analysis of sand dunes over a part of Libya, including susceptibility mapping, hazard and vulnerability assessments, and risk modeling via an expert-based method and geospatial information. The main goal of this paper is to quantify the risk resulting from the sand dunes that migrate from the desert toward the city and assess the effect of sand cover in the agricultural and urban land areas. Therefore, the susceptibility and risk modeling has been performed because of dune migration. The key importance of the analysis is to provide a solution for sand dune risk assessment and replace the traditional methods of assessment, which are not as effective as expected.
Sand dunes cause several types of effects on buildings and transportation. The source of these effects is mostly due to the continuous and permanent wind movement. The effects cause buildings to be covered by sand in a short time. The effects in areas that suffer from droughts are high for an extended period of time. Particularly, Libya is to sand dunes, which has affected its buildings and transportation systems. Sand dunes have also affected the roadways, impeding traffic, and road movement through the continuous movement of sand. Moreover, the effects of sand dunes affect the planning processes. Planning operations require significant efforts and a huge cost to address the phenomenon of sand dunes.

Previous works
Studies on sand dune risk modeling are limited. However, a few existing related studies are essential to understand how sand dune risk can be quantified. This section reviews some of the relevant studies. Early research on sand dunes has shown systematic estimations of wind speeds and sand flux (Heywood 1941;Livingstone 1988). Several other researchers have focused on the standard measurement of wind speed and sand flux using one-of-a-kind techniques, including rotating-cup anemometers, paper flags, and sand traps (Livingstone et al. 2007). Studying sand dunes in African countries started at the end of the 19 th century with the innovative work of (Heywood 1941), while study of Embabi (1995) has evolved as the maximum comprehensive research. The details of sand dune migration, formation and geomorphology can be found in various sand dune ecosystem study (Breckle et al. 2009).
El Gammal and El Gammal (2010) studied the hazard impact and the development of sand dunes in Egypt. The author investigated the sand dune hazard using medium-resolution satellite images (Landsat ETMþ) and field works. Particularly, the researcher studied the shape of sands and their origin, as well as their movements in the study area. Ghadiry et al. (2012) developed a GIS-based model to automate the extraction of sand dunes from SPOT satellite images. They also studied the threats of sand dunes to roads, irrigation networks, water resources, and urban areas. This study is also helpful for urban planning using remote sensing techniques (Salman Aal-Shamkhi et al. 2017). The model was developed in ArcGIS model builder and python scripting. The model took images as input and automatically extracted the sand dunes. The results showed that the rate of sand dune movement ranged between 3 and 9 m/yr. The majority of sand dunes had a rate movement between 0 and 6 m, and extremely few sand dunes had a movement rate between 6 and 9 m. Overall, they found that integrating remote sensing and GIS provided the necessary information to determine the minimum, maximum, mean, rate, and area of sand dune migration. Kaboodvandpour et al. (2015) developed a model to predict sand dunes in Iran using adaptive neuro-fuzzy inference system (ANFIS). They compared this model with another model to predict dust storm occurrences. Accordingly, they used artificial neural network (ANN) and multiple linear regression (MLR) models, which were based on meteorological variables. Their results indicated that the accuracy of the ANFIS model was higher compared with those of MLR and ANN.

Study area
The selected site for the current study is Sabha, which is a city in southwestern of Libya ( Figure 1). The study area is 640 km south of Tripoli and borders the following districts: Wadi Al Shatii in the north, Al Jufrah in the east, Murzuq in the south, and Wadi Al Hayaa in the west. Large sand dunes are resulted due to the desertification process. These sand dunes can be seen at the border region of Algeria up to the Sabha border. Urban settlements and large oases can be found in the desert. The desert near to the Sabha city is popularly known as Awbari sand sea. The red desert near to Sabha is a mixture of terrains, mountains, small valleys, deserts and sand dune in waveforms. The desert is very popular among the locals and it crosses from Daraj to Adiri via Awaynat Wanin. There is a continuation of the Awbari Sand Sea toward the very famous Fezzan lakes. The large civilization at Sabha, Awaynat and Germa have been developed along this lake, which is the home to the ancient people. The climate is a desert in most parts of the district. Dust storms lasting for 4-8 days are common during spring. The surface temperature of the area on average ranges from 20 to 40 C with the lowest and highest temperature in November and May, respectively. The region receives an annual rainfall of 64 mm. No perennial rivers exist in the area; however, the region is abundant with groundwater aquifers. The population density of the area is 1.86 people per sq. km. Sabha is served by an airport and a 800-km railway property. In addition, many irrigation canals exist in Sabha, which are used to provide fresh water for growing crops.

Datasets
Risk assessment of sand dunes using GIS requires several types of spatial information, such as land use, topographic data, climate, and demography. In the current study, several datasets were acquired from various sources. All the data used in this study are open access and freely available. The data used in this research are described in  Table 1. A set of Landsat images was collected from the USGS database to prepare a land cover map of the study area. The Landsat images were checked according to the image quality provided by the data provider. These images were then filtered, and only the best (with less noise and less cloud cover) was selected. A DEM of 30 m spatial resolution was acquired from ASTER database and prepared in GIS, which involved several steps, such as spatial clipping, filtering, and gap filling. In addition, rainfall data were acquired from MODIS databases over a 1-year period (2016), and the average rainfall map was then created in GIS and used for the subsequent analysis. Wind speed and direction are critical information for sand dune modeling; thus, they were collected from online sources and converted into GIS data format. Finally, the thematic layers of wind speed and direction were prepared for the analysis and processing. Population data are essential information required nearly in any risk modeling. Therefore, the current study prepared a population density information for the Sabha region at 250 m spatial resolution in 2016. The population data were provided by the Sabha Municipal Corporation. Finally, a soil map was collected from the USGS website because no local geology map was found. Once these datasets were acquired, they were prepared and assembled in a GIS database and transformed into one unique coordinate system (i.e., WGS 1984). Figure 2 presents the overall workflow for the production of the sand dune risk map for the Sabha region using remote sensing and GIS techniques. First, the input data, which contained six main datasets, namely, Landsat OLI image, DEM, soil map, rainfall data, wind information (speed and direction), and population data, were prepared. These data were then used to derive several controlling factors (i.e., altitude, rainfall, soil texture, wind direction and speed, land cover, and population density) that affect sand dune movements in the study area. In this research, the AHP approach is adopted for risk assessment. This is because AHP performs better than other techniques due to the use of geometric mean that projects the exact value (Pourghasemi et al. 2016). Processing includes four main phases: AHP analysis, hazard analysis, vulnerability assessment, and risk modeling. In the AHP, opinions from five experts were utilized to assign weights to the criteria, namely, altitude, rainfall, soil texture, wind direction and speed, land cover, and population density. These experts are specialists in remote sensing and environmental studies, with familiarity on the study area. On the basis of expert opinions, the AHP method was utilized to estimate the weights of the factors. Subsequently, these factors were overlaid to produce the susceptibility map using the estimated weights. Second, hazard analysis was performed by multiplying sand dune susceptibility information and wind speed. Wind speed is considered the triggering factor for the sand dune in the study area. Vulnerability assessment phase includes the preparation of land cover and population factors, reclassification of each of them into several classes, and assignment of weights to each of the sub-classes of each factor. Finally, the sand dune risk map was produced by multiplying the vulnerability and hazard maps. The final risk map was then reclassified into five categorical classes: very low, low, moderate, high, and very high.

Sand dune controlling factors
Several existing factors affect sand dunes in the study area, including altitude, land cover type, soil texture, rainfall, population, wind speed and direction. Each of these factors are explained in the following subsections.

Altitude
Altitude is one of the most common topographical factors used for many environmental studies, like sand dunes (Tao et al. 2017). Topography is a major factor in activating erosion processes. Consequently, the geological structure of the area may have a clear influence in increasing the amount of sediment. In addition, a clear relationship exists between the altitude and surface temperature, rainfall, and living comfort. Thus, this factor has been considered in many studies and sand dune risk assessments. Figure 3(a) shows the altitude distributions in the study area derived from ASTER DEM at 30 m spatial resolution. The altitude in the study area ranges from 370 to 514 m above mean sea level. The northern part of the region has lower altitudes than the southern part. This difference in topography between the north and south parts of the study area may affect the rates of sediment movement accruing into both parts.

Land cover
Land cover is another important factor used in many applications, e.g., sand dune modeling (Santana-Cordero et al. 2016). It reflects human activities, as well as environmental development or degradation. For example, the development of farming can have positive effects on sand dune activities in a study area. However, moving from forest to bare soil can degrade the local environment of the region and thus lead to sand dunes. In the current study, the land cover map of the survey area was derived from a Landsat OLI image (2016) via object-based image analysis method with 87% accuracy. The satellite image was initially radiometrically corrected by converting the digital numbers into important reflectance values. This step allows the satellite image to reflect approximately the exact ground reflectance, thereby improving the image processing interpretation of the derived products. The image was then geometrically corrected by bringing the pixel values into exact ground locations and properly projected in the local coordinate system of the study area. Once the pre-processing steps were applied, the calibrated image was then used in OBIA to derive the land cover feature found in the scene. First, the image was segmented into homogenous non-overlapping groups using the multiresolution segmentation approach. The parameters of the algorithm were set as 136, 0.6, and 0.3, for scale, shape, and compactness, respectively. These parameters were selected via a trial-and-error method. Second, by using a supervised classification method called support vector machine (SVM), the image objects were classified into four classes: agriculture, bare land, sand, and urban area (Figure 3(b)). The parameters of the SVM algorithm were empirically selected. The penalty parameter# was set to 100, and a linear kernel function was utilized to map the feature into a high dimension. The final land cover map was evaluated via the traditional accuracy assessment approach, overall accuracy, and Kappa statistics. Overall, the SVM model could classify the image objects into the four classes with an accuracy of 87% and a Kappa statistic of 0.85. Land cover is also an important factor for vulnerability and risk assessments of sand dune. Hazard maps are often overlaid with such data to produce maps that demonstrate the effects of potential sand dune events on population and properties.

Soil texture
Rock erosion and fragmentation often result in the formation of several soil types (Blume et al. 2016). Sediments move after the separation from the original rocks and may continue indefinitely. Alternatively, deep sediments are occasionally formed after the movement, like sand dune. Soil type is a major factor that affects plant growth as a response to erosions. In addition, suitable soil wetness and the low water that retains the capability of sand are essential components influencing plant growth. Moreover, vegetation interacts with different soil properties. Thus, soil texture was considered as a factor in dune susceptibility assessment in this work. Figure 3(c) shows the soil map of the study area. The dominant soil types in the study area include Adrar Al Faw (northern part of the area), flat desertic areas (some parts in the middle of the region), and Sarir Al Qattusah (most of the area, especially in the southern and middle parts).

Rainfall
The density and growth of annual plants depend on the amount of rainfall in the study area. The area is covered by annual plants and sparse shrubs because of the relatively low annual precipitation. Figure 3(d) shows the average rainfall distribution over the study area and the dots corresponds to the survey points in 2016. Nevertheless, a typical single year data represent the historical data over a multi-year period. Therefore, typical single year data are the average of historical data over a long period. It is reliable, save time and does not make the process complicated. The amount of rainfall ranged from 0 mm/month to 24 mm/month. The map demonstrates that the survey area faced heavier rainfall in the middle and south parts of the region than the north part.

Population
Many cities worldwide face challenges in a rapidly increasing urbanization and population growth, including the study area. Numerous effects occur for this phenomenon, e.g., human activities, developing farming due to the need for food, and deforestation due to the need for urban areas, where people transform forest and vegetated areas into bare soil or towns. Occasionally, the location of the region within the dry lands implies that erosion activities possibly occur due to this geographical characteristic. Therefore, the area can be utilized to investigate the relative influence of natural processes and human activities. Population data were collected from Diva GIS (http:// www.diva-gis.org/gdata). Diva GIS is the database, which provides free and effective spatial data at varying level such as country level data: administrative boundaries, roads, railroads, altitude, land cover, population density (Hijmans 2009). Agriculture is the main economic activity in the study area. Consequently, people mostly convert forests into farm lands and grasslands. In addition, population information is necessary for sand dune risk assessment, where the hazard map intersects with the vulnerability map (land cover and people). Figure 3(e) shows the population density of the study area at 250 m spatial resolution. The lowest number of individuals per 250 m cell size is 189, which is possibly in the northeast part of the region. Conversely, the highest number of individuals per a raster cell is nearly 3000 people. The maximum population density is found in the middle part of the region.

Wind speed and direction
Wind is a key factor of erosion (Iversen and Rasmussen 1999;Borrelli et al. 2014). A significant portion of the study area is affected by wind activity, and severe winds most likely increase the potential of a sand dune in the study area. A gradual loss of wind speed results in the deposition of sediments and other materials, which create geomorphologic phenomena; the most significance of these phenomena arising from the wind are sand dunes. Wind movement often occurs at different heights, thereby creating a unique form of sediment transport. The transfer by wind deposits is directly linked to wind speed, where sand grains and occasional particles of gravel move through a creep, whereas particles of mud, silt, and dust move by air suspension. This phenomenon occurs particularly when the wind speed is high. Wind speed rates can directly contribute to the growth of new sand dunes because they can carry particles and remain constant for extended periods of time (e.g., in a full day or for several days). Similarly, wind direction has a major role in the direction of sand movement, particularly in the formation of different types of dune. Figure 4(a,b) show the average seasonal wind speed and wind direction over the study area in 2016, respectively. The wind speed ranged from 10 to 17 m/s. The highest wind speeds are found in the southern part of the region. By contrast, the low-speed winds are located in the northern part of the study area. The southern part of the area is nearly a desert; thus, the high speed of winds can increase the probability of sand dune occurrence. In addition, the wind direction map of the study area demonstrates that the wind moves from the southern part toward the northern part of the region.

Sand dune susceptibility mapping using AHP method
Several methods generally exist to produce sand dune susceptibility maps, such as statistical, data mining, and expert-based approaches (Tyoda 2013). The first two approaches often require an adequate number of events or past records. Meanwhile, expert-based methods do not require inventory information; instead, the weights of each controlling factor are determined by expert opinions and the AHP method (Yang et al. 2016). The current study had no available information regarding ancient sand dunes. Therefore, an expert-based method using AHP was utilized. The AHPbased sand dune susceptibility mapping comprises three main steps: gathering of expert opinions, weight determination, and factor overlay. The first phase includes preparation of a questionnaire document and subsequent submission to five experts who are specialists in environmental studies and with familiarity of the study area.

Sand dune risk assessment
The sand dune risk map was generated by three main steps: sand dune hazard mapping, vulnerability assessment, and risk modeling. Hazard maps can be produced by multiplying the sand dune susceptibility map with a triggering factor. In this study, the triggering factor is the wind speed because nearly all the sand dunes in the study area occurred as a result of high wind speed and wind fluctuation (Cooper et al. 2004). The sand dune hazard map was then converted into a probability using a rescale function tool in GIS. The hazard map shows the areas that are prone to sand dunes using wind speed as a temporal component. The second step involves the preparation of the sand dune vulnerability map of the study area. In this step, the land cover and population data were used because in a vulnerability assessment, one attempts to evaluate the elements at risk, which in the case of a sand dune is the land cover (properties) and population density. However, to produce the vulnerability map using these two factors, expert opinions are also needed to assign weights for each class of land cover and population. The land cover of the study area comprises four types: agriculture, urban areas, sand, and bare land. The weights provided by the experts are as follows: sand ¼ 0, bare land ¼ 3, agriculture ¼ 6, and urban ¼ 10. By contrast, the population density was converted into a probability ranging from 0 (less dense) to 1 (highly dense). Subsequently, the reclassified land cover with the assigned weight values and the population probability maps were overlaid, thereby creating the vulnerability map of the study area. Finally, the hazard and vulnerability probabilities were multiplied to produce the sand dune risk map. The risk map was then produced as a probability and as a classified map. The risk classes are very low, low, moderate, high, and very high according to the natural break method implemented in GIS.

Results and discussion
4.1. Results of sand dune susceptibility Table 2 presents the weights of the conditioning factors for each expert and the geometric mean of those weights. The experts provided the importance degree in the pairwise table for the factors, and the AHP method was used to calculate the weights in the table. The main reason for choosing geometric mean as opposed to arithmetic mean is the only aggregation form that makes a valid and approximate value by keeping the reciprocity. Geometric mean holds the axiom while the arithmetic mean do not. However, the geometric mean is convenient and easy to implement in AHP method (Althuwaynee et al. 2014a). The geometric mean is exactly a concave aggregator that provides good performance against arithmetic means. By analyzing the similarities and differences among the expert opinions, the correlation analysis confirms that the lowest R 2 is 0.29 between the second and third experts. By contrast, others opinions are nearly similar to R 2 >0.90. According to the estimated geometric means, land use, and soil seem to be the most influential factors that affect the sand dunes in the study area, whereas rainfall is the least significant factor.
To produce the sand dune susceptibility map of the study area, the geometric means (Table 2) were used to overlay the controlling factors in GIS using the weighted sum function. The susceptibility map was a produced probability (0-1) and was subsequently classified into five certain classes: very low, low, moderate, high, and very high. The latter simplifies the susceptibility map and establishes an easily readable material for decision-making and land use planning. Figure 5(a) presents the sand dune susceptibility as a probability map. In the figure, colors green and red denote low and high probabilities, respectively. The southern part of the study area is prone to sand dunes, whereas the northern part is less prone to sand dunes. The sand dune probability map was reclassified into five categorical classes (i.e., very low,  low, moderate, high, and very high) using the natural break algorithm of ArcGIS ( Figure 5(b)). This map permits easier interpretation and is suitable for decisionmaking. The study area was classified into different zones according to the probability of sand dune occurrence. According to the map (see Figure 11), nearly 25% of the study area located in the south part of the region are at a very high risk zone. By contrast, the central and northern parts of the area are supposed to be at moderate to very low risk zone.

Results of sand dune hazard
Unlike the susceptibility map, the hazard analysis of sand dunes aims to create a map that shows the spatial and temporal (regarding triggering factor) components. Therefore, the hazard map is an indicator when a sand dune may occur in the study area. The temporal component is provided by the wind speed. Sand dune occurrences are expected when the wind rate increases. Figure 6 presents the hazard probabilities over the study area. The central south areas are most likely to experience sand dunes. Generally, the southern part of the region has a higher chance of sand dune occurrence than the northern part.

Results of sand dune risk assessment
Sand dune risk map is based on two major steps: sand dune hazard mapping and vulnerability assessment. The hazard map was prepared as a probability map, the pixel values of which ranged from 0 to 1, indicating low and high hazard probability, respectively. Meanwhile, the vulnerability map was produced by combining land cover and population density maps using expert-based weighting procedure. Regarding land cover, the map was reclassified, and the experts provided each class a specific weight. The sand and bare land were given low weights (0 and 3, respectively), whereas agriculture and urban classes were assigned high weights (6 and 10, respectively). The reclassified land cover map is shown in Figure 7(a). The central part of the area mostly contains urban and agriculture lands, which are expected to be at a high risk zone. In addition, the second important input data required for the vulnerability map is the population density map. The population density was rescaled in the range of 0-1 to be a probability population map. The lowest population is found in the northeast part of the study area. The sand dune vulnerability map was produced once the land cover map is reclassified according to the expert-given weights and the population density was converted into a probability map. The vulnerability map was produced by overlaying the two maps using the raster calculator in GIS. The map indicates that the middle part of the region is mostly prone to sand dunes. Some of northeast part of the region is also at risk of sand dunes. The sand dune vulnerability map was also produced using a probability map, which is needed in the risk analysis. After preparing the sand dune hazard and vulnerability maps, the sand dune risk map was produced by overlaying the two maps in GIS. The highest risks are located in the middle part, where urban and agricultural lands are present. The sand areas are found to have very low risk of sand dunes.
To simplify and make the risk map interpretable and easy for decision-making, it was reclassified into five classes: very low, low, moderate, high, and very high, using the natural break method in GIS. The analysis shows that the highest risk areas are in the middle of the south part of the region. The land cover map of the area indicates that these high-risk areas are agricultural lands. Famers have their method of reducing sand dunes risks (e.g., cultivation of spiny and semi-woody plants, which are characterized by rapid growth and short life span, such as castor bean and cane). However, risks may exist especially when the wind speed is high. Figure 7(b) shows the sand dune risk map overlaid with the average wind speed and direction. From the figure, the high-risk areas are expected to have high-speed winds, and the wind direction is most probably from the south toward the north part of the region crossing the high-risk areas.

Results of sand dune impacts on urban areas and population
Sand dunes have significant impacts on urban areas, agricultural lands, and the population in a particular region. In this study, the effects of high-risk sand dunes were analyzed in GIS. Figure 8(a) shows the urban areas within the study area. The urban areas are mostly located in the middle part of the region. The urban areas were converted into points using the centroid of the polygons to calculate the urban density within the study area. Statistics on the affected cities could be calculated by intersecting the urban density and sand dune risk maps. Figure 8(a) shows the minimum, maximum, and sum of urban densities in each risk class. The risk map contains five classes with various areas. Table 3 lists the results of affected urban densities in the risk map. From the table, >200 km 2 of the study area are under high and very high risks, whereas more than half of the region is at very low risk. The results also show that the highest urban densities are in the high-risk class. The urban densities in the moderate, low, and very low classes are greater than the very high-risk class.
The risk map was intersected with the population density to calculate the affected people by the potential sand dunes in the study area. Figure 8(b) shows that most  population are in the very low class of the risk map. However, a substantial number of individuals who are at high and very high risk zones exist. Table 4 shows that >13,709 people are at high risk zones and 5397 are at very high risk zones.

Validation
Image composition and classification tools of Google Earth Engine were used in determining available images to perform a visual comparison and discover the changes. Figure 9 shows the cultivated crop lands and the Sabha City and its surrounding areas. The changes in Sabha could be understood from various images, as shown in Figure 10. MODIS 16-day NDVI data from January 2015, 2016, and 2017 were used to compare the spatial and temporal changes. Figure 10(a-c) reveals the continuous growth of newly cultivated agriculture fields. Agriculture lands are mostly at the southern part of Sabha City, which is at high risk according to the sand dune risk map. MODIS Aqua Surface Reflectance Global 500 m data were used to understand the changes to the surface of the Sabha region and its surroundings. In Figure 10(d-f), light yellow indicates only sand without water  or moisture, whereas deep brown denotes traces of water in the Sabha City and agriculture lands. Therefore, the moisture or the amount of water is less in Figure 10(f) than those in Figure 10(d,e), which implies the lack of rain in January 2018. Changes can be determined in the desert, which is at northwestern part of the three images. Google Earth imagery was also used to validate the risk map. The three images in Figure 10(g-i) present the changes. The Google Earth images show the significant changes in shape and size of the agriculture lands, as well as the changes on the basis of increasing number of agricultural fields. More agricultural fields are found in January 2016 than in 2014 and 2015. In addition, the images show that the Sabha City is increasing in size. Therefore, the agriculture lands and the city are growing toward the northeastern part because the southwestern part of the Sabha City and the agriculture lands are affected by sand dune migration. MYD09Q1.005 Surface Reflectance 8-Day L3 Global 250 m data were used to identify the changes in the desert. Figure 10(j-l) shows the changes of rough surfaces. The changing surfaces imply the migration of sand dunes. Generally, sand dune migration is toward the southwestern part of the city and the agriculture lands. Hence, some parts of the agriculture lands near the desert are affected by the sand dunes, whereas the central part is safe. However, sand dune migration depends on wind direction and speed, which is higher at the southwestern part of Sabha than the northern part. The value of this comparison is evident. However, some problems still exist in the surface reflectance data. For validation, the resultant sand dune risk map of Sabha, Libya was compared with the updated actual features in Google Earth based on the areas affected by sand dunes. Figure 11(a) shows the south western part of the Sabha city, which is genuinely affected by the sand dunes. Figure 11(b,c) shows that the agriculture lands are affected by sand dunes. In addition, the agriculture lands are at an elevation of >150 m with the highest point at 187 m, whereas the southwest desert is at an elevation of 135 m. Evidently, the urban and agriculture lands act as obstacles and are mostly affected by sand dunes due to high elevation. Figure 11(d) shows the sand dune migration toward the city and agriculture lands.
Therefore, high wind speed from the south part migrates the sand dunes toward the city and agriculture lands. The validation result shows that the prepared risk map Figure 11. Validation of sand dune risk map using Google Earth images.
is relatively comparable and effective. The accuracy of the prepared risk map can be improved by adding parameters and using data of good quality.

Conclusion
In this study, a workflow for sand dune risk modelling in Sabha, Libya is developed using remote sensing and GIS techniques. This study has several significant points, including those that are useful for decision-makers at local administrations in Libya. The main advantage of the study is the proposal of an advanced solution for assessing the risk of sand dunes in Libya, which is a practical replacement for the previously adopted methods. The proposed solution is complete and easy to apply for generating maps that are essential for decision-making and land use planning. The information from the set of susceptibility, hazard, vulnerability, and risk maps of the sand dune in Libya can be effectively used to plan different types of land use for the reduction of deforestation. The impacts of strong wind activities on population and property are also considered.
The main results of this study are as follows. The expert opinions reveal that the land cover and soil are the most significant factors that affect the occurrence of sand dunes in the Sabha region, whereas rainfall was found to be the least critical factor. Moreover, the susceptibility map shows that about quarter of Sabha, which is mostly located in the south part of the region, is under very high proneness for sand dunes. In addition, the hazard map shows that the southern areas are most likely to experience sand dunes. The highest risks are located in the middle part, where the urban and agricultural lands are present. The sand areas are found to have very low risk of sand dunes. More than 200 km 2 of the study area are under high and very high risk zones, whereas more than half of the region is at very low risk zone. Furthermore, a substantial number of individuals are at high and very high-risk zones. The results show that more than 13,709 people are at high risk and 539 are at very high-risk zones.
Although this study has several significant points, it has the following limitations. The derivation of the controlling factors is solely based on freely available data due to the limitation and the availability of commercial data over the study area. In addition, the susceptibility analysis is based on expert-based methods, such as the AHP method; hence, other data mining and statistical techniques are not explored due to limited historical events.

Disclosure statement
No potential conflict of interest was reported by the authors.