Model Builder Project

Model Builder Project: GIS-MCDA for Tourism Center Site Selection in Levis, Quebec

Overview: In this project, I used a form of GIS-Multi Criteria Decision Analysis (MCDA) to derive a layer of suitability for a potential new tourism center in Levis, Quebec.

1. Introduction

1.1 Context

Determining optimal locations for development is an essential query conducted by business and real estate professionals readily, enabled by the power of geographic information systems. In this case, optimal site selection constitutes a form of multicriteria decision analysis (MCDA), where criteria are evaluated, normalised, and aggregated into a final model of suitability, and reflected “most-optimal” locations are indicated as a result (Steele et al. 2009). Thus, this project aimed to answer the following question: where are the most central locations to existing tourism and landmark sites? This question was followed as to source a new location for a tourist centre in Levis, Quebec. Notably, MCDA does have history for land use selection in other fields, inter alia healthcare (Muhlbacher and Kaczynski 2016) and urban development (Fernandes et al. 2018); however, limited work has been done for tourism site selection (see Sahani 2019 for example).

In this instance, existing point layers for notable sites were collected from the Quebec Open Data archive. These layers included a series of vectorized data (e.g., points, lines and polygons), and in efforts to work with both data types (as to include raster), this project used the vectorized multipoint data to create raster layers. The context for this project is that tourism locations have been, and continue to be, centralised or most proximal to points of interaction (Elwin 1989; Eldrandaly and Al-Amari 2014). Thus, locating a centralised hub is of importance to tourism boards, and effectively making tourism programs that allow visitation, within reasonable bounds of how far one might travel, to possible sites of interest.

1.2 Scenario and Questions

In this scenario, the tourism board of Levis was looking to rate, off of existing tourism and recreation sites, where optimal areas for the new tourism center would be. Following a land suitability assessment, the following was queried: Where is the most central locations for a new tourist center, relative to existing tourism and landmark sites?

2.1. Source Data

The open source data collected and used included 1) the locations of campgrounds, 2) Historical landmarks, 3) cemeteries, 4) reservoirs and significant bodies of water, 5) protected wetland locations and 6) picnic sites. All source data layers were located in Levis, Quebec and included in UTM 1983, Zone 19. This projection was selected to reduce geographic distortion as an equal area projection for this geography, a priori to precedent set in the peer-reviewed literature (Fraser et al. 2017)

The data was all point data, and resulting raster layers of euclidean distance were computed for each area, with bounds set to the display extent of the map interface. Euclidean distance calculates a weighted distance from a point, in reference to all points; the result is a matrix of continuous scores of proximity or clustering of possible locations (Danielsson 1980). Because proximity to points is relevant to determining potential sites, no issues in spatial autocorrelation were tested, this belief is a priori to the literature, and followed herein (Behrens et al. 2018). The resulting euclidean distance layers were raster files. Once reclassified to normalised ranges of 0-100, another raster layer existed for each location. After an equal weighted addition of all layers, a final weighted overlay layer was presented. In total, there were six (6) point data sets, and thirteen (13) raster layers used in this analysis.

2.2 Criteria Selection

Criteria selection was determined on the basis of data availability, the number of six sites was used per the assignment requirements of iterating in a function for size (6) instances. In this case, these were the source layers in a euclidean distance calculation and aggregation procedure. Criteria had to be point datasets with more than 1 point to calculate euclidean distances. For example, early criteria with one location were not able to be used in these calculations and returned error messages.

Criteria were evaluated through the euclidean distance function, where the point data were considered the source values and cells increasingly far from the source points were rated as depreciating. The procedure for one layer was imported into model-builder, then edited in python to loop through the length of all input layers (6). This could be automated for any amount of procedures, as long as the file names are specified by the reader for the output. The loop was written with index lists of the two to iterate between the same number of both lists (i.e., when working on the 3rd point layer, it would be saving the euclidean distance layer as the 3rd output name, not the forth or another number).

2.3 Model Builder

The model builder function included three major steps. Firstly, euclidean distance was computed for all point vector layers (6 layers). Secondly, the values were normalized into a 0-100 score in raster calculator. The reason for this was to try and automate the raster calculator procedure, as opposed to following the reclassify procedure where values are truncated; the raster calculator procedure preserves the distribution of the data. Thirdly, normalized distance layers were added together in a weighted overlay function. Equal weights were assigned to all 6 layers (16.677%). Herein, detail is provided on each of these procedures.

2.3.1 Euclidean Distance

Criteria were evaluated through the euclidean distance function, where the point data were considered the source values and cells increasingly far from the source points were rated as depreciating. The procedure for one layer was imported into model-builder, then edited in python to loop through the length of all input layers (6). This could be automated for any amount of procedures, as long as the file names are specified by the reader for the output. The loop was written with index lists of the two to iterate between the same number of both lists (i.e., when working on the 3rd point layer, it would be saving the euclidean distance layer as the 3rd output name, not the forth or another number).

