Comparison of Multicriteria Analysis Techniques for Environmental Decision Making on Industrial Location

Introduction

European legislation calls for a well-planned sustainable development. As such, it has to include a social, economic as well as an environmental dimension. According to Agenda 21 (http://www.un.org/esa/dsd/agenda21/), countries should undertake efforts to build up a comprehensive national inventory of their land resources in order to establish land information systems. The overall objective is to provide information for the improvement or the restructuring of land-use decision processes including the consideration of socio-economic and environmental issues.

In the last decades conflicts caused by competing land uses have increased, particularly in urban areas. Consequently, a lot of research has been done aiming to develop methods and tools that assist complex spatial decision problems. The development of Spatial Decision Support Systems (SDSS) has turned out to be very beneficial in assisting to the solution of complex land-use problems [1, 3].

In addition, any planning process must focus on a mix of hard (objective) and soft (subjective) information. The former are derived from reported facts, quantitative estimates, and systematic opinion surveys. The soft information denotes the opinions (preferences, priorities, judgments, etc.) of the interest groups and decision makers. The idea of combining the objective and subjective elements of the planning process in a computer based system lies at the core of the concept of SDSS [1, 3].

SDSS can be defined as an interactive, computer-based system designed to support a user or a group of users in achieving a greater degree of effectiveness in decision making when solving a semi-structured spatial decision problem [3]. SDSS also refers to the combination of GIS and sophisticated decision support methodologies, e.g. in terms of multicriteria analysis techniques [3, 6], and are therefore suitable to manage sustainable development of urban areas.

Although the development of multicriteria analysis began mainly in the '70s (the first scientific meeting devoted entirely to decisionmaking was held in 1972 in South Carolina) its origins can be dated back to the eighteenth century [4]. Reflections on French policies in the action of judges and their translation into policy (social choice), led people like Condorcet to deepen in decision taken supported in several criteria [4].

In the last two decades of the twentieth century there was an increased trend of integration of Multicriteria Evaluation techniques (MCE) and Geographic Information Systems (GIS), trying to solve some of the analytical shortcomings of GIS “For example see [4, 7, 15]“. Wallenius et al. [16], made a study of the evolution in the use of MCE techniques from 1992 to 2006, showing that the use of multiattribute techniques has increased 4.2 times during this period. In recent years, there has also been a great effort in the integration of MCE and GIS techniques on the Internet “For example see [17, 20]”.

Since we consider land-use decision making in general as an intrinsic multicriteria decision problem, in our opinion these are valid methodologies to support the land-use decision process by means of a land-use suitability analysis.

Land-use suitability analysis aims to identify the most appropriate spatial pattern for future land uses according to specified requirements or preferences [3, 21, 22]. GIS-based land-use suitability analyses have been applied in a wide variety of situations, including ecological and geological approaches, suitability for agricultural activities, environmental impact assessment, site selection for facilities, and regional planning [3, 6, 11, 17, 21,23, 28].

Different attempts to classify Multicriteria Decision Making (MCDM) methods by diverse authors exist in the literature [4, 6, 7, 11, 26, 29]. The majority of them agree that additive decision rules are the best known and most widely used Multiattribute Decision Making (MADM) methods in GIS based decision making. Some of the techniques more commonly described in literature are: Simple Additive Weighting (SAW), Ordered Weighting Averaging (OWA) technique, the Analytical Hierarchy Process (AHP), ideal point methods (e.g. TOPSIS), concordance methods or outranking techniques (e.g. PROMETHEE, Electre).

Nevertheless, the integration of these techniques continues to pose certain problems or difficulties at the time of developing specific applications. Among the most notable drawbacks are [4]:

• The impracticality of applying pairwise comparison techniques as PROMETHEE with long series of data due to limitations posed by existing informatics systems.

• The difficulty on the implementation of some MCE methods, thereby leading to a difficult analysis of the results, as well as an ignorance of the internal procedure of the methods by non-specialist users.

• The need to generate data processing software attached to the GIS, based on algorithms that describe MCE methods, which naturally implies that many users of these systems cannot access these methods.

In this chapter, we compare the results obtained by the application of two distinctive land-use suitability analyses to the location of industrial sites, applying two different multicriteria analysis techniques. The multicriteria analysis employed has been performed in a raster environment and been used for two objectives. During the site search analysis each pixel was considered a potential location alternative. This analysis used a SAW method which signifies a weighted summation. It can thus easily be performed in GIS [24, 25, 30, 31]. A site selection analysis then used the PROMETHEE-2 methodology [32] and a set of predefined alternatives [30, 33]. All of the techniques used in the project were coded and integrated within ArcGIS by Marinoni [5, 34].

A problem in the application of multicriteria analysis is the definition of weights for a given set of criteria. A variety of approaches does exist, see for example [26], and the probably best known weight evaluation method is the AHP [35], which we have used in our case as well.

Another problem is the specification of the criteria performance scores which are often subjective in their determination. Data which have been measured directly will certainly be regarded as more reliable than data which have been estimated, interpolated, taken from a map or simply interpreted. Thus, the method of criteria data collection plays a central role [5]. A stochastic approach which takes account of the uncertainty of input values and which is presented at a last step in this chapter could be a way out of this dilemma.

1. Background and Methodology

1.1. Study area and project background

Zaragoza city and its surroundings are located in the Ebro corridor, a highly dynamic economic area within the Iberian Peninsula. The climate in this area is semi-arid with mean annual precipitation of about 350 mm and a mean annual temperature of about 15° C.This city is crossed by the cited Ebro river and two of its main tributaries, the Gállego and Huerva rivers (Figure 1). Geologically, Quaternary alluvial terraces of the Ebro river were deposited above Tertiary gypsum formations, forming a covered karst area with intense karstification processes. The Quaternary materials are an important source of sand and gravel which are needed for civil engineering purposes. In addition, it hosts important groundwater reservoirs, used for domestic, industrial and agricultural purposes.

The availability of these resources has been one of the reasons of the fast development of the city in the last decades. But this fast development has also led to negative interactions with the environment and man-made infrastructure. Intense irrigation triggered land subsidence which in turn caused costly damage and/or destruction of infrastructure such as roads, buildings, gas and water supply networks [36]. Many infrastructures that have been built occupy areas where soils of high fertility had naturally developed, making these areas inaccessible to agriculture. Also, many ecologically important areas have been harmed and an increased contamination of the aquifer has been observed [37].

Figure 1.

Location and geomorphology of study area.

Based on the above, the area surrounding Zaragoza, which represents a rapidly growing urban area, merits closer investigation in terms of geoscientific factors. Thus, a research project was initiated to develop a methodological workflow which will facilitate the sustainable development in the surroundings of a growing city. Our main objective was to perform a land-use suitability analysis to identify the most appropriate future land-use patterns. Therefore a variety of tasks needed to be performed such as:

• Characterization of the study area and collection, analysis and processing of the available information for its introduction into a GIS environment.

• Geo-hazards and geo-resources detection, description and modelling with the help of GIS and 3D techniques.

• Land-use suitability analysis by means of SDSS.

Here, we report on the land-use suitability analysis to find most suitable locations for industrial facilities. As mentioned above, we compare the results obtained by the application of two distinctive multicriteria analysis techniques for environmental decision making on industrial location. For more details on the general project workflow and geo-resources and geo-hazards modellingsee [24, 25, 30, 31, 33, 37, 39].

1.2. Methodology

It is important to differentiate between the site selection problem and the site search problem. The aim of site selection analysis is to identify the best location for a particular activity from a given set of potential (feasible) sites. Where there is no predetermined set of candidate sites, the problem is referred to as site search analysis [3].

In terms of the MCE methods applied, the main advantage of the SAW approach can be considered its low degree of complexity as which made it attractive to be used for the site search analysis in this project. It is precisely this simplicity that makes weighted summation actually quite widely applied in real-world settings [8, 40, 42].

The site selection analysis has been performed by the implementation of PROMETHEE-2 which belongs to the ‘family’ of outranking techniques. Since the mentioned techniques require pairwise or global comparisons among alternatives, these methods become impractical for applications where the number of alternatives ranges in the tens or hundreds of thousands (Pereira and Duckstein, 1993). For a more detailed description of both methodologies see [24, 25, 30, 33, 35, 36].

Figure 2.

Workflow of the land-use suitability analysis.

In order to perform both site search and site selection, several steps needed to be covered. These included (Figure 2):

• Definition of alternatives (decision options): feasible location areas.

• Definition of constraints: areas with land-use restrictions.

• Definition of important factors in the decision process: identification of criteria.

• Determination of criteria weights

The criteria weights were determined with the AHP. This technique represents another MCE method and involves pairwise comparison of criteria where preferences between criteria are expressed on a numerical scale usually ranging from 1 (equal importance) to 9 (strongly more important). This preference information is used to compute the weights by means of an eigenvalue computation where the normalized eigenvector of the maximum eigenvalue characterizes the vector of weights. Empirical applications suggest that this pairwise comparison method is one of the most effective techniques for spatial decisionmaking approaches based on GIS [15, 43]. There exist many well-documented examples of application of this method with success [44, 46].

It is well known that the input data to the GIS multicriteria evaluation procedures usually present the property of inaccuracy, imprecision, and ambiguity. In spite of this knowledge, the methods typically assume that the input data are precise and accurate. Some efforts have been made to deal with this problem by combining the GIS multicriteria procedures with sensitivity analysis [47] and error propagation analysis [48]. Another approach is to use methods based upon fuzzy logic [3].

In many situations it is hard to choose the input values for multicriteria analysis procedures, since the criteria values for the different alternatives usually do not have a single realization, but can obtain a range of possible values [5]. Performing a multicriteria analysis with the mean values produces some kind of mean result, but the uncertainty in either the input values or the result cannot be quantified. A solution to this dead-end is a stochastic approach, which utilizes probability distributions for the input parameters instead of single values. A stochastic multicriteria analysis implies that the analysis is performed multiple times with varying input values for the criteria involved. These criteria input values (or performance scores) are drawn from probability distributions that are inferred from empirical criteria populations (e.g pixels on a map, expert knowledge). Such an approach uses the whole range of possible criteria value outcomes and extreme events are according to their low outcome probabilities realistically represented as rare events. In a last step we explored the influence of criteria weightsby conducting a sensitivity analysis.

1.2.1. Site search analysis

Within the site search analysis, every pixel was considered a decision alternative. Constraints depict the areas where industry is and will not be allowed. These restrictions are generally characterized by the existence of other land uses (e.g. urbanareas), the protection of natural areas and land management planning. These restrictions are (Figure 3):

• Natura 2000 network areas: natural reserve of the oxbows in La Cartuja (map provided by the Aragon Government).

• Urbanized areas: obtained from the topographic map scale 1:25,000 from the National Geographical Institute (IGN, Instituto Geográfico Nacional), imported to ArcGIS and updated.

• Infrastructures (roads, rail roads, canals) and their area of protection: also extracted from the topographic maps. The area of protection of roads and train rails was delineated as defined by the Spanish Roads Law and according to the Spanish Railway Sector Law, respectively.

• Other restrictive planning: Zaragoza Land Management Planning (PGOUZ, Plan General de Ordenación Urbana de Zaragoza), mapping provided by Zaragoza Council, and natural resources planning of the thickets and oxbows of the Ebro river, provided by the Aragon Government.

• Cattle tracks: tracks traditionally used by the seasonal migration of livestock which are protected by law, provided by the Aragon Government.

• Industrial areas where no space is left for new industries. Provided by the Aragon Institute of Public Works (IAF,Instituto Aragonés de Fomento) from the Aragon Government.

Figure 3.

Industrial restrictions.

A variety of social, economic and environmental factors were taken into consideration. Figure 4 shows the mapping of all the variables that were considered relevant for industrial development. Areas considered less suitable are kept in red while a higher suitability is shown in green.These variables are:

• Important areas from the environment point of view: natural areas included in the Natura network 2000 as SPAs (Special Protection Areas for birds), and the SACs (Special Conservation Areas), habitats, points of geological interest and other areas which mapping has been provided mainly by the Aragon Government.

• Doline (sinkhole) susceptibility: model developed within the project using a quantitative method, a logistic regression technique [38].

• Groundwater protection: a model developed also within the project, performed with Gocad [37] and applying a methodology by the German Geological Survey [49].

• Flooding hazard: a flooding hazard mapping developed along the Ebro river [50] was digitised and introduced in the land-use suitability analysis. This model shows the different periods of return of flood events.

• Agricultural capability of the soils: mapping developed within the project [39] applying the Cervatana Model [51].

• Slope of the terrain: developed from the DEM (resolution 20x20 m) from the Ministry of Agriculture (SIG oleícola).

• Geotechnical characteristics of the subsoil: different geomorphological units with better or worse geotechnical characteristics, described according to the PGOUZ. This classification has been applied to the geomorphological units derived from the geological map, scale 1:50,000, from National Geological Survey (IGME, Instituto Geológico y Minero de España).

Figure 4.

Variable mapping. 1) natural protected areas, 2) doline susceptibility, 3) groundwater protection, 4) flooding hazard, 5) agricultural capability of the soils, 6) slope percentage, 7) geotechnical characteristics.

