Urban growth modelling in Qazvin, Iran: an investigation into the performance of three ANFIS methods

ABSTRACT Determining the rules that lead to the expansion of urban areas has always been a challenging factor in urban modeling. To overcome this issue, an ANFIS model is proposed to enhance the simulation of urban growth through automatic production of transition rules. Hence, 22 ANFIS models were trained using different division methods and a Cellular Automata-based Markov Chain (CA-MC) was developed to examine their efficiencies. The results of accuracies reveal that the ANFIS accompanied by subtractive clustering (ANFIS-SC) simulating urban growth with a Kappa coefficient of 0.76 and overall accuracy of 93.41% is superior to other ANFIS and CA-MC models.


Introduction
The proportion of the world's population that resides in urban areas has grown from approximately 5% to over 50% in recent centuries (Ezeh et al. 2017), and it is expected that this proportion will rise to 70%, or even more, by 2050 (Kabisch and Haase 2011, Broere 2016, Hegazy et al. 2017).This issue causes a more serious problem to emerge, that is, the accelerating rate of urban expansion in most parts of the world, insofar as it is predicted that the world's urban areas will have tripled between 2000and 2030(McDonnell and Hahs 2015, Southon et al. 2017, Khamis et al. 2018).Rapid urban expansion has led to a variety of complex problems, such as increasing informal settlements, air pollution, deforestation, soil erosion, and climate change (Kamusoko and Gamba 2015, Xu et al. 2019, Wang et al. 2020).Additionally, rapid growth of cities has strained their capacity to provide services, such as education, health care, and transportation (Ewing et al. 2002).Moreover, the level of urbanisation in Iran increased from 33.7% in 1960to 74.9% in 2018(World Bank 2018), and this high rate of urbanisation has led to the fast expansion of built-up areas in recent years (Bihamta et al. 2015, Taravat et al. 2017).Such issues can be prevented from occurring by modelling and studying future urban expansion (Xu et al. 2019).
Over the last three decades, many models have been developed to simulate urban growth processes and land use/land cover (LULC) changes (van Vliet et al. 2016, Karimi et al. 2019).Among these models, cellular automata (CA) is a popular (Liu et al. 2017, Musa et al. 2017) and strong technique for simulating the pattern of built-up areas (Aburas et al. 2016, Mustafa et al. 2018).Nevertheless, determining the transition rules, which leads to the calibration of the CA model, has always been a key challenge in this model because of the complex relationship between urban expansion and its multiple driving factors (Blecic et al. 2015, Berberoğlu et al. 2016, Li et al. 2017).To handle this challenge, several calibration techniques have been operated alongside CA model, such as multicriteria evaluation (MCE) (Wu and Webster 1998), spatial autoregressive model (SAR) (Feng et al. 2021), and logistic regression (LR) (Wu 2002).More recently, some evolutionary algorithms such as particle swarm optimisation (PSO) (Wang et al. 2021) were also utilised.Although these methods overcome the calibration challenge of CA, they also introduce some drawbacks into the model.First, in some statistical models such as regression models, a fitness function (linear or nonlinear) must be determined and applied for mathematical modelling of observed data.Therefore, determining the fitness function and its performance in dealing with complex and spatial urban systems are challenging constraints.Second, knowledge-based methods such as MCE may incorporate errors due to human judgements into weighting and consequently model results.In contrast, data-driven methods reduce the dependence on experts and are more reliable on condition that they utilise valid data.Third, modelling nonlinear and complex relationships in urban growth phenomenon is too difficult for models such as LR due to low modelling capabilities (Triantakonstantis and Mountrakis 2012).And fourth, some methods such as PSO can easily be involved in local optimum in high-dimensional spaces (Li et al. 2014).This is due to lack of a solution such as validation data that evaluates error function and prevents overfitting of model function to data.In recent years, some machine learning methods have been mainly able to overcome these problems.For instance, Artificial Neural Networks (ANN) have been extensively utilised in urban growth simulation (Li and Yeh 2001, 2002, Li et al. 2013) and have overcome the local optimum problem, the human wrong judgements, and spatial, nonlinear, and complex modelling.Nevertheless, neural networks have some disadvantages, such as difficult calibration (Kamusoko and Gamba 2015) and incompatibility with the uncertainty of transition rules (Azari et al. 2016).Azari et al. (2016) utilised fuzzy set theory to tackle the uncertainty of transition rules.However, fuzzy models require a justifiable approach for calibrating parameters of MFs (Al-Ahmadi et al. 2009).
Several works have demonstrated that the weaknesses of fuzzy sets and ANN can be overcome by using an Adaptive Network-based Fuzzy Inference System (ANFIS).An ANFIS is a hybrid intelligent system that carries the benefits and learning potential of neural networks and fuzzy logic (Parvinnezhad et al. 2019), and eliminates many limitations of these two algorithms (Ahmadlou et al. 2019).Mohammady et al. (2013) utilised ANFIS for urban growth modelling from 2000 to 2006 in Sanandaj, Iran.They employed figure of merit and percent correct match for model evaluation, which obtained 64% and 94%, respectively.Mohammady and Delavar (2016) and Mohammady (2016) have demonstrated that the ANFIS model has more accurate results in the simulation of urban growth compared with the ANN and LR.
The ANFIS model employs and trains a Fuzzy Inference System (FIS) with distinct division methods, such as Grid Partitioning (GP), Subtractive Clustering (SC), and Fuzzy C-means (FCM) to divide the input data space.It has been reported in some scientific studies that clustering type can play an important role in ANFIS model performance, such as Adedeji et al. (2020) and Moradi et al. (2019), whereas there is a paucity of high-quality research investigating this subject in urban development modelling.It is explicit that calibrating various kinds of ANFIS with urban parameters, comparing the accuracy of trained models, and finally, employing the most well-trained ANFIS model, leads to more precise transition rules and, as a result, more accurate forecasting of LULC changes.This study, therefore, aims to implement ANFIS accompanied by GP (ANFIS-GP), SC (ANFIS-SC), and FCM (ANFIS-FCM) models to assess the effect of clustering methods and their parameters on the accuracy of ANFIS in urban growth modelling.Besides, for urban planners around the world, it intends to provide a potent model that overcomes uncertainty over complex urban systems with no need for expert knowledge.Additionally, it aims to simulate future built-up areas and warn about the possible harmful effects of spontaneous urban growth on the environment.To fulfill these aims, we executed multiple scenarios for the ANFIS training to achieve the most accurately calibrated model.Additionally, we developed a popular (Wang et al. 2020) existing model, namely CA-MC (Rimal et al. 2017, 2018, R. Wang et al. 2018), to provide a benchmark against which to interpret the ANFIS outcomes.The reason why we chose CA-MC over other methods is that this algorithm has been extensively utilised in urban growth modelling and has had relatively acceptable simulations (Nugroho and Al-Sanjary 2018).Afterwards, we compared the simulation results of CA-MC with the most acceptable ANFIS model for the 2019 simulated map.We finally simulated the upcoming urban expansion of Qazvin city for the year 2029, utilising the ANFIS model, for the use of urban planners and policymakers.

