Monitoring and forecasting of land use/land cover (LULC) in Al-Hassa Oasis, Saudi Arabia based on the integration of the Cellular Automata (CA) and the Cellular Automata-Markov Model (CA-Markov)

ABSTRACT

This study sought to integrate the Cellular Automata-Markov Model (CA-Markov) and the Cellular Automata (CA) using sensing data for land cover maps for the years: 1988, 2000, 2013 and 2020 to monitor, detect, and predict the spatial and temporal of Land Use/Land Cover (LULC) change in Al-Hassa Oasis, Saudi Arabia. The maximum likelihood classifier (MLC) method was used as a supervised classification algorithm to super control classification coefficients. Training areas were extracted to classify green and water from Landsat images using The Natural Vegetation Index (NDVI) and Normalized Difference Water Index (NDWI) indicators, respectively. The Markov form has been used Cellular Automata (CA)To create a Transition probability Matrix (TPM)), Transition Area Matrix (TAM) and Transition Suitability Maps (TSM), which is the main factor for obtaining to simulate the prediction of the future of land use/land cover in the years 2030 and 2050, while was used Cellular Automata-Markov Model (CA-Markov) To validate and Forecasting LULC. The results showed an expansion in urban areas by 63.5% equivalent to about 105.7 km2, and an increase in agricultural land by 6.6% equivalent to 10.7 km2 despite the loss of about 3.2% of the old agricultural land during the period 2000–2013 AD. Also, the area of marshes and water increased by about 41.7% equivalent to about 8.7 km, which led to the loss of 5.9% or the equivalent of about 139 km2 of arid lands during the past three decades. The projected LULC maps for 2030 and 2050 indicate that these patterns of arid lands that change into built-up areas and agricultural lands will continue over the next thirty years due to urban growth and the development of agricultural activities. The main strength of this study is the use of a relatively long period that allowed us to conduct a very good assessment of LULC in the study area, as well as to predict the different patterns for 2030 and 2050. Moreover, the results obtained in this study provide information to help decision makers to carry out effective practices for future planning and management of urban land growth and land use planning, especially in the Saudi Vision 2030.    

KEYWORDS:

• Image classification
• change prediction
• change detection
• land use/land cover-spatial suitability
• urban growth
• environmental planning

1. Introduction

Knowledge related to the monitoring and control of past and present land use/land cover plays an important role in the decision-making process (Moghadam & Helbich, 2013; Patino & Duque, 2013) and its prediction simulation can provide input data for various environmental and planning models (Dimitrios et al., 2019; Islam et al., 2018; Li, Chen, et al., 2017). Monitoring and observing land use/land cover help advance our understanding of past trends, while modeling, based on prediction simulation, provides insight into potential future developments (Appiah et al., 2017; Chim et al., 2019; Mishra et al., 2020). Both complementary approaches are considered necessary strategies to formulate better land/land cover use policies (Taubenböck et al., 2009).

A variety of land-use models have been developed to simulate and predict the future of land use/land cover change primarily by coupling remote sensing (RS) and geographic information system (GIS) data (Sohl et al., 2016). Land cover change models are tools used to analyze and provide a better understanding of the causes and consequences of land use/land cover changes, they can also be used to investigate the complex driving forces that influence spatiotemporal patterns of land use/land cover changes (Verburg et al., 2004)

Models that are characterized by simplification and abstraction play a very important role in monitoring, observing, and predicting spatial patterns of land use/land cover, and in recent years many models of monitoring and prediction simulations of land use/land cover change have been applied, the most important of which are as follows: Markov (Guan et al., 2019; Mohamed & El-Reay, 2018; Mondal et al., 2019). The SLEUTH model (a cellular automata-based computer simulation model) (Clarke & Gaydos, 1988). Multi-criteria analysis and the Analytical Hierarchy Process (AHP) (Svoray et al., 2012; Zakerinejad & Marker, 2014). Multivariable statistical methods (MSMs) such as: Linear regression (LiR) (Chaplot et al., 2005). Logistic regression (LR) (Conoscenti et al., 2014; Vanwalleghem et al., 2008). Machine Learning Methods (SVM) (Kheir et al., 2007; Märker et al., 2011; Pourghasemi et al., 2017). Artificial Neural Networks (ANN) (Märker et al., 2011).

In this study, an innovative approach has been integrated for three computer software models, two of which are within the IDRISI Selva v.17.0 program, these models are the Land Change Model (LCM) and the CA-Markov Model, whereas the third model is the multi-criteria analysis model based on geographic information systems (GIS-MCA). The Analytic Hierarchy Process (AHP) was also used to determine the weights of the criteria affecting the spatial suitability map, and the land change model (LCM) was used to monitor the features of change in the land use/land cover map during the period 1988–2020. The Land Change Model (LCM) is an integrated software environment for land use/land cover analysis and prediction. To validate the results, LCM has been included in the IDRISI program (Tewolde & Cabral, 2011), where only thematic bitmaps that have the same land cover classes are listed in the same sequential order in the land use/land cover analysis (Roy et al., 2014), then the LCM assesses land cover changes between two different times, calculates the changes, and presents the results with various graphs and maps. This model was adopted because it proved effective in predicting land use/land cover change in the study area (Megahed et al., 2015; Shooshtari & Gholamalifard, 2015).

In recent decades, the cellular automata (CA) model has been of interest to researchers because of its simplicity, flexibility, and ability to combine spatiotemporal dimensions (Chen et al., 2016; Santé et al., 2010). According to various studies, the cellular automata model is a useful tool for modeling land use/land cover changes. It is also applicable to various topics, such as urban growth prediction, traffic simulation, and land use/land cover changes (Rezazadeh & Mirahmadi, 2009; Ward et al., 2000).

The cellular automata model – Markov (CA-Markov) is based on the main hypothesis, which is the extrapolation of future growth by tracking its behavior in the past and assuming the stability of the trends of this growth and the values of change in the future. Successive satellite visualizations of the city are an important basis for this stability.

Among all the spatiotemporal simulation techniques, the CA-Markov model is the most popular and applied in most prediction and simulation studies of urban growth (Aburas et al., 2016; Li, Chen, et al., 2017; Santé et al., 2010; White & Engelen, 2000). It has been used to model land-use changes and to predict future urban growth since the 1990s (Batty & Xie, 1994; Wu, 1998, 2002), where prediction simulations are based on the assumption that “historic urban areas will influence patterns of future growth through interactions between different types of land use” (Clarke et al., 1997).

Studies related to monitoring land use changes and predicting them in the arid and semi-arid environments in general, and the Kingdom of Saudi Arabia in particular, can be divided into two main sections. The first section is concerned with studies that dealt with detecting land-use changes only without predicting. The most important thing to note in this study is that it is numerous and widespread across the Kingdom, and most of it is based on different Landsat visuals, whose archives cover different dates or visuals of the type Sentinel-2. The classification process was based on the supervised classification method or the controlled classification or the uncontrolled classification using some indicators, such as NDVI and NDWI. The maximum likelihood classifier (MLC) mostly captured the vast majority in most of the studies reviewed, and the resulting classification category ranges from four to seven categories: The most important of which are bare lands, desert lands, urban areas, agricultural lands, marshes, water, and roads. The results of the validity of the classification are verified through the kappa indicator. The most prominent studies that represent this trend are studies (Chouari 2021; Alharthi et al. 2020; Alqurashi and Kumar 2016; Abdelkarim, et al., 2020; Alqurashi & Kumar, 2014; Alsharif et al., 2022; Shrahily et al., 2022);

As for the second type of studies, they are studies that dealt with detecting changes while predicting the future of land uses. It is noted that these studies are very few in the dry environment of the Kingdom of Saudi Arabia. Moreover, these studies used the Cellular Automata-Markov Model (CA-Markov) mainly in the process of predicting land uses (Abdelkarim et al., 2022; Al-Alola et al., 2021; Alsharif et al., 2022). The research methodology, methods, and tools do not differ much in studies conducted outside the Kingdom of Saudi Arabia, with the availability of various methods for producing maps of spatial suitability (driving force), which help in the process of predicting land uses, as the methods and techniques vary in producing these maps, such as: Frequency rate, fuzzy logic, artificial neural network, MCE, AHP, decision tree, and multivariate regression lines. In this regard, many studies can be found, such as studies (Mansour et al., 2022; Pakawan & Saowanee, 2019; Solomon et al., 2019; Wang et al., 2021, 2022).

Through reviewing previous studies, it is clear that no independent and in-depth study was concentrated on detecting and predicting land uses in Al-Hassa Oasis, which clearly lacks such studies. All that was found about it is nothing more than references within a number of geographical studies that dealt with the Kingdom in general or some studies of small parts of Al-Hassa Oasis, such as the study of Chouari (Chouari, 2021) about the detection of land cover change using multi-spectral Landsat data: A case study of the wetlands in Al-Asfar Lake east of Al-Hassa city. The study focused on a very small part located east of Al-Hassa Oasis known as Al-Asfar Lake. The other study is Guma and Aljughaiman (Biro Turk & Aljughaiman, 2020) about the LULC assessment as associated with soil salinity and irrigation water above an oasis in an arid environment. The study focused on the agricultural eastern part of Al-Hassa Oasis and not the entire Oasis. In addition to producing a map of the land cover for a specific year only as in the study of Salih (Abdelatti et al., 2017; Almadini et al., 2019; Alqarni et al., 2018; Salih, 2018). It is about producing a map of the land cover in Al-Hassa Oasis according to the classification of a satellite visual of the type landsat-7 ETM for the year 2017 AD. The researcher applied the supervised classification method using the maximum likelihood classifier (MLC). Despite the increasing importance of urbanization in Al-Hassa Oasis due to oil revenues in the region, the increasing migration from rural to urban areas, recruitment, undertaking works in the oasis at the expense of the countryside, the high environmental sensitivity of the Oasis, which overcomes the sandy deserts from all sides, and the increasing pressures on the map of land uses in the large Al-Hassa Oasis during the past years, there is now a comprehensive study to detect the changes associated with the map of land uses during a clear time series to reveal the features of the effects accompanying this imposed urbanization.