Many multicriteria methods, as the SAW methodology, require criteria standardization to bring all of them to a common scale. The classification ensures that the weights properly reflect the importance of a criterion.The standardization method used here may be classified as a subjective scales approach [26] since the variables are classified in subjective ranges. These ranges can be selected following standards, legal requirements, or the classes already determined in the geo-resources and geo-hazards models used as criteria in the decision process. Six categories were selected considering the adaptation of these classes to the variables to be introduced. For more details on the standardization approach see [25, 30, 31].

Weights for criteria are assigned with the help of the AHP. An AHP extension was specifically developed for the ArcGIS environment at the Institute of Applied Geosciences of the Technische Universität Darmstadt [34]. This tool can be downloaded from the ESRI web page (http://arcscripts.esri.com/). For more details on the AHP performance see [25, 30, 31]. In a last step all classified raster files (criteria) are multiplied by its corresponding weight and summed up.

1.2.2. Site selection analyst

The main objective of a site selection analysis is the ranking of feasible alternatives. Generally, outranking methods, such as PROMETHEE-2, require pairwise or global comparisons among alternatives. Here location alternatives are represented by industrial areas, as defined in the Aragon Institute of Public Works (IAF) database, which signify spaces for the establishment of new industries. Geometrically, these alternatives represent a polygon each. A total of twenty seven industrial areas were evaluated for the site selection analysis.

As alternatives are directly compared along their criteria values, the application of outranking methods does not require a transformation or standardization of criteria values. The restrictions (constraints) and criteria are the same used for the site search analysis. Alternatives located completely in restrictions areas were eliminated from the analysis. However, there exist some industrial polygons, representing one alternative, located partially in restricted areas, as these polygons are partially occupied or crossed by a road or a cattle track. It has implied the inclusion of the constraints as an additional criterion in the decision process. The criterion representative of use restrictions was then reclassified into two different values; zero in the area where industry is forbidden or not possible due to the presence of other uses, and one in areas where this use is permitted or feasible.

It is important to define whether a higher value of a particular criterion leads to an improvement or to a decrease in land-use suitability. In the case of industrial development, an increase in the value of all criteria, with the exception of groundwater protection and geotechnical characteristics, implies a suitability decrease. For example, a higher groundwater protection value implies an increase in suitability to industrial use location while an increase in doline development susceptibility implies a decrease in industrial use location suitability.

Geometrically, every alternative is a polygon so that within each polygon a variety of criteria values (pixels in the criteria layers) are to be found. The question then arises which of the multiple criteria realization to use for the multicriteria analysis evaluation. Therefore, a multicriteria GIS extension was developed to draw site specific values (minimum, maximum, mean etc.) for raster cell populations that lie within the polygonal outline of a location alternative. For our analysis the mean value was used for all criteria since, in our opinion, this value better symbolizes all alternative values. Minimum and maximum values are usually rare events with a low probability of occurrence.

PROMETHEE-2 methodology uses preference function, which is a function of the difference between two alternatives for any criterion [32]. Six types of functions based on the notions of criteria, are proposed. For more details on preference functions see [5, 30, 32].

We exclusively used the “usual criterion” preference function that is based on the simple difference of values between alternatives as this function helps to discriminate best between available alternatives which we wanted to achieve.

The pair comparison of alternatives produces a preference matrix for each criterion (Figure 5). Having calculated the preference matrices along each criterion, a first aggregation is performed by multiplying each preference value by a weighting factor w (expressing the weight or importance of a criterion), and building the sum of these products [5]. This results in a preference index, Π (see Figure 5). The AHP has also been integrated in this tool and used for criteria weighting.

Figure 5.

Schematic calculation of the preference index Π. Source [5].

The final ranking of alternatives is performed by calculating the net flow Φ (a1) for every alternative, a, which is a subtraction between the leaving flow and the entering flow. The higher the net flow is, the higher is the preference of an alternative over the others (Table 1).The leaving flow Φ+ (a1) represents a measure of the outranking character of a1 (how a1 is outranking all the other alternatives). Symmetrically, the entering flow Φ- (a1) is giving the outranked character of a1 (how a1 is dominated by all the other actions).

1.2.3. Stochastic PROMETHEE-2

The stochastic PROMETHEE-2 approach requires the assignment of theoretical distribution types to every criterion of the available alternatives. Distribution models were inferred based upon the criteria value populations (pixel values) within each location alternative (polygon) along all criteria. The software used to fit distribution types and to perform distribution fitting test was @Risk [52]. In a next step the distribution models were used within a Monte Carlo Simulation (MCS). The number of iterations nwas set to 5000.

Starting a MCS with niterations for the specified distributions produces nrealizations for every cell of the input matrix [5]. Figure 6 shows the principle of one iteration cycle.Values are randomly drawn between 0 and 1 and input values (criteria performance scores) are determined using the inferred theoretical model distribution. With nbeing 5000, the multicriteria analysis is repeated 5000 times. The results may then be used to establish a rank distribution for a specific alternative or a distribution of alternatives for a specific rank (see Figure 7 for a four hypothetical scenarios demonstration).

Π	a1	a2	a3	Φ+(ax)	Φ(ax)	Rank
a1	-	0.25	0.75	1.0	0	2
a2	0.75	-	0.75	1.5	1	1
a3	0.25	0.25	-	0.5	-1	3
Φ-(ax)	1	0.5	1.5

Table 1.

Example of possible preference indices, leaving, entering and net flow calculations and final ranking.Slightly modified after [5].

However, the alternative possessing the highestnumber of first ranks may not necessarily be the best [5]. Therefore, it was suggested calculating a dimensionless mean stochastic rank MSR for every alternative.

$$MSRAj,m=\frac{1}{n}{\displaystyle \sum \begin{array}{l}n\\ i=1\end{array}}\left(Ri*i\right)\forall j=\mathrm{1,}\dots ,n$$E1

where:

m: number of iterations

A_j: j^th alternative

n: number of available alternatives

R_i: rank count for the i^th rank

In order to compare mean stochastic ranks of simulations with different iteration counts, the MSR value must be standardized which leads to the stochastic rank index SI [5]:

$$SIAj,m=\frac{MSRAj,m-MSRmin,m}{MSRmax,m-MSRmin,m}$$E2

where:

m: number of iterations

SI_Aj: stochastic rank index for the j^th alternative

MSR_Aj: MSR for the j^th alternative

MSR_min: the lowest possible MSR value

MSR_max: the largest possible MSR value

The more the SI value approaches 0, the better the alternative.

Figure 6.

PROMETHEE input value determination for one iteration cycle. Source: [5].

Figure 7.

Left: example distribution of 4 scenarios (s1,..., s4) for rank 1. Right: rank distribution for scenario 1. Source: [5].

2. Results and validation

2.1. Site search analysis

The criteria preference values have been assigned by the authors after discussions with experts from different stakeholder groups from the Zaragoza City Council and the Ebro River Authority (CHE, Confederación Hidrográfica del Ebro). The highest preference values (and therefore the highest weights) were given to the groundwater protection and environmentally high value areas (Table 2). Hazard criteria were considered less important as some of the encountered geological hazards (land subsidence and sinkhole development) can be mitigated or avoided by applying more suitable (but more costly) construction techniques.

The validation of a model consists in checking whether the structure of the model is suitable for the purpose and if it achieves an acceptable level of accuracy in predictions. In the case of explanatory or predictive models, validation is usually carried out by checking the degree of agreement between the data produced by the model and data from the real world [4]. In the case of our project in order to validate the model has been verified that the result follows the preferences in the assignation of the weights to the criteria.

Preference matrix	A	B	C	D	E	F	G	Weight
A	1.00	2.00	3.00	1.00	5.00	8.00	6.00	0.2288
B	0.50	1.00	2.00	0.50	3.00	7.00	4.00	0.1736
C	0.33	0.50	1.00	0.33	2.00	6.00	3.00	0.1131
D	1.00	2.00	3.00	1.00	5.00	8.00	6.00	0.2288
E	0.20	0.33	0.50	0.20	1.00	4.00	2.00	0.0678
F	0.13	0.14	0.17	0.13	0.25	1.00	0.50	0.0251
G	0.17	0.25	0.33	0.17	0.50	2.00	1.00	0.0427

Table 2.

Pairwise comparison matrix, criteria weights for site searchanalysis.A) Groundwater protection, B) Doline susceptibility, C)Flooding hazard, D)Location of natural areas, E) Agricultural capability of soils, F) Slope percentage, G) Geotechnical characteristics.