Study area
Qazvin Province, with a top urbanisation rate of 73.1%, is located in the northern half of Iran, as one of the 31 provinces of the country.Qazvin city, as the central city of this province, lies at an altitude of 1278 m (4193 ft.) above sea level.The study area of this research includes Qazvin city and its suburbs, located between 36.2266° and 36.3604°N, 49.9614° and 50.0555°E (128.63 km 2 ) (Figure 1).The shapefiles and satellite image presented in Figure 1 were acquired from GIS bureau of Qazvin municipality and Landsat 8 satellite, respectively.
The most significant reason for choosing this city as our study area was that the city is suffering from rapid growth in both population and urban areas.According to the census of the Statistical Center of Iran, the population of Qazvin city was 402,748 in 2016, whereas in 2006this number had been 349,821 (Statistical Center of Iran 2019).This increase in population has led to rapid outward city growth.Low population-density zones and some formal and informal settlements have appeared surrounding the city in the last decades.Emergence of these areas causes some serious damage to the environment and natural resources.Consequently, an investigation of the rapid growth of this city appears to be vital for urban planners.

ANFIS architecture
In fuzzy logic, the MF specifies the grade of membership of each element (membership degree) to a fuzzy set.However, determining the optimal location and figure of an MF in a fuzzy inference system is a challenging problem.Jang (1993) developed an integrated neuro-fuzzy model which adjusts MF parameters through a learning process like ANN algorithm, consequently determining the optimal location of MFs.To put it more simply, this model employs a fuzzy inference system and trains it by a specific training algorithm.In fact, ANFIS carries the advantages of both ANN in terms of training and fuzzy reasoning in terms of using FIS and Linguistic terms (such as: young, middle-aged, old, . . .for age variable) at the same time (Zhou et al. 2017).This model is composed of five layers (Fig. S1 in Supplementary Materials).
The parameters of the first and fourth layers are recognised as premise and consequence, respectively.Unlike the parameters of the other layers, they are adjusted throughout the training process.The first layer has come to be known as fuzzification layer.Nodes in this layer take the numeric inputs and express the membership grade of the inputs set as the output by fuzzy Linguistic terms.For example, for the input x and fuzzy set Ai, the output of the first layer is defined as follows.
where μ A i x ð Þ represents the membership grade of the input x in fuzzy set Ai.The second layer is known as the rule layer.The general structure of a fuzzy If-Then rule in the ANFIS network is as follows.
where x and y are ANFIS inputs, A and B are fuzzy sets, Fi is an output function, and pi, q i , and r i are consequence parameters belonging to the output function.Any T-norm (a function that can be utilised for interpretation of the intersection of fuzzy sets in fuzzy logic) that performs the AND part of fuzzy rule can be deployed in the second layer as the fuzzy implication function.Hence, this layer utilises the minimum T-norm and multiplies the membership grades coming from the first layer in order to obtain the weight of each rule.Eq.3 represents the functionality of this layer for two x and y inputs.
where w i is the weight of ith rule and μ A i x ð Þ and μ B i y ð Þ represent the membership grades of the x and y.The third layer computes the normalised weights of ith rule ( � w i ) using Eq.4.
The fourth layer, which is denoted as the defuzzification layer, multiplies the normalised weights from the previous layer by Then part of fuzzy If-Then rule (output function).Eq.5 exhibits the operation of this layer.
The only node in the fifth layer finally adds the outputs of the previous layer and generates the conclusive output of ANFIS (Karaboga and Kaya 2019).For more information on the structure and functionality of ANFIS, please refer to Jang (1993).