The researcher has benefited greatly from previous studies in terms of enriching the scientific background, and shaping the approach followed in studies of detecting change and prediction of land uses through many available scientific references. The current study of the researcher is in line with the previous studies in terms of theoretical, scientific, methodological, and methods used in the studies of detection and prediction of changes. The most important thing that distinguishes this study from previous studies is its integration, comprehensiveness, and concentration on the land use map of a dry desert oasis that has undergone great environmental and urban changes during the past 40 years. This is important in terms of planning for political and planning decision-makers in the Kingdom of Saudi Arabia in light of the vision of the Kingdom 2030 AD in general and the increasing economic importance of the first palm oasis in the Kingdom. In addition to the future prediction of the map of land uses until the year 2050 AD through its integration of remote sensing techniques and GIS.

The problem of the study lies in the following: (1) This study area has high environmental sensitivity as a result of the characteristics of the fragile dry desert environment, and its susceptibility to rapid successive changes, as the Al-Hassa Oasis is one of the areas threatened by desertification due to the encroachment of sand, drought, soil erosion by wind, the deterioration of vegetation cover, urban growth, and shrinking agricultural area. Moreover, ignoring the environmental sensitivity of the oasis, like the removal of vegetation cover, palm trees, and change of land uses may lead to the deterioration of the natural environment of the oasis. Therefore, ecological and environmental balance considerations of the natural and biological systems represent a fundamental factor in formulating the map of future land uses in the region. (2) The continuous change in the patterns of land use/land cover and continuous urban growth represents one of the most prominent urban challenges. Besides, dealing with the challenges of land use/land cover growth in a fast-paced, large-scale, comprehensive, dimensional, and heterogenous manner, burdens any government or local administration regardless of its nobility and longevity. This planning defect posed great challenges in the urban development and distribution of services during the current and future periods.

This study aims to: (1) Monitor and observe the most important features of land use/land cover changes in Al-Hassa Oasis during the period 1988–2020. (2) Analyze the most important areas, trends, and factors of the driving forces for land use/land cover during the period 1988–2020 and their role in guiding the map of the future land use/land cover in the years 2030 and 2050. (3) Produce future land use/land cover maps for the years 2030 and 2050 for the Al-Hassa Oasis.

2. Materials and methods

2.1. Study area

Al-Hassa Oasis is located in the eastern region, about 125 km southwest of the city of Dammam, about 75 km from the Arabian Gulf, with longitude from 25°12”18.05” to 25° 38”30.01’‘North, and latitude from 49° 42’59.96” to 49° 32’32.14” East (Figure 1). Al-Hassa Oasis is bordered on the north by the city of Abqaiq and the Abu al-Hamam desert, on the west by the Al-Dahna desert and the Summan hill, on the east by the Al-Jafurah desert and Al-Asfar lake along with the Arabian Gulf, and from the south is the Al-Jafurah desert and beyond it is the Empty Quarter desert (the Rub al Khali desert).

Figure 1. The location of the study area in Saudi Arabia in 2020.

The unique geographical location of Al-Hassa, which links the Kingdom of Saudi Arabia to the countries of the Gulf Cooperation Council, on the one hand, and is between the cities of Riyadh and Dammam, on the other hand, led to the ease of commercial, industrial, cultural, and civilizational exchange, which paved the way for the emergence of ancient civilizations in the region. According to the satellite visual classification for the year 2020, Landsat 8/OLI, and its source is the United States Geological Survey (USGS), the area of Al-Hassa Oasis reached about 2729.3 km² divided as follows: The urban areas of Al-Hassa Oasis constitute about 166.62 km² with a rate of 6.1%. The agricultural lands are about 162.31 km long², with a percentage of 6%. The water areas are about 20.77 km² with a percentage of 0.8%. And the vast majority are desert lands, which are about 2379.6 km² with a percentage of 87.1%.

The population of Al-Hassa Oasis is about 1,272,633 people, representing about 95.6% of the total population of Al-Hassa Governorate, which is about 1,331,009 people. Most of the population is concentrated in the municipalities of Al Hofuf and Al Mubarraz with a percentage of 63.8%. Additionally, two of the Al-Hassa governorate’s cities, Al Hofuf and Al Mubarraz, which are among the 10 largest cities in the whole Kingdom, are located in the Al-Hassa Oasis. Al-Hassa Oasis consists of six municipalities: Al Hofuf, Al Mubarraz, Al Jafr, Al Umran, Al Uyun, and Jawatha. The largest population is in the municipality of Al Hofuf with 438,563 people, and these municipalities include about 291 neighborhoods and the population density is concentrated in the center and north of the oasis. The Oasis’s terrain is flatter and less rugged than the rest of the surrounding areas, except for some scattered hills that protrude above the general level of the oasis surface. These hills in the Al-Hassa Oasis are usually called Qarat, as is the case with Jabal Qarah in the middle of agricultural lands, and Jabal Abu Hasis. There, the levels of the land’s surface range from 152 to 156 m, while the western parts rise to reach 268 m. The slope of the surface of the Al-Hassa Oasis, in general, is directed towards the east and north-east and indicates the direction of the slope in the oasis and the adjacent areas towards the valleys in the drainage basins.

2.2. Study procedures and data processing

The study procedures and data processing, which are illustrated in Figure 2, included the main first, monitoring land use/land cover changes in Al-Hassa Oasis during the period 1988–2020, second, using the cellular automata-Markov model (CA-Markov) for predicting land use/land cover in the years 2030 and 2050.

Figure 2. Study procedures and data processing.

2.2.1. Monitoring land use/land cover changes

Planning and monitoring changes in land use/land cover (LULC) are vital for sustainable development, planning, and management (Saadat et al., 2011). The process of monitoring changes in land use/land cover for the Al-Hassa Oasis during the period 1988–2020 has gone through a set of basic stages that can be reviewed in detail as follows:

2.2.1.1. Characteristics of satellite images

Remote sensing data, such as Landsat images, is freely available on the internet and covers almost the entire Earth, which has allowed for continuous monitoring of land use/land cover change since 1988 (USGS, 2013). The first stage is to download four satellite images from the United States Geological Survey (USGS) to monitor the changes in land use/land cover that had occurred in the Al-Hassa Oasis during the period 1988–2020. The satellite images for the years 1988 and 2000 from the Thematic Mapper (TM) sensor carried on the American satellite Landsat 4–5, while the satellite images for the years 2013 and 2020 were obtained from the Operational Land Imager (OLI) sensor on the Landsat 8 satellite. The Enhanced Thematic Mapper (ETM) sensor on the Landsat 7 satellite was not used due to technical problems with the sensor, which was evident in its output images. Figures 3 and 4 and Table 1 show the characteristics of satellite images during the period 1988–2020.

Figure 3. False colour composite (FCC) of Land sat 4–5 (TM)+ bands 7(red), 4(green), 2(blue) for the years 1988, and 2000.

Figure 4. False colour composite (FCC) of Land sat 4–5 (TM)+ bands 7(red), 5(green), 3(blue) for the years 2013 and 2020.

Table 1. Characteristics of satellite images of Al-Hassa Oasis during the period 1988–2020.

Image capturing date	Satellite orbit	Standard	Satellite	Spatial resolution (m)	Sensor type	Number of bands
16/10/1988	164	42	Land sat 4-5	30	TM	7
14/10/2000	164	42	Land sat 4-5	30	TM	7
16/10/2013	164	42	Land sat 8	30	OLI	11
19/10/2020	164	42	Land sat 8	30	OLI	11

2.2.1.2. Image processing

Spectral enhancement is performed to increase spectral accuracy and ease of processing. This is done to avoid aliasing in the spectrum “spectral overlapping” in Landsat images, which might result in interference in the spectral fingerprints, so spectral enhancement of the visual is performed to increase the spectral accuracy and reduce spectral overlap through the Principal Component.

The preprocessing of the present research was divided into three main stages. (1) Radiometric correction was performed to reduce, eliminate, and correct errors in the digital numbers of the image. The radiometric correction (also known as) was performed by calculating the top-of-atmosphere (TOA) reflectance using Equations (1) and (2) (Pinto et al., 2020):

(1) ${\mathit{\rho }}^{{ }^{\mathbf{\prime }}}\mathit{\lambda }={\mathit{M}}_{\mathit{\rho }}\times {\mathrm{Q}}_{\mathit{c}\mathit{a}\mathit{l}}+{\mathit{A}}_{\mathit{\rho }}$(1)