Figure 8 shows the final results of the land-use suitability analysis for new industrial development. The lefthand side of figure 8 shows the suitability map under sustainability. The grey sections indicate the areas where industrial locationis not possible due to the constraints. Although the suitability analysis sometimes presents good values, the constraints imply that these areas cannot be exploited due to any restriction.The most suitable locations for industrial development are on the pediments or glacis (Figure 1) and Tertiary materials outside environmental protected areas where the groundwater vulnerability and flood risk is lower. The least suitable locations are the floodplains with high groundwater vulnerability and flood risk, environmentally protected areas around the river bed and other areas in the higher terraces which present more susceptibility to doline development.

To test the robustness of the results, a sensitivity analysis of the model has been performed where higher weights were given to economic aspects. The highest weights were assigned to doline susceptibility and flooding hazard, which might cause the destruction of future industrial sites (Table 3). Slope and geotechnical characteristics of the soils were also assigned high values as a more technically difficult terrain will increase the construction budget. Figure 8 shows the results of this last approach where the best locations for industry were identified to be the pediments and slopes in Tertiary sediments. The least favorable locations are on the flood plain and low river terraces, where sinkhole susceptibility shows higher values. In fact, in order to measure the correlation between both results the Pearson coefficient of correlation between both raster images has been calculated giving a value of 0.874, significant at a 0.01 level, implying a high agreement between both results.