Markov chain analysis
Markov chain analysis (MCA) describes a sequence of events, where the probability of an event occurring at a time interval influenced only by its state at the previous time interval (Fan et al. 2008).Generally, according to the Markov chain theory, if the sequence of the random variables of x n takes the distinct values of a 1 , a 2 , . .., a n , then the Markov chain can be represented by Eq.6 (Fan et al. 2008).
With the assumption that the values of the variables were changed from time t to time t-1, the transition probability in time t can be computed by the following relation (Balzter 2000).
Term p in Eq.7, which is named the transition probability matrix, represents the probabilities of transitioning.Assuming that the states are 1, . .., m, the p can be displayed as the following matrix.
where the p 1m represents the probability of transition from state 1 to state m.As a stochastic model, the Markov chain can express the probability of the LULC changing in simulating land use changes (Wu et al. 2019); however, it cannot obtain the spatial distribution of projected LULC (Nguyen et al. 2019).In order to resolve this issue, the CA-MC model, which integrates the CA with the MCA model, has been employed in various urban growth and LULC changes studies.This integrated model includes a spatial dimension considering the effects of neighbouring cells by CA neighbourhood.

Methodology
Satellite imageries are qualified primary data sources for studying LULC changes.Hence, for this research, we adopted Landsat images (from 2000, 2010, and 2019) due to their accessibility and appropriate spatial resolution.Figure 2 provides a general overview of the principal performed steps of this study.Image pre-processing was performed on the data because the received images were accompanied by errors.Afterwards, images were classified and the six driving factors of urban growth as input layers of the ANFIS and CA-MC models were generated.This study developed and implemented a two-step model that comprised two phases of calibration for 2000 to 2010 and prediction for 2010 to 2019.The calibration phase consisted of three parts.In the first and second parts, ANFIS-GP, ANFIS-SC, and ANFIS-FCM models were calibrated and compared in terms of performance and accuracy of training in pairs using Mean Squared Error (MSE) and Root Mean Square Error (RMSE), and the most appropriate model was selected.In the third part, a CA-MC model with WLC for the generation of transition suitability images and AHP for factor weighting was developed.In the prediction phase, the built-up areas of 2019 were predicted using the CA-MC and most well-trained ANFIS models.Eventually, the predicted map of CA-MC was evaluated using overall accuracy, Kstandard (Kappa standard), Kno (Kappa for no information), Klocation (Kappa for location), and KlocationStra (Kappa for stratum-level location), and the predicted map of ANFIS was assessed utilising overall accuracy and Kappa.