(2) $\mathit{\rho }\mathit{\lambda }=\hspace{0.17em}\frac{{\mathit{\rho }}^{{ }^{\mathbf{\prime }}}\mathit{\lambda }}{sin\left(\mathit{\theta }\mathit{S}\mathit{E}\right)}$(2)

where $\rho { }^{\mathrm{\prime }}\lambda$ is the TOA reflectance, Up is the band multiplicative value, Ap is the band additive value, $\mathrm{Q}$^cal is the digital numbers of the satellite bands, $\rho \lambda$ is the TOA reflectance, and $\theta SE$ is the sun elevation angle. (2) The spectral bands of each Landsat datasets were combined through layer stacking to get a full color composite image comprising the spectral bands from band 1 to band 7. (3) The subsequent step was clipping the four images independently to match the boundary of the study area.

2.2.1.3. LULC classification and accuracy

In this study, it was relied upon the supervised classification method using the maximum likelihood classifier method (MLC) in the classification of Landsat satellite images for years in 1988, 2000, 2013, and 2020. Although several methods have been devised to implement supervised classification, the maximum likelihood classifier, also known as Bayes’ decision rule (Weng, 2002), is one of the most widely used supervised classification algorithms due to its availability and simple training process (Jensen & Cowen, 1999). In this method, satellite image pixels are classified based on their probability of belonging to a particular LULC category. The concept of the method assumes equal probabilities for all classes and a normal distribution for all input bands, and the input bands used in the study to produce composite maps with false colors consist of bands 7 ، 4 ، 2 ، for Land sat 4–5 (TM) and 7 ، 5 ، 3 ، for Landsat 8 (OLI), the spectral signature of each image pixel was matched to the training samples of the study area and satellite images classified into 4 broad land use/land cover types (LULC). They are as follows (water, agriculture, desert land, and desert) (Table 2).

Table 2. LULC classes in the study area and the corresponding description.

LULC Classes	Description
Lakes	All forms of water surfaces (clean-turbid), including rain water swamps, sewage estuaries, open surface water tanks, rain marshes and salt marshes.
Agricultural land	Palm trees, cultivated lands, grass, agriculture, and croplands
Desert lands	Desert, sand, mountains, vacant land, bare exposed rocks, and dry salt areas
Urban and built-up	Residential and commercial settlements, industrial areas, urban agglomerations, roads, bridges, transportation networks, infrastructure, and urban amenities, educational and health services

Satellite classification and visualization for each year were added together in two Viewer spaces. The classification accuracy process was performed through Accuracy Assessment, the classification path to be tested for its accuracy was determined, and random points were added that are used to test the classification accuracy (number of random sampling points = 1000 points). Accordingly, fieldwork data was collected in January 2020 to identify land categories based on prior knowledge of the study area during field visits. A GPS device was used to collect and locate ground facts and ground control points (GCPs). The selected reference pixels were representatives of the four land cover categories: Water, agriculture, arid land, urban land, and buildings. Overall, 75% of the digital reference pixels were used for training and 25% were used for validation (Table 3).

Table 3. Description of digitized reference pixels.

Classes				Sampling
Lakes	Agricultural land	Desert lands	Urban and built-up	Training	Validation	Total
60	240	300	400	750	250	1000

The analysis of the Landsat images showed that each class has a distinctive spectral signature. However, the Landsat medium-resolution images have lots of mixed pixels. Consequently, the differentiation between them, especially in urban areas, was sometimes not accurate and induced confusion. Two critical processing techniques were performed on the pre-processed images, the Natural Vegetation Index (NDVI), which is an indicator that measures the difference between cultivated and uncultivated land and provides a measure of plant density and condition. Land cover is affected by vegetation, vegetation density, and vegetation greening. This indicates that the photosynthetic power of the Earth’s surface cover, and Normalized Difference Water Index (NDWI) was used to extract the vegetated areas and waterbodies, respectively, from the successive Landsat images. The NDVI and NDWI were extracted and analyzed individually for all Landsat data (using thresholds specific to each image) As the Landsat-8 data have narrower spectral bands. (Figures 5 and 6) (Taufik et al., 2016), according to the following Equations (3)(3) $\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}=\hspace{0.17em}\frac{\mathit{N}\mathit{I}\mathit{R}-\mathit{R}\mathit{e}\mathit{d}}{\mathit{N}\mathit{I}\mathit{R}+\mathit{R}\mathit{e}\mathit{d}}$(3) and (4(4) $\mathrm{N}\mathrm{D}\mathrm{W}\mathrm{I}=\hspace{0.17em}\frac{\mathit{N}\mathit{I}\mathit{R}-\mathit{S}\mathit{W}\mathit{I}\mathit{R}}{\mathit{N}\mathit{I}\mathit{R}+\mathit{S}\mathit{W}\mathit{I}\mathit{R}}$(4) ) (Holtgrave et al., 2020).

(3) $\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}=\hspace{0.17em}\frac{\mathit{N}\mathit{I}\mathit{R}-\mathit{R}\mathit{e}\mathit{d}}{\mathit{N}\mathit{I}\mathit{R}+\mathit{R}\mathit{e}\mathit{d}}$(3)

(4) $\mathrm{N}\mathrm{D}\mathrm{W}\mathrm{I}=\hspace{0.17em}\frac{\mathit{N}\mathit{I}\mathit{R}-\mathit{S}\mathit{W}\mathit{I}\mathit{R}}{\mathit{N}\mathit{I}\mathit{R}+\mathit{S}\mathit{W}\mathit{I}\mathit{R}}$(4)

Figure 5. NDVI of Landsat imagery from 1988 to 2020.

Figure 6. NDWI of Landsat imagery from 1989 to 2020.

where NIR is the near-infrared band (band 4 in Landsat-5, Landsat-7 and band 5 in Landsat-8), Red is the red band (band 3 in Landsat-5, Landsat-7 and band 4 in Landsat-8), and SWIR (band 5 in Landsat-5, Landsat-7, and band 6 in Landsat-8) is the short-wave infrared band.

The last step in the classification process was to verify the accuracy of the image classification by comparing the results of the supervised classification with the basic fact information from producers and users, and overall accuracy, as well as the Kappa accuracy statistic. The accuracy of the classification was estimated using the stratified random sampling method via the accuracy assessment tool in ERDAS Imagine software, through stratified randomly sampled GTD (Wang et al., 2021). Stratified random sampling method is a common method of assessing accuracy. This sampling method allows for the possibility of equivalent sampling for all LULC categories in the entire study area without bias (Dey et al., 2021; Kafy et al., 2021; Wang et al., 2021). The accuracy of the product and user was calculated as suggested by (Fung & Ledrew, 1988) to obtain the average accuracy of the visualizations used in the study from 1988 to 2020, the classified image was compared with the ground fact data to verify and calculate the accuracy of the classification. The accuracy confusion matrices showed that the overall accuracy was 91.23%, 93.65%, 95.32%, and 98.66% and the Kappa coefficients were 89.64, 92.41, 93.66, and 96.98 for the four Landsat images in 1988, 2000, 2013, and 2020, respectively (Table 4). The kappa coefficient was calculated as follows: Equations (5(5) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}=\frac{\mathit{N}{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}\mathit{i}-{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}+\left(\mathit{X}+\mathit{i}\right)}{{\mathit{N}}^2-{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}+\left(\mathit{X}+\mathit{i}\right)}\hspace{0.17em}\times 100\mathbf{\%}$(5) ,6(6) $\mathbf{O}\mathbf{v}\mathbf{e}\mathbf{r}\mathbf{a}\mathbf{l}\mathbf{l}\mathbf{A}\mathbf{c}\mathbf{c}\mathbf{u}\mathbf{r}\mathbf{a}\mathbf{c}\mathbf{y}=\frac{{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}\mathit{i}}{\mathit{N}}\hspace{0.17em}\times 100\mathbf{\%}$(6) ,7(7) $\mathbf{U}\mathbf{s}\mathbf{e}{\mathbf{r}}^{\prime }\mathbf{s}\mathbf{A}\mathbf{c}\mathbf{c}\mathbf{u}\mathbf{r}\mathbf{a}\mathbf{c}\mathbf{y}=\frac{\mathit{X}\mathit{i}\mathit{i}}{\mathit{X}+\mathit{i}}\hspace{0.17em}\times 100\mathbf{\%}$(7) and 8(8) ${\mathbf{P}\mathbf{r}\mathbf{o}\mathbf{d}\mathbf{u}\mathbf{c}\mathbf{e}\mathbf{r}}^{\prime }\mathbf{s}\mathbf{A}\mathbf{c}\mathbf{c}\mathbf{u}\mathbf{r}\mathbf{a}\mathbf{c}\mathbf{y}=\frac{\mathit{X}\mathit{i}\mathit{i}}{\mathit{X}+\mathit{i}}\hspace{0.17em}\times 100\mathbf{\%}$(8) ) (Faichia et al., 2020; Khawaldah, 2016).

Table 4. Overall classification accuracy and Kappa statistics.

Year	LULC category	Producer accuracy (%)	User accuracy (%)	Overall accuracy (%)	Khat (%)
1988	Lakes	92.85	93	91.23	89.64
	Agricultural land	84.65	86
	Desert lands	98.88	87
	Urban and built-up	98.12	95
2000	Lakes	90.70	92	93.65	92.41
	Agricultural land	86.77	85
	Desert lands	99.18	88
	Urban and built-up	99.45	94