Figure 8.

Result of a) site search analysis under sustainability and b) sensitivity analysis for industrial development site search analysis (objective economic development).

2.2. Site selection analysis

As a consequence of the introduction of a new criterion depicting restriction of use or constraints as explained in section 2.2.2., the criteria weights used for the site search analysis are not valid implying a new calculation of them using the AHP. Table 4 shows the preference matrix and criteria weights of the site selection analysis for industrial development. The criteria preference values are the same as for the site search analysis (Table 2) under the sustainability scenario, however the constraints obtained the highest preference values and as a consequence the highest weight, in order to avoid the outranking of alternatives located partially in forbidden areas.

The preference indices and the leaving and entering flow generated after the application of PROMETHEE-2 methodology are presented in Table 5. Figure 9 shows the location of the alternatives of the site selection analysis. The best alternatives are generally located south of Zaragoza city, outside the alluvial sector (i.e. alternatives 25, 26 and 27). In contrast, the worst locations are the alluvial areas in the surroundings of El Burgo de Ebro (i.e. alternatives 4 and 18), the industrial areas in the north of Zaragoza city (i.e. alternatives 16 and 17), and the Logroño Road Corridor, upstream of Zaragoza (i.e. alternative 8).

Preference matrix	A	B	C	D	E	F	G	Weight
A	1.00	0.25	0.33	2.00	5.00	0.33	0.33	0.0735
B	4.00	1.00	4.00	5.00	8.00	2.00	3.00	0.3492
C	3.00	0.25	1.00	4.00	7.00	1.00	2.00	0.1774
D	0.50	0.20	0.25	1.00	4.00	0.25	0.25	0.0505
E	0.20	0.13	0.14	0.25	1.00	0.14	0.14	0.0224
F	3.00	0.50	1.00	4.00	7.00	1.00	2.00	0.1892
G	3.00	0.33	0.50	4.00	7.00	0.50	1.00	0.1378