Data sources and imageries preparation
Landsat images have been largely used for urban planning due to aspects of appropriate spatial resolution and shorter revisit time (S.Wang et al. 2018).Hence, for this research, we adopted the Landsat 7 ETM+ (Enhanced Thematic Mapper) images on 22 May 2000, and 3 June 2010, and Landsat 8 OLI/TIRS (Operational Land Imager/Thermal Infrared Sensor) image on 20 June 2019.The images were downloaded from the United States Geological Survey (USGS) database (https://earthexplorer.usgs.gov/)as the primary data sources of this study.As each sensor had a distinct spatial resolution, all images were resampled to 30 m resolution, and then projected to the Universal Transverse Mercator (UTM) 39 N coordinate system.
Before utilising the satellite images, pre-processing steps must be performed on the imageries.Pre-processing refers to those processes that apply to the raw satellite images to reduce or eliminate the errors of sensor and platform and includes radiometric and geometric corrections (Asubonteng et al. 2018).The obtained images have georeferenced and orthorectified with high geometric accuracy (fewer than 6.8 m RMSE), and no additional geometric corrections were required.
Atmospheric correction, as one of the most important steps in pre-processing, eliminates the influences of atmospheric phenomena, such as absorption and scattering (Jha et al. 2021).Hence, we applied an absolute atmospheric correction to all imageries by using the fast line-of-sight atmospheric analysis of hypercubes (FLAASH) algorithm in the ENVI platform.
The Scan Line Corrector (SLC) of Landsat 7 ETM+ has failed since 31 May 2003, which made a loss of about 22% of data in the SLC-off images (Masoumi et al. 2017).Therefore, in order to fill the gap areas of the image of 2010, we used the Local Linear Histogram Matching (LLHM) approach by adding an add-on to the ENVI software and adopting another Landsat 7 image (image on 19 June 2010) from the same area of the initial image.
We selected the suitable three-band combination by applying the Optimum Index Factor (OIF) technique (Table S1 in Supplementary Materials).Afterwards, we carried out a supervised method, namely maximum likelihood, to classify the imageries in the ArcMap software.Urban, agricultural, barren, and transport were selected as the land use classes based on study objectives.For assessing the classification performance, we collected the land use of random sample points by GPS and generated a confusion matrix.An overall accuracy of 80% and a kappa coefficient of 0.72 were received for the classification.According to Pijanowski et al. (2005), the accuracy of classification was acceptable.

Driving factors and input layers
Up to this time, various driving factors of urban expansion, such as distance to primary and/or secondary roads, distance to previous urban, slope, elevation, population, and access to jobs have been considered in diverse investigations (Pijanowski et al. 2014, Mustafa et al. 2017, Shafizadeh-Moghadam et al. 2017a, 2017b, Radwan et al. 2019).In determining driving variables, it is highly important to pay attention to the particular characteristics of the local area.Qazvin, as our study area, possesses several colleges and universities at the present.Most of these institutions have been built outside of the city limits, which makes the value of nearby land grow and the chance of LULC change increase.Additionally, past research has not focused on the distance from educational institutions as a driving force in urban growth or sprawl modelling (Santé et al. 2010, Aburas et al. 2016, Hosseini and Hajilou 2019), except for a very limited number of studies such as Triantakonstantis et al. (2014).Hence, distance to the universities by five other spatial variables, namely slope, aspect, distance to the urban areas, distance to the urban transportation network, and distance to the urban squares, were selected as the driving factors of this study.
Accordingly, six raster layers related to driving factors were created in the ArcMap environment.Two layers of slope and aspect were achieved using Aster Digital Elevation Model (DEM) v.3 (acquired from NASA data resources: https://earthdata.nasa.gov/).Four distance raster layers were generated by Euclidean distance ArcMap tool using shapefiles produced by Google Earth and LULC map.Furthermore, the layer of areas with unchangeable LULC, such as military bases and parks, were generated by the LULC map.The value of these pixels remained 0 (unchangeable land use) in all study stages.All accomplished raster layers for 2019 with 30 m pixel resolution and UTM 39 N coordinate system are presented in Fig. S2 (Supplementary Materials).
Since pixels in every layer had a distinct range of quantities and units, we normalised all of them in the range of 0 to 1 by Eq.8 in the MATLAB (version 2016b) software.
where X(i, j) is the normalised value of each pixel, min(X(i, j) is the minimum value of all pixels in a specific layer, and max(X(i, j) is the maximum value of them.

ANFIS models implementation
ANFIS can be categorised into four types, namely ANFIS-GP, ANFIS-SC, ANFIS-FCM, and ANFIS accompanied by a context-based FCM (ANFIS-CFCM) according to the division method of inputs (Yeom and Kwak 2018).The division methods (such as GP and SC) cluster the input data to fuzzy clusters based on pixel attributes, regardless of pixel regional positions.Afterwards, an MF is allocated to each cluster and a FIS is generated (Adeleke et al. 2020).As the common ANFIS methods, we implemented ANFIS-GP, ANFIS-SC, and ANFIS-FCM using GenFIS1 (generates initial FIS), GenFIS2, and GenFIS3 MATLAB functions, respectively.All models were accomplished with various input parameters and Gaussian MF.Fig. S3 (Supplementary Materials), which achieved by ANFIS editor Graphical User Interface (GUI), represents our ANFIS-SC structure with the six driving factors of this study as inputs of the model and one output representing the development probability of each pixel.
The number of generated rules in ANFIS-GP is equal to the number of MFs to the power of the number of input layers, whereas the number of generated rules in ANFIS-SC and ANFIS-FCM is equal to the number of MFs.Hence, if we had only five MFs, GenFIS1 would have generated 15,625 fuzzy if-then rules with our six input layers.It takes weeks or months for a common computer to run an ANFIS model with this amount of rules.Therefore, we calibrated both ANFIS-GP and ANFIS-SC in a separate part decreasing the number of MFs (for ANFIS-GP), input layers, and training dataset to 3, 4, and 50,000, respectively.Consequently, the ANFIS-GP model was implemented and the ability to compare the results of these two models was provided.Due to the better performance of ANFIS-SC in the first part, we implemented ANFIS-SC and ANFIS-FCM as the second part of the calibration.We accomplished ten ANFIS-SC and eight ANFIS-FCM considering six input layers, zero error tolerance, and 200 epochs in ANFIS Architecture.It is worth mentioning that we examined various radii of influence for the SC method in ANFIS-SC and the number of generated clusters for FCM clustering in ANFIS-FCM using both hybrid and back-propagation optimisation methods.

Implementation and validation of CA-MC model
To implement the CA-MC model, we accomplished the MCA model in the Idrisi software by introducing land use maps of 2000 and 2010, and two time intervals from 2000 to 2010 and 2010 to 2019 into the model.As a result, the transition probability matrix was obtained.One of the inputs of CA-MC is the transition suitability image that reveals the suitability of converting to another land use.We used weighted linear combination (WLC) as one of the MCE methods to obtain these maps.To achieve this aim, first, we obtained transition potential images (express the potential of change in the specific time) for three urban, agricultural, and barren land uses and then converted all three maps to transition suitability images in Idrisi.It is worth mentioning that we normalised all pixels between 0 and 255 by fuzzy reasoning due to differences in units and pixel values, and utilised the AHP for factor weighting in the MCE method.
We assessed the CA-MC model performance by comparing the output map of the model and the land use map of 2019 and obtaining the kappa coefficients in Idrisi.To fulfill this aim, we computed various Kappa coefficients, namely Kstandard, Kno, Klocation, and KlocationStra to validate the model in both quantification and location.For a detailed description of these coefficients refer to Supplementary Materials.

Results and discussion
All input parameters and execution results of ANFIS-GP and ANFIS-SC models, executed in the first part of the calibration section, are provided in Table 1.According to this table, although the RMSE of ANFIS-GP was approximately equal to the RMSE of ANFIS-SC, its MSE was greater than the MSE of ANFIS-SC.Besides, it is apparent that the execution time of the ANFIS-GP was longer than ANFIS-SC due to a larger number of generated rules.On the other hand, ANFIS-SC was able to approximately obtain an equal accuracy to ANFIS-GP with more MFs and fewer rules.More MFs and fewer generated rules lead to more coverage of the input space and less execution time of the model, respectively.
Additionally, the execution results and characteristics of ANFIS-SC and ANFIS-FCM in the second part are presented in Tables 2 and 3, respectively.Maximum execution time of 50,000 seconds was considered for all ANFIS models due to long execution time of some programs.All models were successfully executed within the time limit, except for the ANFIS-SC with the radius of 0.1.By comparing all RMSE and MSE errors of the training and test data, it was found that the ANFIS-SC with a radius of influence of 0.2 (the hybrid one) was the most accurate for training data.The ANFIS-SC with a radius of 0.4 had the least MSE for test data; however, it had more RMSE and MSE for training data in comparison with the ANFIS-SC with a radius of 0.2.The results generally indicated the superiority of the ANFIS-SC to the ANFIS-FCM and ANFIS-GP in urban growth modelling.Although we are unable to compare this finding with similar studies due to lack of urban development studies in this field, it is consistent with the results of Adedeji et al. (2020), Moradi et al. (2019), andAkkaya (2016) studies accomplished in different areas.However, some studies such as Adeleke et al. (2020) (a study on modelling waste generation) have represented ANFIS-GP as a more precise model.This contradiction could have arisen because of differences in the inputs, MFs type and number, and ANFIS training quality.
The line charts in Figure 3 reveal the relation between the RMSE error of model training and the number of clusters/centre's range of influence in both back-propagation and hybrid methods.From these line charts it can be observed that the accuracy of training in ANFIS-FCM has not necessarily increased with the increase of clusters; however, the accuracy of ANFIS-SC has increased with decreasing radius (increasing clusters).Furthermore, it is quite obvious that the models that have used the hybrid optimisation algorithm have always had less error than the back-propagation algorithm.Hence, utilising an ANFIS-SC with a hybrid algorithm and smaller radius of influence can lead to more accurate modelling of urban expansion.
As it is indicated in Figure 4, the target-output plot, error plot, and error histogram chart of the ANFIS-SC with the radius of 0.2 were generated in order to further examine the outputs of the ANFIS model.The target-output plot represents the ideal output of ANFIS (0 and 1) as target and the predicted output as output at the same time, whereas the error plot reveals the difference between them, and the error histogram chart illustrates the distribution of errors.Looking at the error histogram chart, it is apparent that more errors were close to zero, which means the predicted outputs of the ANFIS-SC have been close to the ideal output.Additionally, the target-output and error plots reveal that the ANFIS had an appropriate recognition of land use changes due to the proximity of outputs to targets.Consequently, the results of presented ANFIS have been satisfactory.According to the results and implementation of proposed model, the major model drawbacks and limitations are: 1 -Accuracy of ANFIS-SC depends largely on the radius of clusters.The accuracy and run-time increases with decreasing radius, and a smaller radius would mean more MFs and heavier computations.Hence, we need for computers with more powerful hardware to run more precise ANFIS models.Additionally, type of MFs is a challenging part in ANFIS implementation.2 -We are not able to utilise all driving factors influencing the expansion of the city due to the lack of data on some factors like historical prices of land parcels in our study area.To our knowledge, in developing countries such as Iran, there is no comprehensive historical urban data.3 -The use of high-resolution DEMs and satellite imageries enhance the accuracy of classification and  simulation, however the images of some high-resolution sensors are not available in our study area.On the other hand, key strengths of proposed model are: 1 -ANFIS model overcomes observed uncertainty in spatial phenomena such as urban expansion, due to the use of interpretable Linguistic terms.Therefore, ANFIS can have precise results not only in urban development but also in other spatial phenomena.2 -ANFIS operates well in dealing with nonlinear problems.Therefore, as it is clear from the results, it has been well executed in modelling urban expansion as a nonlinear problem. 3 -There is an accelerating rate of urban expansion in developing countries like Iran.In such cases, the CA neighbourhood can be a barrier to accurate forecasting.Our model overcomes this limitation of CA and has good results in dealing with rapid or sprawl expansions, especially in developing countries.
The predicted map of the ANFIS-SC (radius of 0.2) model for 2019 had a Kappa coefficient of 0.76 and overall accuracy of 93.41, whereas the CA-MC model had a Kappa of 0.66 and overall accuracy of 77.96.The kno of 0.66, Klocation of 0.79, and KlocationStra of 0.79 for the CA-MC model indicated that this model can have satisfactory results; however, the ANFIS model had much better results than CA-MC due to better Kappa and overall accuracy.Besides, a Kappa between 0.6 and 0.8 is considered a very good kappa (Pijanowski et al. 2005).Hence, as it is presented in Figure 5, the predicted map for 2029 was achieved by ANFIS-SC model.Comparing this map with the 2019 land use map, it is apparent that the projected development will occur in the north more than in other directions.A possible explanation for this may be the existence of two large universities by enrolment in the north of the area in addition to the suitability of other driving factors.These universities with more than 36,000 enrolled students cause a large number of trips to the north of the city during the day.The examination of the road land use changes over the last few years reveals that urban decision-makers have improved the road transport network in the northern areas to meet the requirements of the students.The built-up roads to the universities have improved access into the northern barren lands; consequently, they have increased the probability of LULC change of these lands.Additionally, if we eliminate the factor of proximity to universities in our study, other changeable factors including proximity to roads, urban areas, and squares would gain greater weight in training of the model.Consequently, the forecast would be associated with more errors.Hence, large educational institutions can generally play an important role in increasing residential land use in the region, especially in small and medium-sized cities.This result accords with Triantakonstantis et al. (2014) observations, which showed that proximity to universities is a reason for rapid urban growth.They also represented that urban development at a distance of fewer than 500 metres to universities has increased by 360%.Therefore, this factor can be considered more in future studies on urban development.

Conclusion
This study was established to develop a powerful model based on ANFIS, and to investigate the role of clustering method on the accuracy of ANFIS model training and simulation power in urban expansion modelling.The most obvious finding to emerge from this study was that the proposed ANFIS model had a high accuracy in modelling the complex phenomenon of urban expansion and was able to be an acceptable alternative to the traditional model of cellular automata.Besides, according to the results of this study, the ANFIS model offered a more precise simulation of LULC in comparison with the CA-MC model.According to the results, we can infer that the ANFIS-SC model will lead to more precise urban expansion modelling in comparison with the ANFIS-FCM.Moreover, although the ANFIS-GP has precise and acceptable results, it is not a great option for modelling urban expansion due to longer execution time and inability to execute and high complexity in handling high-dimensional problems (Adedeji et al. 2020).
An important aim of the present research is to warn about the harmful effects of spontaneous urban expansion on the environment, using a functional model and predicting the land use changes from a ten-year perspective.The valuable finding of our study is that city expansion will only occur to the north towards the barren lands.The sprawl development in these barren lands can pose a significant threat to the environment of this area that is a habitat for various animals.Contrary to expectations, the predicted map of 2029 indicates that agricultural lands that are located in the south, west, and southwest of the city will not be damaged in the future due to land use changes.This finding is likely to be attributed to the lower price of barren land compared with agricultural land, the fertility of agricultural land, and supportive government policies protecting agricultural land, and it will be of broad use to policymakers.
Urban planners and decision-makers need to pay more attention to sustainable urban development in poor and developing countries due to increasing trend of urbanisation.To overcome this challenge, more powerful urban models would be an efficient solution for planners.According to the points raised in the discussion, the presented model is able to overcome some shortcomings of CA, ANN, and fuzzy logic and has the capability of expressing uncertainty in urban growth.It is also a qualified option for studying various types of development, such as sprawl, edge-expansion, and outlying in all growing cities with rapid growth rates.Therefore, the results of this study can pave the way for urban planners to implement effective support policies to control the rapid growth of cities around the world.
ANFIS models were implemented based on cellular spaces in this study.Hence, it might be possible to use a parcel-based land use classification in further investigation, and instead of using cells in the model, parcels can be used as the smallest unit for studying land use changes.Additionally, the centre's range of influence between 0.2 and 0.8 in ANFIS-SC was investigated due to the mentioned limitations, whereas it revealed that the model results improved with decreasing radius of clusters.Hence, further studies need to be carried out to achieve more accurate training of the model using a radius smaller than 0.2, and cautiously investigate the relationship between decreasing radius and increasing model accuracy.Moreover, a further study on implementing ANFIS-CFCM in urban growth modelling and comparing the results with the other three ANFIS models are suggested.

Figure 1 .
Figure 1.The map and satellite image of the study area: (a) location of Qazvin province in Iran, (b) location of Qazvin city in the county, and (c) satellite image of Qazvin city.

Figure 2 .
Figure 2. Principal performed steps of the study.

Table 3 .
Clustering method Centre's range of influence Number of MFs Number of rules Optimisation method Execution time (second) Execution results of ANFIS-FCM in the second part of calibration.Clustering method Number of clusters Number of MFs Number of rules Optimisation method Execution time (second)

Figure 3 .
Figure 3. RMSE error based on the number and radius of the clusters: (a) for ANFIS-FCM and (b) for ANFIS-SC.

Figure 4 .
Figure 4. Target-output plot, errors plot, and the error histogram chart of ANFIS-SC: (a) target-output plot, (b) errors plot, and (c) error histogram chart.

Table 1 .
Execution results of ANFIS-GP and ANFIS-SC in the first part of calibration.

Table 2 .
Execution results of ANFIS-SC in the second part of calibration.