2013	Lakes	91.66	95	95.32	93.66
	Agricultural land	83.23	88
	Desert lands	100.00	86
	Urban and built-up	99.70	97
2020	Lakes	94.66	96	98.66	96.98
	Agricultural land	87.45	90
	Desert lands	100.00	88
	Urban and built-up	100.00	98

(5) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}=\frac{\mathit{N}{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}\mathit{i}-{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}+\left(\mathit{X}+\mathit{i}\right)}{{\mathit{N}}^2-{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}+\left(\mathit{X}+\mathit{i}\right)}\hspace{0.17em}\times 100\mathbf{\%}$(5)

(6) $\mathbf{O}\mathbf{v}\mathbf{e}\mathbf{r}\mathbf{a}\mathbf{l}\mathbf{l}\mathbf{A}\mathbf{c}\mathbf{c}\mathbf{u}\mathbf{r}\mathbf{a}\mathbf{c}\mathbf{y}=\frac{{\sum }_{\mathit{i}=1}^{\mathit{r}}\mathit{X}\mathit{i}\mathit{i}}{\mathit{N}}\hspace{0.17em}\times 100\mathbf{\%}$(6)

(7) $\mathbf{U}\mathbf{s}\mathbf{e}{\mathbf{r}}^{\prime }\mathbf{s}\mathbf{A}\mathbf{c}\mathbf{c}\mathbf{u}\mathbf{r}\mathbf{a}\mathbf{c}\mathbf{y}=\frac{\mathit{X}\mathit{i}\mathit{i}}{\mathit{X}+\mathit{i}}\hspace{0.17em}\times 100\mathbf{\%}$(7)

(8) ${\mathbf{P}\mathbf{r}\mathbf{o}\mathbf{d}\mathbf{u}\mathbf{c}\mathbf{e}\mathbf{r}}^{\prime }\mathbf{s}\mathbf{A}\mathbf{c}\mathbf{c}\mathbf{u}\mathbf{r}\mathbf{a}\mathbf{c}\mathbf{y}=\frac{\mathit{X}\mathit{i}\mathit{i}}{\mathit{X}+\mathit{i}}\hspace{0.17em}\times 100\mathbf{\%}$(8)

where N is the number of all pixels used for observation, $\mathit{r}$ is the number of rows in the error matrix (number of classes), $Xii$ is the diagonal value of the contingency matrix of row $i$ and column $i$, $X+i$ is column pixel number $i$, and $Xi$ + row pixel number $i$.

2.2.2. Monitor changes using the Land Change Model (LCM)

Geographical information systems and satellite-based image analysis can be used as effective tools for urban planning. The satellite visuals provide updated information on the land cover (Piekarczyk, 2014). At this stage, the results of classifying the visuals in 1988, 2000, 2013, and 2020, which had been processed and produced in the previous stages, were adopted and entered directly into the Land Change Model (LCM). Land Change Modeler (LCM) built in TerrSet software version 18.31 was used to simulate future LULC changes. It is an empirically driven stepwise process performing change analysis, transition potential modeling, and change prediction (Eastman, 2015) based on the historical.

The land change model (LCM) is an innovative land planning and decision support system integrated into the IDRISI Selva v. 17.0 program, which was developed by Clark Labs to analyze land cover change, and experimentally model the relationship between land change and explanatory variables. The LCM tool was applied to analyze land use/land cover change for the Al-Hassa Oasis during the period 1988–2020. The “Validate” tool was used in the IDRISI program to verify the accuracy of the model, which eventually showed very high accuracy for different Kappa coefficients.

2.2.3. Markov Model Cellular Automata (CA)

The Markov form is used Cellular Automata (CA) to create a Transition Probability Matrix (TPM), Transition Area Matrix (TAM), and Transition Suitability Maps (TSM) (Ozturk, 2015) Which is the main factor and is the key to obtaining simulation of the land use in Cellular Automata-Markov Model (CA-Markov) (Eastman, 2015; Mishra et al., 2014), and the TPM is obtained by scheduling two multi-time images, each category has the potential to change ground cover to another within a predetermined period of time, and the TAM. On the estimated number of pixels that may change during a specific number of time units from one specific ground cover category to another, while reflecting the map of the convenience of the transition (TSM) about the possibility of each land cover class pixel in a place with the designated category after a specified time (Nadoushan et al., 2015). In this study, land use data for the four time periods used in the study were entered in ASCII format within the environment of the IDRISI selva v.17.0 program. In addition, the Cellular Automata (CA) model was used to obtain (TSM), (TPM), and (TAM) for the time period 2000–2010 for the purpose of obtaining a simulation of the spatial distribution of land uses in 2020. It was also used to obtain (TSM), (TPM), and (TAM) for the time period 2010–2019 for the purpose of obtaining a simulation of the spatial distribution of land uses in the year 2030, 2050. The Cellular Automata (CA) Model can be explained as follows in the Equation (9)(9) $\mathit{S}\left(\mathit{t},\mathit{t}+1\right)=\mathit{f}\left(\mathit{s}\left(\mathit{t}\right),\mathit{N}\right)$(9) (Sang et al., 2011).

(9) $\mathit{S}\left(\mathit{t},\mathit{t}+1\right)=\mathit{f}\left(\mathit{s}\left(\mathit{t}\right),\mathit{N}\right)$(9)

Where S is the set of limited and discrete cellular states, N is the Cellular field, t and t + 1 indicate the different times, and f is the transformation rule of cellular states in local space.

2.2.4. Markov Cellular Automation Model

Accurate simulations and land use/land cover predictions are necessary to manage urbanization and explicitly address trends and spatiotemporal distributions, and the CA-Markov model is one of the most widely used models for this purpose (Tingting et al., 2019). The process of future prediction of land use/land cover in Al-Hassa Oasis was carried out using the IDRISI selva v. 17.0 program by applying the famous CA-Markov prediction model, which was developed at Clark University laboratories to predict a certain phenomenon during a particular period. Because the prediction process requires certain information in different periods for the same area under study, multi-temporal data are very effective and successful in achieving the prediction process by tracking the stages of development of phenomena at different times. The result of the Markov model depends on the transition probability (Coppedge et al., 2007). And the transition probability (Pij) between state (j) and (i) is the probability that the land cover (i in pixels) in time (x) changes to the land cover class (j) in time (x + 1), the change dynamics for any study area depends on the past or present land cover state (Coppedge et al., 2007) and is mathematically calculated as approved by (Subedi et al., 2013) using the two Equations (10)(10) ${\mathit{L}}_{\left(\mathit{x}+1\right)}={\mathit{P}}_{\mathit{i}\mathit{j}}\times {\mathit{L}}_{\left(\mathit{x}\right)}$(10) and (11(11) $P_{ij}=\left[\begin{array}{c}{P}_{11}\\ P_{21}\\ P_{31}\\ P_{m1}\end{array}\begin{array}{c}{P}_{12}\\ P_{22}\\ P_{32}\\ P_{m2}\end{array}\begin{array}{c}{P}_{13}\\ P_{23}\\ P_{33}\\ \cdots \end{array}\begin{array}{c}\cdots \\ \cdots \\ \cdots \\ \cdots \end{array}\begin{array}{c}{P}_{1m}\\ P_{2m}\\ P_{3m}\\ P_{mm}\end{array}\right]$(11) ).

(10) ${\mathit{L}}_{\left(\mathit{x}+1\right)}={\mathit{P}}_{\mathit{i}\mathit{j}}\times {\mathit{L}}_{\left(\mathit{x}\right)}$(10)

Transition probabilities are obtained from the transition samples that occur during a specified period ((Adhikari & Southworth, 2012)^{) and are presented in the transition matrix (P).}

(11) $P_{ij}=\left[\begin{array}{c}{P}_{11}\\ P_{21}\\ P_{31}\\ P_{m1}\end{array}\begin{array}{c}{P}_{12}\\ P_{22}\\ P_{32}\\ P_{m2}\end{array}\begin{array}{c}{P}_{13}\\ P_{23}\\ P_{33}\\ \cdots \end{array}\begin{array}{c}\cdots \\ \cdots \\ \cdots \\ \cdots \end{array}\begin{array}{c}{P}_{1m}\\ P_{2m}\\ P_{3m}\\ P_{mm}\end{array}\right]$(11)

Where:

L (x) and L(x+1) are the land use/land cover conditions at (x + 1) and (x), respectively. $0<P_{ij}<1and\sum _{j=1}^m{P}_{ij}=1,(i,j=1,2,3,\dots ,m\mathrm{\_}$ is the transition probability matrix.

The kappa values of the CA – Markov model for LUCC simulation ranged from −1 to 1, where positive values are a sign of agreement, and negative values illustrate a lack of agreement Kappa _ 0.5 shows high agreement, 0.5 _ Kappa _ 0.75 marks a moderate level of agreement, and 0.75 _ Kappa <1 (Faichia et al., 2020; Omar et al., 2013; Ruben et al., 2020) indicates a high level of agreement. Equations (12(12) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{n}\mathbf{o}\mathbf{i}\mathbf{n}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{m}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}{\mathbf{K}}_{\mathbf{n}\mathbf{o}}=\frac{\left({\mathbf{M}}_{\left(\mathbf{m}\right)}{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}{\left({\mathbf{P}}_{\left(\mathbf{p}\right)}-{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}$(12) ,13(13) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{n}\mathbf{o}\mathbf{i}\mathbf{n}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{m}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}{\mathbf{K}}_{\mathbf{l}\mathbf{o}\mathbf{c}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}}=\frac{\left({\mathbf{M}}_{\left(\mathbf{m}\right)}{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}{\left({\mathbf{P}}_{\left(\mathbf{m}\right)}-{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}$(13) and 14(14) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{n}\mathbf{o}\mathbf{i}\mathbf{n}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{m}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}{\mathbf{K}}_{\mathbf{s}\mathbf{t}\mathbf{a}\mathbf{n}\mathbf{d}\mathbf{a}\mathbf{r}\mathbf{d}}=\frac{\left({\mathbf{M}}_{\left(\mathbf{m}\right)}{\mathbf{M}}_{\left(\mathbf{n}\right)}\right)}{\left({\mathbf{P}}_{\left(\mathbf{p}\right)}-{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}$(14) ) for the summary statistics are:

(12) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{n}\mathbf{o}\mathbf{i}\mathbf{n}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{m}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}{\mathbf{K}}_{\mathbf{n}\mathbf{o}}=\frac{\left({\mathbf{M}}_{\left(\mathbf{m}\right)}{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}{\left({\mathbf{P}}_{\left(\mathbf{p}\right)}-{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}$(12)

(13) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{n}\mathbf{o}\mathbf{i}\mathbf{n}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{m}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}{\mathbf{K}}_{\mathbf{l}\mathbf{o}\mathbf{c}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}}=\frac{\left({\mathbf{M}}_{\left(\mathbf{m}\right)}{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}{\left({\mathbf{P}}_{\left(\mathbf{m}\right)}-{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}$(13)

(14) $\mathbf{K}\mathbf{a}\mathbf{p}\mathbf{p}\mathbf{a}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{n}\mathbf{o}\mathbf{i}\mathbf{n}\mathbf{f}\mathbf{o}\mathbf{r}\mathbf{m}\mathbf{a}\mathbf{t}\mathbf{i}\mathbf{o}\mathbf{n}{\mathbf{K}}_{\mathbf{s}\mathbf{t}\mathbf{a}\mathbf{n}\mathbf{d}\mathbf{a}\mathbf{r}\mathbf{d}}=\frac{\left({\mathbf{M}}_{\left(\mathbf{m}\right)}{\mathbf{M}}_{\left(\mathbf{n}\right)}\right)}{\left({\mathbf{P}}_{\left(\mathbf{p}\right)}-{\mathbf{N}}_{\left(\mathbf{n}\right)}\right)}$(14)

where no information is N(n), medium grid cell-level information is M(m), and perfect grid cell-level information across the landscape is P(p).

2.2.4.1. Through the following steps, it is possible to predict the future land use employing the CA-Markov model

1. Entering the land use/land cover data in the six time periods shown in ASCII format within the environment of IDRISI selva v.17.0 software and using the CA-Markov model to obtain maps of appropriate transition, transition probability matrix, and transition area matrix for the period 2000–2013 to obtain a simulation spatial distribution of the land use map in 2020.

2. Predicting the map of land use/land cover for the year 2020 using the cellular automata model – Markov (CA-Markov) on the map of land use/land cover for the year 2013 and the matrix of transition areas resulting from the analysis of the maps of the years 2000 and 2013 in the previous step. The number of cycles was determined by 7 according to the number of years to be predicted from 2013 to 2020. The default filter 5 × 5 was used and a default error rate of 0.15 was determined.

3. Comparing the map of land use/land cover for the year 2020 resulting from the prediction with the existing map, which was obtained by the oriented classification of the satellite visual for the year 2020, using the “Validate” simulation accuracy tool in IDRISI program, and this step was to verify the accuracy of the simulation so that we can be assured of the results that will be obtained for the simulations in 2030 and 2050.

4. Using CA-Markov within the environment of the IDRISI selva v.17.0 program to obtain appropriate transition maps, the transition probability matrix, and the transition area matrix for the period 2013–2020 to obtain a simulation of the spatial distribution of land uses in the years 2030 and 2050.

5. Predicting of land use/land cover maps for the years 2030 and 2050 using the cellular automata model – Markov (CA-Markov). The number of cycles was determined by 10 and 30 cycles according to the number of years to be predicted, starting from 2020 until 2030 and 2050. The default filter 5 × 5 was used, and a default error rate of zero was determined.

3. Results

3.1. Monitoring land use/land cover changes

3.1.1. Monitoring the total changes of the land use/land cover map

The map of Transition of land use/land cover changes in Al-Hassa Oasis during the period 1988–2020 shown in Figure 7 revealed that the lands that were transformed into urban areas reached an area of 105.7 km², and the largest percentage of these changes included the desert areas with a percentage of 93.33% of its area, which included the desert lands adjacent to the center of the oasis with a north–south extension with the presence of urban scattered to the west such as the National Guard and the industrial zone, in addition to the expansions at the expense of raising agricultural lands in the four main oases around the oasis, and the lands that became agricultural reached an area of 10.7 km².

Figure 7. Total losses and gains in land use/land cover from 1988 to 2020.

The largest percentage of these changes was at the expense of the desert areas around the four oases by 98.62%, in addition to the agricultural expansion towards the south of the oasis, and the lands that turned into salt lakes reached an area of 8.7 km², and the largest percentage of these changes led to increasing the Al-Asfar salt lake, northeast of the oasis, by almost 100%, and the area of lands that turned into roads reached 15.8 km², which included the largest percentage of these changes at the expense of desert areas by 89.36%. Table 5 shows the changes in land use/land cover during the period 1988–2020, and Table 6 shows Statistics of LUCC classification for 1988, 2000, 2013, and 2020, and Figure 8 shows Gains and losses in each type of land-use in 1988–2000, 2000–2013, 2013–2020, and 1988–2020, Figure 9 shows Spatial trend surface (fourth-order trend) of LULC dynamics in 2020.

Figure 8. Gains and losses in each type of land-use in 1988–2000, 2000–2013, 2013–2020, and 1988–2020.

Figure 9. Spatial trend surface (4^th-order trend) of LULC dynamics in 2020: (a) Urban and built-up to Agriculture, (b) Desert lands to Urban and built-up, (c) Agriculture to Urban and built-up, (d) All to Urban and built-up.

Table 5. Changes in land use/land cover during the period 1988–2020.

Classification number	Land Use Categories	The area of change during the period 1988–2020 (km^2)
1	Desert	2476.15
2	Agriculture into urban	2.6
3	Desert into urban	103.5
4	Roads into urban	0.49
5	Urban into agriculture	0.012
6	Desert into agriculture	10.7
7	Agriculture into water	0.002
8	Desert into water	8.7
9	Agriculture into desert	2.2
10	Water into desert	0.29
11	Roads into desert	1.08
12	Urban into roads	0.02
13	Agriculture into roads	1.4
14	Desert into Roads	15.8

Table 6. Statistics of LUCC classification for 1988, 2000, 2013, and 2020.

LU Types	1988		2000		2013		2020
	Area (km^2)	Area (%)	Area (km^2)	Area (%)	Area (km^2)	Area (%)	Area (km^2)	Area (%)
Lakes	12.11	0.44	21.18	0.77	18.72	0.68	20.77	0.76
Agricultural	151.63	5.53	162.36	5.92	157.19	5.73	162.31	5.92
Urban and built-up	60.89	2.22	69.1	2.52	128.79	4.71	166.62	6.1
Desert lands	2519.41	91.81	2491.4	90.79	2439.34	88.9	2394.34	87.26
Total	2744.04	100	2744.04	100	2744.04	100	2744.04	100

3.1.2. Monitoring changes in urban areas during the period 1988–2020

Al-Hassa Oasis witnessed a great boom and an unprecedented development of the urban space, as urban construction increased during that period approximately 2.7 times its area, from 60.89 km² in 1988 to about 166.6 km² in 2020, and the urban increase amounted to about 105.71 km², and the area of increase was distributed over 32 years, with a growth rate of 3.2%, and an average annual increase of 3.3 km² over that period. The urban area in 1988 represented about 60.89 km², due to the urban expansions and extensions of Al Hofuf, especially to the south and west, such as Khalidiya in the south, the Industrial City, and National Guard in the west at the expense of desert areas. This increase continued during the period between 1988 and 2000, until it reached 69.1 km² in 2000, with an increase of 13.5% compared to what it was in 1988, due to the continued increase in urban growth to the south and east at the expense of desert areas and the noticeable increase in rural urbanization villages, especially Al-Muazil to the south and At Taraf to the east at the expense of agricultural lands (Figures 10, 11, and 12)

Figure 10. Change in land use/land cover for Al-Hassa Oasis between 1988 and 2000.

Figure 11. Change in land use/land cover for Al-Hassa Oasis between 2013 and 2020.

Figure 12. A comparison of urban growth and Agricultural land between 1988 and 2000.

This increase also continued during the period 2000–2013 until it reached 128.79 km2 in 2013, with an increase of 86.4% compared to 2000, due to the steady urban growth that occurred during this period and the connection of some new urban communities to Al Hofuf, such as Khalidiya in the south and Mubarraz in the west at the expense of desert areas and the significant increase in the area of rural-urban villages such as At-Taraf and Al Umran at the expense of agricultural lands. Likewise, the increase continued during the period between 2013 and 2020, until it reached 166.62 km² in 2020, with an increase of 29.4% compared to what it was in 2013, because of urban growth and extensions, especially in the southern extension and in the area between Al Hofuf and the National Guard in the west and Al-Muaizil in the north. Figure 13 clarifies the relationship between urban growth and population growth at the level of the municipalities of Al-Hassa Oasis in 2020, whereas the municipality of Al Hofuf ranked first in terms of urban growth by 32.9%, and the population thereof reached 438,563 people, and the municipality of Al Mubarraz ranked second in terms of urban growth, which reached 30.9% with population reached 411,220 people, whereas the municipality of Juatha ranked last with 4.9% of urban growth and its population reached 65,318 people.

Figure 13. Relationship between urban growth and population growth in 2020.

3.1.3. Monitoring changes in agricultural lands during the period 1988–2020

Al-Hassa Oasis witnessed many changes in the area of agricultural lands, as agricultural land was added during the period 1988–2020, approximately 1.07 times its area, from 151.63 km² in 1988 to about 162.31 km² in 2020, and the increase in agricultural lands amounted to about 10.67 km². The area of increase was distributed over 32 years, with a growth rate of 0.2%, and an average annual increase of 0.33 km² over that period.

The area of agricultural lands in 1988 was about 151.63 km², the increase continued during the period between 1988 and 2000 until it reached 162.36 km² in 2000, with an increase of 7.1% compared to what it was in 1988, to increase the agricultural area, especially in the northern parts at the expense of the desert areas. However, this area declined during the period between 2000 and 2013, reaching 157.19 km² in 2013, with a decrease of 3.2% compared to what it was in 2000. The reasons for this are due to the lack of rain and high salinity, especially in the eastern oasis in the area east of Jabal Qarah, Al Umran, and Al-Asfar regions. And in some places of the northern oasis that surround the villages of Al Shuqaiq and Al Jalila, and also because those agricultural lands were affected by the phenomenon of sand encroachment, as well as the urban encroachment on agricultural lands. The area of agricultural lands increased again during the period between 2013 and 2020, until it reached 162.31 km² in 2020, with an increase of 3.3% compared to what it was in 2013, due to the increase in the agricultural scattered spread in the northern and southern parts at the expense of the desert areas (Figures 10, 11 and 12).

3.1.4. Monitoring changes in water (marshes)

Marshes surround the Al-Hassa Oasis from all sides and spread across various parts of the Oasis land. The reason for the formation of marshes is the flat surface of the oasis and the surrounding areas from the north, east, and south, and the slight slope to the east and northeast. The existence of depressions inside the Oasis or the adjacent areas in which rainwater is collected led to the formation of some marshes. The area of water in 1972 represents about 4.11 square kilometers (topographic maps scale 5000:1) concentrated in marshes Al-Asfar northeast of Al-Hassa, which is a drainage area for adjacent agricultural lands. The area of water began to increase between 1988 and 2000, reaching 21.14 square kilometers in 2000, an increase of 57.3% over 1988, owing to an increase in the expansion of water in marsh. The area of water again decreased between 2000 and 2013, reaching 18.72 km² in 2013, a decrease of 11.4% from 2000, owing to the decline in water again in marsh. Furthermore, the water area increased again during the period (2013–2020) to reach 20.77 square kilometers in 2020, with an increase of 11% from what it was in 2013, due to the increase in the water area in the marsh (Figures 14 and 15). The reasons for water decline and expansion in marsh during time are attributed to the amount of rainfall and water discharge, as they generally expand in the winter and recede in the summer. It also expands in rainy years and recedes in dry years. Yellow Lake or Yellow marsh is exposed to a range of geomorphological hazards, the most important of which is the change of the lake’s shore line from year to year, the creep of sand from year to year, the increase in salinity, the change of chemical properties, and the exposure of the lake to salt weathering, which exposes the beer to significant environmental changes.

Figure 14. Detection of changes in land use/land cover for Al-Hassa Oasis between 1988 and 2000.

Figure 15. Change in land use/land cover for Al Asfar lake between 2013 and 2020.

3.2. Predicting the future of land use/land cover for Al-Hassa Oasis until the year 2050

3.2.1. The actual and predicted land use/land cover map for the year 2020

Land use modeling and simulation are among the main methods of land use/land cover study, and they are important tools for exploring the drivers of land use, supporting urban planning, policy formulation, and evaluating land use impacts on the ecological environment (He et al., 2017; Kuang, 2011). The cellular automata model – Markov (CA-Markov) was used on the layers of land use/land cover in the years 2000 and 2013, which are agricultural lands, urban areas, desert areas, water, and roads. And the two layers of land use/land cover in the years 2000 and 2013 were entered and the value of the interval between the two images was set to 13 indicating the difference of years between the land use/land cover in both images, and the number of years was set to 7, referring to the difference of years between 2013 and 2020, so that we can predict the map of land uses for the year 2020 and then compare it with the map of the actual land use/land cover for the year 2020, which was obtained from the oriented classification of the satellite visual of the captured area in 2020, to verify the accuracy of the model so that we can be assured of the prediction results for the year 2030 and 2050, in this context the transition probability matrix was obtained (Table 7) and the matrix of transition areas, Table 8 of land use for the years 2000 and 2013.

Table 7. Matrix of the transition probability between different land uses/land cover during the period 2000–2013.

	2013
Year	Criterion	Agricultural land	Urban and built-up	Lakes	Desert lands	Roads Networks
2000	Agricultural land	0.9893	0.0025	0.0000	0.0082	0.0000
	Urban and built-up	0.0231	0.9152	0.0000	0.0617	0.0000
	Lakes	0.0000	0.0000	0.9791	0.0209	0.0000
	Desert lands	0.0144	0.0046	0.0008	0.9791	0.0010
	Roads Networks	0.0003	0.0000	0.0000	0.0012	0.9985

Table 8. Matrix of transition values (number of Pixels) for different land use/land cover during the period 2000–2013.

	2013
Year	Criterion	Agricultural land	Urban and built-up	Lakes	Desert lands	Roads Networks
2000	Agricultural land	0.0000	0.0000	0.0000	326158	1521
	Urban and built-up	0.0000	0.0000	437	385171	3711
	Lakes	0.0000	0.0000	4349	91335	7
	Desert lands	344674	411182	91251	10954860	105181
	Roads Networks	0.0000	0.0000	51	110706	677

The CA-Markov model was implemented, and the number of cycles was set to be 7 years equal to the number of years to be predicted. The simulation map for the year 2020 was obtained (Figure 17) and compared with the actual land use/land cover map for the year 2020.

Table 9 and Figure 16 show the accuracy of the prediction model, which ranged between 77.05% and 100%, with an average accuracy of about 91.96% for different land use/land cover, as the area of agricultural lands resulting from the forecast for 2020 was about 151.77 km², while the actual area resulting from the classification was 162.31 km² with an accuracy of 93.51%, as for the urban mass, the area resulting from the prediction was 128.38 km², while the actual area was 166.62 km² with an accuracy of 77.05%, and the desert areas reached an area resulting from the prediction of 2425.98 km², while the actual area was 2379.6 km² with an accuracy of 98.09%, and the area of the water resulting from the prediction was 23.17 km², while the actual area was 20.77 km² with an accuracy of 89.64%, and finally, the area of the road network in the prediction map was 38.36 km², and the actual area was 38.36 km² with an accuracy of 100%.

Figure 16. The comparison land use/land cover between real map 2020 and simulated map 2020.

Figure 17. (a) the actual map of the land use/land cover in 2020. (b) the map of the land use/land cover in 2020 resulted from the prediction.

Table 9. A comparison between the statistics of the prediction map and the actual map for the year 2020 to evaluate the accuracy of the model.

Criterion	Areas in square kilometers
	Actual (2020)	Predicted (2020)	Prediction accuracy %
Agricultural land	162.31	151.77	93.51
Urban and built-up	166.62	128.38	77.05
Lakes	20.77	23.17	89.64
Desert lands	2379.6	2425.98	98.09
Roads Networks	38.36	38.36	100
Sum	2767.66	2767.66	2767.66

3.3. Model validation

We did not just verify the accuracy of the model by comparing the areas of the different types of land uses in the prediction map for the year 2020 and their analogies to the actual land use map for the same year as shown. The Kappa value shows the agreement between the actual and reference map (equal to an agreement of chance), and the upper and lower bounds for kappa are 1.00 (occurs when there is total agreement) and −1.00 (represents agreement less than chance) (Congalton, 1991).

The “Validate” tool was used in the IDRISI program to verify the accuracy of the model, which showed very high accuracy for the various Kappa coefficients (Table 10) and (Figure 18) where the value of the coefficient Kno equals (0.9421), the Klocation coefficient equals (0.8734), the KlocationStrata coefficient equals (0.8734), and the coefficient of Kstandard equals (0.8091), after checking the accuracy of the model and using kappa coefficients, we were confident to use the model to predict land uses for the years 2030 and 2050.

Figure 18. Kappa parameters from the comparison of the simulated LULC map and the reference map.

Table 10. Kappa indicators for evaluating the accuracy of a prediction model.

Kappa index	Kappa value for the CA-Markov model
Kno	0.9421
Klocation	0.8734
KlocationStrata	0.8734
Kstandard	0.8091

According to (Zadbagher et al., 2018) a model is valid if the overall Kappa score (K_standard) exceeds 70% or (0.7), a K_standard score near 90% is a very strong indicator of overall accuracy and model performance, and the remaining K scores, which all exceed 85%, that there are almost no errors or quantification and a very small location between the predicted and actual land cover map for the year 2020, and therefore the simulation has a strong ability to predict both the amount and locations of change.

3.4. Predicting the future of land uses/land cover for the years 2030 and 2050

Creating an accurate simulation of future land use/land cover prediction is one of the most important challenges in urban studies involving spatial modeling (Aburas et al., 2017), and future land use/land cover prediction is simulated using the Markov model (Balzter, 2000). The CA-Markov model has been used on the layers of land use/land cover for the years 2013 and 2020, represented by agricultural lands, urban mass, desert areas, water, and roads. The value of the time interval between the two images was adjusted to 7 years, indicating the difference of years between the land cover in both images, and the number of years after the second image was set to 10, indicating the difference of years between 2020 and 2030. This is to be able to predict the land use/cover map for the years 2030 and 2050, and in this context, the transition probability matrix, (Table 11) and the transition areas matrix, (Table 12) and appropriate maps for the transition from land use/land cover for the years 2013–2020 were obtained. The cellular automata model – Markov (CA-Markov) was implemented and the number of cycles was set to be 10, 30 equals to the number of years to be predicted, and a simulation map was obtained for the years 2030 and 2050, assuming that the rates of development of all different types of land use/land cover are maintained.

Table 11. Matrix of the transition probability between different land use/land cover during the period 2020–2030.

	2030
Year	Criterion	Agricultural land	Urban and built-up	Lakes	Desert lands	Roads Networks
2020	Agricultural land	0.9893	0.0025	0.0000	0.0082	0.0000
	Urban and built-up	0.0231	0.9152	0.0000	0.0617	0.0000
	Lakes	0.0000	0.0000	0.979	0.0209	0.0000
	Desert lands	0.0144	0.0046	0.0008	0.9791	0.0010
	Roads Networks	0.0003	0.0000	0.0000	0.0012	0.9985

Table 12. Matrix of transition values (number of Pixels) for different land use/land cover during the period 2020–2030.

	2030
Year	Criterion	Agricultural land	Urban and built-up	Lakes	Desert lands	Roads Networks
2020	Agricultural land	420457	1069	1	3476	1
	Urban and built-up	9279	367789	5	24786	15
	Lakes	7	2	97921	2088	0
	Desert lands	170122	54075	9753	11540393	12024
	Roads Networks	31	0	2	143	117831

3.4.1. Predicting the future of land uses/land cover in 2030

The CA-Markov model has been widely used in predicting land-use/land cover changes, the CA-Markov model chain generates a transition probability matrix, while cellular automation controls the evolution and changes in cells (Chang & Chang, 2006; Kamusoko et al., 2009), and it is clear from the Figures 19 and 20 and Table 13 that the land uses/land cover of the Al-Hassa Oasis in the year 2030 vary among themselves in terms of the amount of change expected for each type of different land uses/land cover, either increasing or decreasing it in the year 2020.

Figure 19. Map of future land uses/land cover of Al-Hassa Oasis: (a) in 2030, (b) in 2050 using CA-Markov within IDRISI selva v.17.0 program.

Figure 20. Losses and gains between land use/land cover from 2020 to 2030.

Table 13. The expected areas for each type of land use/land cover in 2030 for Al-Hassa Oasis compared to 2020.

Phenomenon	Areas in km^2
	2020	2030 all uses	Predicted change
Agricultural land	162.31	167.14	+4.83
Urban and built-up	166.62	217.97	+51.35
Lakes	20.77	24.37	+3.6
Desert lands	2379.6	2319.82	−59.78
Roads Networks	38.36	38.36	-
Sum	2767.66	2767.66	-

The increase in the urban mass came at the expense of the desert area, which, on the contrary, was the most likely to lose space, as the total area of 2319.82 km2 was recorded in 2030 with a decrease of 59.78 km2 compared to the year 2020 when the total area of the desert was 2379.6 km2. While the expected increase in agricultural lands in 2030 compared to what it was in 2020 amounted to 4.83 km2 to record a total area of 167.14 km2 instead of 162.31 km2 in 2020, and the expected increase in water for 2030 was 3.6 km2 to record a total area of 24.37 km2 instead of 20.77 km2 in 2020..

Where we find that the urban mass is the most likely type in terms of area increase, as the expected increase in the urban mass area in 2030 compared to what it was in 2020 is 51.35 km², to record a total area of 217.97 km² instead of 166.62 km² in 2020.

3.4.2. Predicting the future of land uses/land cover for the year 2050

The future land uses/land cover prediction model contributes to understanding future changes in land use and developing the right policies to support the growth of a sustainable city (Yu & Du, 2011). It is clear from Figures 19 and 21 and Table 14 that the land uses/land cover of the Al-Hassa Oasis in 2050 compared to 2020 vary with each other in terms of the amount of change expected for each type of land uses/land cover, either by increasing or decreasing.

Figure 21. Losses and gains between land use/land cover from 2020 to 2050.

Table 14. The expected areas for each type of land use/land cover in 2050 for Al-Hassa Oasis compared to 2020.

Phenomenon	Areas in km^2
	2020	2050 all uses	Predicted change
Agricultural land	162.31	167.93	+5.62
Urban and built-up	166.62	220.31	+53.69
Lakes	20.77	24.37	+3.6
Desert lands	2379.6	2316.69	− 62.91
Roads Networks	38.36	38.36	-
Sum	2767.66	2767.66	-

Where we find that the urban mass is the most likely type in terms of area increase, as the expected increase in the urban mass area for the year 2050 compared to what it was in 2020 is 53.69 km2, to record a total area of 220.31 km2 instead of 166.62 km2 in 2020, and the increase in the urban mass came at the expense of the desert areas, which, on the contrary, was the most likely to lose area, as it recorded a total area of 2316.69 km2 in 2050, with a decrease in area equal to 62.91 km2 compared to the year 2020, when the total area of desert areas was 2379.6 km2. While the expected increase in agricultural lands in 2050 compared to what it was in 2020 amounted to 5.62 km2, to record a total area of 167.93 km2 instead of 162.31 km2 in 2020, and the expected increase in water for 2050 was 3.6 km2 to record a total area of 24.37 km2 instead of 20.77 km2 in 2020.

4. Discussion

According to the authors’ knowledge and the reference inventory, this is the first study that aims to analyze the spatiotemporal dynamics of land use and urban growth and predict them in the future in the Al-Hassa Oasis in the Kingdom of Saudi Arabia. The results of this study can be the basis for several applications (Alqurashi & Kumar, 2014; Alsharif et al., 2022; Rahman, 2016) that have not been followed up so far due to the lack of spatiotemporal data over long historical periods (Abdelkarim et al., 2022). The study monitored and predicted spatial and temporal land use changes in Al-Hassa Oasis, Saudi Arabia, using land cover maps for the years 1988, 2000, 2013, 2020. LULC maps were classified into four main categories: (Water, Agriculture, Desert lands, and Desert) through the supervised classification method using the maximum likelihood classifier (MLC) method, Maximum likelihood classifier (MLC), it is a robust and popular supervised classification algorithm, to derive mapping rules for all other pixels into the class values (Dey et al., 2021; Kafy et al., 2021; Wang et al., 2021).

The maximum likelihood classifier (MLC) method was used as a supervised classification algorithm to super control classification coefficients. Training areas were extracted to classify green and water from Landsat images using NDVI and NDWI indicators, respectively. The classified image was compared with ground fact data to verify and calculate the accuracy of the classification. The accuracy confusion matrices showed that the overall accuracy was 91.23%, 93.65%, 95.32%, and 98.66% and the Kappa coefficients were 89.64, 92.41, 93.66, and 96.98 for the four Landsat images in 1989, 2002, 2013, and 2021, respectively. These statistics reflect that the different ranked maps meet the minimum recommended overall accuracy of 85%. They also show that there is a remarkable symmetry between the reference data and the land use/land cover classification according to the selected indices (Foody, 2002). These results are consistent with the results obtained by (Abdallah et al., 2019; Rahman et al., 2017; Salih, 2018).

The LCM tool was applied to analyze the LULC change of Al-Hassa Oasis during the period 1988–2020. The “Validate” tool was used in the Idrissi software to check the accuracy of the model, which ultimately showed very high accuracy for different Kappa coefficients. The results revealed significant changes in the map of land use for the Oasis – where the urban areas increased by 63.5%, equivalent to about 105.7 km² during that period, which is about 2.7 times its area in 1988 AD, from 60.89 km² in 1988 to about 166.6 km² in 2020. The urban increase was about 105.71 km², and the area of the increase was distributed over 32 years, with a growth rate of 3.2%. The largest percentage of these changes was at the expense of desert areas with 93.33% of their area. This significant increase can be explained by the fact that there is a positive relationship between urban growth and oil discovery in the Kingdom of Saudi Arabia between 1974 and 2019, especially in the eastern region, which includes the study area. This has been accompanied by major shifts in infrastructure and urbanization. Al-Hassa Oasis has gotten the largest share of laborers due to its proximity to the oil areas. This was confirmed by the population growth data in Al-Hassa Oasis, which increased to 61.9% in the official census of the General Authority for Statistics during the period between 1974 and 2020 (Saudi General Authority for Statistics [SGS], 2020), as the total population of the Kingdom of Saudi Arabia witnessed a huge increase by 500% between 1974 and 2019 (Saudi General Authority for Statistics, 2020). These families who came to Al-Hassa for work began to pump new population, and this population growth witnessed by Al-Hassa Oasis played a major role in the desert urbanization. There is no doubt that population growth plays an important role in expanding the urban area (Singh et al., 2015; Ye & Bai, 2008). It is expected that this urban transformation will continue in the composite land use map during the coming decades, and it is expected that the city’s population will reach 1,563,541 people by 2050 (SGS, 2020). This urban transformation in the land use map will not only lead to an increase in energy consumption for the population and an increase in burdens the new infrastructure and utilities, but it will also create more environmental challenges and problems, especially air, soil, and water pollution, constricting the agricultural lands in the east of the Oasis, lack of groundwater, well water depletion, rising wastewater, and soil salinity. Accordingly, sustainable desert urban life has been affected, and thus nature and resource flows are diminished. In similar studies conducted elsewhere, the change in land patterns and urban growth across desert areas led to many environmental challenges, such as water scarcity along with environmental degradation and groundwater shortage. Studies (Mansour et al., 2022). The agricultural land increased during the study period by 6.6%, equivalent to about 10.7 km². The largest percentage of these changes was at the expense of the desert areas around the four oases by 98.62%, in addition to the agricultural expansion towards the south of the Oasis. However, this increase is considered a deceptive and unreal increase in its entirety, as it is an increase in new places away from the old oasis, which has already lost about 3.3% of the agricultural land during the study period, specifically during the period 1988–2000. These results are consistent with the results obtained by (Abdelatti et al., 2017; Almadini et al., 2019; Alqarni et al., 2018). The land that turned into water marshes reached an area of 8.7 km², with an increase of about 41.7%. The largest percentage of these changes was at the expense of desert areas, with an increase in the Sabkhat Al-Asfar area by almost 100% northeast of the Oasis.

To predict the dynamics of LULC change for Al-Hassa Oasis for the years 2030 and 2050, we used the Cellular Automata-Markov Model (CA-Markov) and the Cellular Automata (CA). The model has been used efficiently in recent studies to simulate spatial and temporal changes in land uses due to its combination of dynamic spatial projections and probabilistic theories (Aneesha Satya et al., 2020; Huang et al., 2020; Wang & Zhang, 2001; Zheng et al., 2015); The projected LULC maps for 2030 and 2050 indicate that these patterns of arid lands that change into built-up areas and agricultural lands will continue over the next 30 years due to urban growth and the development of agricultural activities. We find that the urban mass is the most likely type in terms of the area increase, as the expected increase in the urban area reached 51.35 km² and 53.69 km², respectively, for 2030 & 2050 AD. This is in line with the expected population growth scenarios of the General Authority for Statistics during 2030–2050, while the expected increase in agricultural lands for the years 2030–2050 AD was 4.83 km² and 5.62 km².

The stage of classifying the satellite images of the Al-Hassa region faced major spectral problems related to: First, the overlap of the spectral range between the phenomena of different land uses, such as the overlap of the spectral range between agricultural lands and water, and the overlap of the spectral range between urban areas and some parts of the desert regions dispersing the spectral range. Second, the dispersion of the spectral range of one of the land use phenomena, such as desert lands, in which the spectral range is dispersed in a large extent due to the different geomorphological forms and surface geological compositions. These reasons made it difficult to derive the phenomena of land uses in separate categories.

Therefore, various indicators such as the NDVI and NDWI were used and multiple treatments were performed to separate the anomalous spectral fingerprints from the different categories and further improve the classification categories, which increases their accuracy and reliability. Many studies have encountered complexity in the process of classification and mapping in arid and semi-arid environments and have overcome these challenges by using various specified image approaches. These methods have mitigated errors that are wrongly classified and better distinguished between urban and arid areas as in the studies of (Abdallah et al., 2019; Alsharif et al., 2022). Another common problem found when mapping land cover for dry and semi-arid areas is the mixed pixel and spectral confusion between impervious surfaces and bare soil. This study also faced another challenge, as the study of the change detection in the phenomena of land use relies on a specific time period and a number of satellite visualizations suitable and selected for the study area during this time range in a way that achieves the best possible results. The longer the time span of the study and the greater the number of studied satellite visualizations, the more accurate the final results of the study will be. Satellite visuals of the study area were chosen with semi-rigid intervals starting from the available satellite visuals suitable for the study area in 1988 until the satellite visualization of the study year in 2020, with the exclusion of distorted or inappropriate ones, whether due to problems in the sensor itself, as in the ETM sensor, or problems related to weather conditions, as in 1990 and 2010, so they were excluded.

Even though the Cellular Automata-Markov Model (CA-Markov) and the Cellular Automata (CA) Model can be used as an effective tool for mapping and monitoring spatial and temporal changes and predicting LULC. Generally, the weaknesses inherent in these methods are that they do not take into account the environmental, social, or demographic conditions (Yang et al., 2022), for example, combining (CA-Markov), (MCE), and (AHP), which would improve the display of dynamic growth in urban areas and land use changes (Wang et al., 2022). An investigation into urban growth and land use change will be conducted later, focusing on combining the CA-Markov, Multi-Criteria Evaluation (MCE), and Analytical Hierarchy Process (AHP). Placing transition rules in the form of suitable maps based on the MCE allows us to consider different components commonly used in land use planning and decision-making such as changes in population density, terrain (slope), proximity to network routes, proximity to central business area water resources, restricted military zones, and oil. The Analytical Hierarchy Process (AHP) is an easy-to-use, powerful, and useful decision-making tool that allows collective decision makers to compare and select options as part of collective decision-making.

5. Conclusions

The study monitored and predicted spatial and temporal land use changes in Al-Hassa Oasis, Saudi Arabia using land cover maps for the years 1988, 2000, 2013, and 2020 by integrating the Cellular Automata-Markov Model (CA-Markov) and the Cellular Automata (CA). LULC maps were classified into four main categories (Water, Agriculture, Desert lands, and Desert). Through the supervised classification method using the Maximum Likelihood Classifier (MLC) method, the training areas were extracted to classify the green and water from Landsat images using NDVI and NDWI indices, respectively. A high classification accuracy result was achieved to ensure more reliability in the results. This was to understand the LULC changes in Al-Hassa Oasis comprehensively. The LCM tool was applied to analyze the LULC change of the Al-Hassa Oasis during the period 1988–2020. The “Validate” tool was used in the Idrissi software to check the accuracy of the model, which ultimately showed very high accuracy for different Kappa coefficients. The results revealed significant changes in the map of land uses for the Oasis – where the urban areas increased by 63.5%, equivalent to about 105.7 km² during that period, which is about 2.7 times its area in 1988 AD, from 60.89 km² in 1988 to about 166.6 km² in 2020. The urban increase was about 105.71 km², and the area of the increase was distributed over 32 years, with a growth rate of 3.2%. The largest percentage of these changes was at the expense of desert areas with 93.33% of their area. The agricultural lands increased during the study period by about 6.6%, equivalent to about 10.7 km². The largest percentage of these changes included the desert areas around the four oases by 98.62%, in addition to the agricultural expansion towards the south of the oasis. However, this increase is considered a deceptive and unreal increase in its entirety, as it is an increase in new places away from the old oasis, which already lost about 3.3% of agricultural land during the study period. The land that turned into water marshes reached an area of 8.7 km², with an increase of about 41.7%. The largest percentage of these changes was at the expense of desert areas, with an increase in the Sabkhat Al-Asfar area by almost 100% northeast of the Oasis. To predict the dynamics of LULC change for Al-Hassa Oasis for the years 2030 and 2050, we used the Cellular Automata-Markov Model (CA-Markov) and the Cellular Automata (CA). The model has been used efficiently in recent studies to simulate spatial and temporal changes in land uses due to its combination of dynamic spatial projections and probabilistic theories. The projected LULC maps for 2030 and 2050 indicate that these patterns of arid lands that change into built-up areas and agricultural lands will continue over the next 30 years due to urban growth and the development of agricultural activities. We find that the urban mass is the most likely type in terms of the area increase, as the expected increase in the urban area reached in 2030 & 2050 AD51.35 km² and 53.69 km², respectively. This is in line with the expected population growth scenarios of the General Authority for Statistics during 2030–2050, while the expected increase in agricultural lands for the years 2030–2050 AD was 4.83 km² and 5.62 km². Environmental sustainability is one of the necessary requirements in the map of future urban land uses in 2030 AD and 2050 AD in Al-Hassa Oasis. This will only be achieved through preserving agricultural lands east of the Oasis, controlling urban expansions, preserving valleys and environmental rehabilitation for them, benefiting from their water and resources, and stopping the drilling of wells randomly in a way that contributes to the preservation of water, food, and housing security in order to establish the rules for sustainable urban development and reduce the negative environmental impacts. The results obtained from this study provide information to assist decision makers to undertake effective practices for future planning and management of urban land growth, especially in the vision of the Kingdom of Saudi Arabia 2030 for various sectors.

Disclosure statement

No potential conflict of interest was reported by the author.