Table 3.

Pairwise comparison matrix, criteria weights for site search analysis under economic aspects.A) Groundwater protection, B) Doline susceptibility, C)Flooding hazard, D)Location of natural areas, E) Agricultural capability of soils, F) Slope percentage, G) Geotechnical characteristics.

Preference matrix	A	B	C	D	E	F	G	H	Weight
A	1.00	2.00	3.00	1.00	5.00	8.00	6.00	0.50	0.197
B	0.50	1.00	2.00	0.50	3.00	7.00	4.00	0.33	0.121
C	0.33	0.50	1.00	0.33	2.00	6.00	3.00	0.25	0.087
D	1.00	2.00	3.00	1.00	5.00	8.00	6.00	0.50	0.197
E	0.20	0.33	0.50	0.20	1.00	4.00	2.00	0.17	0.048
F	0.13	0.14	0.17	0.13	0.25	1.00	0.50	0.11	0.020
G	0.17	0.25	0.20	0.17	0.50	2.00	1.00	0.14	0.030
H	2.00	3.00	4.00	2.00	6.00	9.00	7.00	1.00	0.300

Table 4.

Pairwise comparison matrix, criteria weights for site selection analysis.A) Groundwater protection, B) Doline susceptibility, C)Flooding hazard, D)Location of natural areas, E) Agricultural capability of soils, F) Slope percentage, G) Geotechnical characteristics, H) Constraints.