2.3.2 Normalization in Raster Calculator

Subsequently, values were normalised following a distribution centralization procedure in Raster Calculator. As raster calculator arguments are immutable in python for looping (ESRI, n.d.), this procedure was manually conducted for the six raster euclidean layers by accessing the feature from the catalogue pane. An example in the script of how it would work is provided below. The logical maths to rescale values from 0-100 is:

(layer - layer minimum value)/(layer maximum value - layer minimum value)*100

The layers were outputed and saved as the file names they would have been in outputNormalize would the procedure have worked automatically. In an example of a write to text, the average values for each layer pre-normalized were written to a text file. The following text is included here.

2.3.3. Weighted Overlay

Subsequently, a weighted overlay was used to join the maps into one single image. Using the weighted sum function, all layers (as a list) were added together and renamed into a new raster file called ‘weightedOverlay’. No specifications were made for weights, and as optional parameters, these were then assumed to be equal across all layers (i.e., 100/6 = 16.667%). The resulting map layer shows a heat map of land suitability where more preferable locations for tourist locations should be built.

2.4 Python Script

One loop was used in the python script. This for loop used indexed lists to iterate euclidean distance over the length of the amount of layers in the imput array. It was written this way by range(len() so that any array can be used in future, and it will iterate for all possible arrays of imput layers, not just the six in this instance. The loop is on lines 31-33, with the the indexed lists shown in lines 24-26.

3. Results

Following the alterations to the model, Euclidean distances and normalized values were computed for all points. The final result is a land suitability map which ranks suitable areas based off of their proximity to key landmarks or other tourism features. The first example of euclidean distances is presented in Figure 4 and the weighted overlay is presented in Figure 5. Note, how following the raster calculation procedure for normalization that the weighted overlay is expressed in values between 0-100 (which it would not be if values were not normalized or reclassified.

4. Discussion

The advent and use of Geographic information systems for decision making shows much promise. Especially in instances where multiple criteria must be considered, and their consideration non-numerically would lead to erroneous conclusions, GIS-MCDA shows great utility. In addition, the use of python to modulate and alter base functions in ArcMap makes the application of different decision support models unique and highly applicable to individual scenarios. Specific to tourism site selection, or site selection in general, smarter and more data-informed decisions directly contribute to smart city development, and in turn, are valuable tools to be considered as part of sustainable development.

4.1 Next steps:

• feasible parcels for development should be used to constrain feasible alternatives (e.g., developable sites) for future tourism locations. Parcels with the highest average rating from the weighted overlay should be considered for development.

• A hot-spot-cold-spot (e.g., Geo-Getis-ordination) should be run to compute feasible sites for future development of tourism activities, or to identify areas that are deprived of tourism and landmarks and of possible redevelopment potential.

5. Limitations and Challenges

In this tutorial I learned that Euclidean distance requires multiple points to compute distance from and cannot be one central point, as distances are calculated as relative between points. I also learned how to index lists and call upon them in functions so that the right variables will be called upon, in proper order that the function iterates on them. In future, I would run an aggregation test for a hot-spot cold spot analysis to source actual development locations as opposed to solely just computing a land suitability layer. Also, the way i wrote the file with these indexed lists meant that, as long as it was within the function, it would write for any amount of outputs to a text file, and would not have to be done for all, then aggregated later - this is just a simple example, but one thing it showed me is how inside and outside of functions can be used to logically develop calculations.

6. References

Behrens, T., Schmidt, K., Viscarra Rossel, R. A., Gries, P., Scholten, T., & MacMillan, R. A. (2018). Spatial modelling with Euclidean distance fields and machine learning. European journal of soil science, 69(5), 757-770.

Danielsson, P. E. (1980). Euclidean distance mapping. Computer Graphics and image processing, 14(3), 227-248.

Eldrandaly, K. A., & AL-Amari, M. A. (2014). An expert GIS-based ANP-OWA decision making framework for tourism development site selection. International Journal of Intelligent Systems and Applications, 6(7), 1.

Elwin, J. (1989). Site selection in tourism. Site selection in tourism., 403-406.

Fernandes, I. D., Ferreira, F. A., Bento, P., Jalali, M. S., & António, N. J. (2018). Assessing sustainable development in urban areas using cognitive mapping and MCDA. International Journal of Sustainable Development & World Ecology, 25(3), 216-226.

Fraser, C., Bernatchez, P., & Dugas, S. (2017). Development of a GIS coastal land-use planning tool for coastal erosion adaptation based on the exposure of buildings and infrastructure to coastal erosion, Québec, Canada. Geomatics, Natural Hazards and Risk, 8(2), 1103-1125.

Malczewski, J., & Rinner, C. (2015). Introduction to GIS-mcda. In Multicriteria decision analysis in geographic information science (pp. 23-54). Springer, Berlin, Heidelberg.

Mühlbacher, A. C., & Kaczynski, A. (2016). Making good decisions in healthcare with multi-criteria decision analysis: the use, current research and future development of MCDA. Applied health economics and health policy, 14(1), 29-40.

Sahani, N. (2019). Assessment of ecotourism potentiality in GHNPCA, Himachal Pradesh, India, using remote sensing, GIS and MCDA techniques. Asia-Pacific Journal of Regional Science, 3(2), 623-646.

Steele, K., Carmel, Y., Cross, J., & Wilcox, C. (2009). Uses and misuses of multicriteria decision analysis (MCDA) in environmental decision making. Risk Analysis: An International Journal, 29(1), 26-33.

Figure 6 (below). Copy of script, organized by different functions. Link available here

# -*- coding: utf-8 -*-

# ---------------------------------------------------------------------------

# py2try.py

# Created on: 2022-02-21 09:32:32.00000

# (generated by ArcGIS/ModelBuilder)

# Description:

# ---------------------------------------------------------------------------

# 1. Import processing module:

import arcpy

import numpy

# 2. Import Local variables from Geodatabase:

camp = "021L14_campgro_p"

historical = "021L14_hi_site_p"

cemetery = "021L14_cemeter_p"

resevoirs = "021L14_u_reser_p"

wetland = "021L14_wetland_a_Layer"

picnic = "021L14_picnic_p"

textfile = open('C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\text.txt', "w")

# 3. Specify List objects for for function, includes lists of file paths, layers for imput, and layers for outputs (per the rescaled euclidean distance layers, respectively.

paths = ['C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file1', 'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file2', 'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file3', 'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file4', 'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file5', 'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file6']

layers = [camp, historical, cemetery, resevoirs, wetland, picnic]

outputs = ['file1', 'file2', 'file3', 'file4', 'file5', 'file6']

# declare as an array, write into it in python

MeanDistanceValue = [0, 0, 0, 0, 0, 0]

# 4. Process to run through both lists, which are indexed

for x in range(len(layers)):

arcpy.gp.EucDistance_sa(layers[x], paths[x], "", "", "", "PLANAR", "", "")

# 5. Raster Calculator Automation Attempt

# "In Python, Map Algebra expressions should be created and executed with the Spatial Analyst module, which is an extension of the ArcPy Python site package. They are non computable functions within python, and should be excecuted within the function pane of ArcMap."

# This is what I would use if there was a way in ArcMap to loop or python code running rater math operations (THERE IS NOT).

inputNormalize = ["(\"%file1%\"-\"%file1%\".minimum)/(\"%file1%\".maximum - \"%file1%\".minimum)*100",

"(\"%file2%\"-\"%file2%\".minimum)/(\"%file2%\".maximum - \"%file2%\".minimum)*100",

"(\"%file3%\"-\"%file3%\".minimum)/(\"%file3%\".maximum - \"%file3%\".minimum)*100",

"(\"%file4%\"-\"%file4%\".minimum)/(\"%file4%\".maximum - \"%file4%\".minimum)*100",

"(\"%file5%\"-\"%file5%\".minimum)/(\"%file5%\".maximum - \"%file5%\".minimum)*100",

"(\"%file6%\"-\"%file6%\".minimum)/(\"%file6%\".maximum - \"%file6%\".minimum)*100"]

outputNormalize = ['C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file1_raster',

'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file2_raster',

'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file3_raster',

'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file4_raster',

'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file5_raster',

'C:\\Users\\bserre\\Documents\\ArcGIS\\Default.gdb\\file6_raster']

arcpy.gp.RasterCalculator_sa(inputNormalize[x], outputNormalize[x])

# In turn, Normalization was then conducted individually using the catalogue pane to run rasterCalculator. Would look like something below >>

# the outputNormalize list is the same file directory of the files once they were normalized following the raster math normalization procedure (layer-layerminimum/layermaximum-layerminimum).

# 6. Combining Rasters into single weighted sum

WeightedSum_sa(outputNormalize, "weightedOverlay")

# 7. Summarizing Raster Values, and writing to file.

MeanDistanceValue[x] = arcpy.GetRasterProperties_management(outputNormalize[x], "MEAN")

print(MeanDistanceValue[x])

temp = " File Name: " + outputs[x] + ", Mean Value : " + str(MeanDistanceValue[x]) + "km \n"

textfile.write(temp)

textfile.close()

##END##