2.3. Stochastic PROMETHEE-2

In the stochastic approach, distribution types have to be assigned to every alternative and criterion and a MCS is performed over a MCE method (here PROMETHEE-2) meaning that the multicriteria analysis is performed a specified number of times (here: 5000; hence stochastic PROMETHEE-2). It should be noted that, due to local/regional variability, the local distributions of a criterion are highly likely to be different for each location alternative. Thus, although it seems reasonable, at first sight, to determine one distribution type for one criterion, if location dependent statistical analyses indicate varying distribution types, then varying types should be assigned to one criterion.

Alt	Φ-	Φ+	Φ	Rank PROMETHEE	SI	Rank stochastic	SAW; Mean value
1	13.69	12.22	-1.48	14	0.41	10	3.95
2	13.41	8.74	-4.67	20	0.75	22	3.09
3	13.65	8.62	-5.03	22	0.65	19	3.26
4	15.76	8.50	-7.25	25	0.50	13	3.95
5	12.77	9.28	-3.49	16	0.70	21	3.25
6	7.48	12.93	5.45	8	0.27	7	4.60
7	6.42	17.08	10.66	4	0.06	1	5.38
8	13.31	6.91	-6.40	24	0.69	20	3.55
9	11.97	8.25	-3.72	17	0.56	14	3.75
10	12.30	8.11	-4.18	18	0.63	17	3.32
11	9.49	10.92	1.44	10	0.42	11	3.69
12	13.30	8.96	-4.34	19	0.63	18	3.79
13	8.99	11.42	2.44	9	0.56	15	3.57
14	8.61	15.29	6.69	7	0.19	5	4.66
15	10.74	9.97	-0.77	12	0.60	16	4.76
16	15.17	4.93	-10.25	26	0.95	27	3.61
17	13.40	7.00	-6.40	23	0.88	26	3.75
18	19.37	6.63	-12.75	27	0.86	25	2.55
19	11.97	8.74	-3.24	15	0.45	12	4.32
20	12.74	7.78	-4.96	21	0.76	24	3.9
21	10.83	9.58	-1.25	13	0.76	23	3.83
22	5.25	15.16	9.91	5	0.33	8	4.68
23	4.90	15.82	10.92	3	0.21	6	5.68
24	10.13	9.98	-0.15	11	0.39	9	5.61
25	5.78	17.63	11.85	1	0.11	4	4.88
26	6.18	17.32	11.14	2	0.10	3	4.80
27	7.22	17.04	9.83	6	0.07	2	5.32

Table 5.

Leaving floe, entering flow, net flow and rank for site selection analysis, stochastic rank index and final rank for stochastic approach, and mean value in SAW methodology for every alternative of location.

Table 6 show the distribution types assigned to every alternative and criterion for the suitability analyses. Distribution fitting tests were performed to confirm/reject a modeled distribution type. The software used was @Risk [52]. If physical properties can only have non-negative values distribution types can (and should) be selected such that this feature is reflected, for example by choosing an exponential distribution.

Figure 9.

Location of alternatives for site selection analysis.

The more commonly used distributions for continuous variables (i.e. slope percentage) are normal and lognormal, but also logistic and exponential distributions are present in some variable (i.e. groundwater protection and/or doline susceptibility).

A binomial distribution was selected for categorical variables having two possible outcomes.If there are more than two categories (possible outcomes) the use of a categorical distribution can be problematic, implying the inclusion in the decision process of categories not present in the alternative. For example, if one alternative presented values 1 and 4 in agricultural capability criterion, the distribution selected by the fitting test would have given values 2 and 3 to this alternative, which are not present in the real world. Thus, instead of assigning a distribution, the percentage of cases (p value in Table 6) in every category was calculated and used as the probability of occurrence of every category. This was also the case for some continuous variables, which presented few different values, thus complicating the distribution selection (i.e. alternative 3 in doline susceptibility criterion). In these cases, the percentage or probability of occurrence of every value was introduced in the analysis.

In the case of the criterion “susceptibility to doline development”, difficulties were experienced as some alternatives showed continuous values close to value 0 (see Figure 10). Since it was not possible to apply the percentage of values in these cases, a decision was made to apply an exponential distribution in order to avoid the introduction of negative values in the suitability analysis, even though the adopted solution was not absolutely satisfactory. Finally, some alternatives presented the same value for the whole alternative (unique value in the tables). Some representative examples of the selected distribution types can be seen in Figures 10. The bars symbolize the original (empirical) values retrieved from the pixels within each location alternative; the solid line the fitted theoretical distribution model.

Al.

A

B

C

D

E

F

G

H

Al.

A

B

C

D

E

F

G

H

1

p

ln

p

p

p

ln

u

b

15

p

ln

u

u

p

ln

p

b

2

l

n

p

p

p

ln

u

b

16

e

p

u

u

u

ln

p

u

3

e

p

u

u

u

n

u

b

17

l

n

u

u

u

ln

p

b

4

l

n

u

p

p

ln

p

b

18

l

ln

u

p

p

e

p

b

5

l

ln

p

u

u

ln

u

b

19

l

e

u

u

p

ln

p

b

6

l

ln

u

u

p

ln

u

b

20

l

ln

u

u

u

n

p

b

7

u

u

u

p

iu

ln

p

b

21

p

e

u

u

p

n

u

b

8

l

e

u

u

u

ln

u

b

22

l

n

u

u

p

n

u

b

9

p

e

u

u

u

ln

u

b

23

p

p

u

u

p

ln

p

b

10

l

e

u

u

p

ln

u

b

24

p

p

u

u

u

n

u

u

11

e

u

u

u

p

ln

u

b

25

u

u

u

p

p

n

u

u

12

l

p

u

p

u

l

u

b

26

u

u

u

p

p

ln

p

b

13

p

e

u

u

p

n

u

b

27

p

e

u

p

p

ln

p

b

14

p

u

u

p

p

ln

p

b

Table 6.

Distribution types for every alternative and criterion for industrial settlements suitability analysis. Al.) Alternative, A) Groundwater protection, B) Doline susceptibility, C) Flooding hazard, D) Location of natural areas, E) Agricultural capability of soils, F) Slope percentage, G) Geotechnical characteristics, H) Constraints. p) percentage, ln) lognormal, u) unique value, b) binomial, l) logistic, n) normal, e)exponential, iu) Intuniform

The results of the site selection suitability analysis based on stochastic PROMETHEE-2 can be seen in Table 5. In general, there are few differences in the SI values and total flows between the first rankings: alternatives 7, 25, 23, 26 and 27 (Figure 9). In addition, all these alternatives are located in the areas with higher suitability values in the site search analysis (SAW mean value in Table 5). Nevertheless, alternative 24 presents a high mean value in the SAW methodology but is ranked 11 and 9 in the PROMETHEE-2 and the stochastic approach. This is due to the fact that the major part of this polygon is located in a restricted area. In fact, the worst rankings, alternative 16, 17 and 18, are located partially inside restricted areas. However, in the case of alternative 24 the weight assigned to the constraint factor in PROMETHEE-2 and the stochastic approach it was not enough to rank this alternative in the last positions. Besides, the first rank changes from alternative number 25 to alternative number 7, in the PROMETHEE-2 and the stochastic approach, respectively. This is the consequence of assigning a unique mean value in the PROMETHEE-2 approach to the alternatives. Figure 11 saws the values of the SAW methodology for both alternatives. It can be observed how, although both alternatives present similar mean value, alternative 7 present homogeneous high values, implying more percentage of high values in its distribution, while alternative 25 present a variety of suitability values, implying a less percentage of high suitability values. The stochastic approach overcomes this handicap by simulating values along the whole range of values inside the distribution assigned to the alternatives.

Figure 10.

Examples of distributions assigned to alternatives.

3. Discussion and conclusions

The industrial use suitability map developed with the SAW and AHP methods integrated in a GIS for the surroundings of Zaragoza, is a substantial aid in the land-use management of this city. Besides, an additional benefit is achieved by integrating geoscientific aspects in the land-use decision process, as demanded by Agenda 21.

Figure 11.

SAW values for alternatives 7 and 25.

A fundamental problem of decision theory is how to derive weights of criteria. One disadvantage of the AHP method is the inherent subjectivity of assigning preference values between criteria. The weights derived from these preference values have usually a profound effect on the results of the suitability analysis. However, in our particular case, in the industrial suitability analysis, there were no strong differences between the results of the site search analysis performed under the concept of sustainable development or the site search analysis performed under the concept of economic development, although different weights were assigned to the criteria in both approaches.

If differences are greater, a possible solution is to establish a set of suitability maps and to combine these to select the most suitable areas.

After some talks with different managers in the administration and following the approach under sustainability aspects, our results suggest that the best location for new industries is on the pediments and Tertiary sediments outside the natural protected areas, where the groundwater vulnerability and flood risk is lower, although the geotechnical characteristics of the terrain are less favorable, according to the PGOUZ. The least favorable location embodies the floodplain with high groundwater vulnerability values and the natural protected areas around the river bed, and other areas in the higher terraces which are more susceptible to doline development.

An advantage of outranking methods as PROMETHEE-2 is the fact that criteria do not need standardization or transformation processes which reduces subjectivity. However, in a spatial multicriteria analysis decisions still need to be made, as for example what characteristic value (from the population of pixels within a location alternative polygon) to use for a subsequent multicriteria analysis (e.g.maximum, minimum, mean, etc.). If using PROMETHEE-2 more decisions needs to be made in regards to the selection of the preference function as well as which set of criteria weights to use.

It is important to notice the similarity of the results after applying the site search analysis and the site selection analysis. In general, the highest rank positions are present in alternatives located in areas where the site search analysis also presented the highest suitability values. Some differences can be observed in alternatives located in areas with restrictions, as in the site selection analysis constrains are included as criteria. This is the case of alternative 24 which presented a high mean value in the SAW methodology but present a rank 11 and 9 in the PROMETHEE-2 and the stochastic approach. In this case, the weight assigned to the constraint factor in PROMETHEE-2 and the stochastic approach it was not enough to rank this alternative in the last positions. Thus, a higher weight should be given to the constraint factor in the site selection approaches.

Performing a PROMETHEE-2 with the mean values produces a mean result, but the uncertainty in either the input values or the result cannot be quantified. The stochastic approach helps approaching this problem by using probability distributions for the input parameters, instead of single values. For spatial multicriteria analysis in a variable data environment it is our recommendation to use stochastic approaches although, in this case, the process was not absolutely integrated in the GIS, and as a consequence it is very time consuming.

Acknowledgments

This research was funded by the Deutsche Forschungsgemeinschaft (DFG, Ho 804/7-1+2). We also greatly acknowledge support of the Ebro River Authority (CHE), the Aragón Region Authority and the Zaragoza Council.