During precipitation events, water which is not absorbed into the ground flows over the ground surface eventually contributing to surface water features such as rivers and lakes. This phenomenon is called run-off and it is a major contributor of non-point source pollution, which can contaminate surface water features. Run-off can also cause erosion damage to property and infrastructure and in extreme cases, risk to human life. From May 1996 to December 2017, the National Climatic Data Center (NCDC) reported a total of 1292 flood/flash flood events caused by stormwater run-off with 32 deaths, 15 injuries, over $115 million in property damage, and over $5 million in crop damage  . Urban runoff is a major contributor to elevated levels of suspended sediments and of harmful pollutants in streams, lakes, rivers, and oceans. Common runoff pollutants include nutrients, organic materials, pathogens, hydrocarbons, pesticides, metals, chlorides, grease, trash, and toxic substances  . These pollutants pose a serious threat to surface water quality and risk contamination of surrounding areas. Additionally, arid and semi-arid watersheds are especially vulnerable to run-off, due to infrequent but high-intensity episodic rainfall events which precipitate high volumes of water in a short period of time.
The management and control of stormwater run-off are of principle concern to community planners and local governments, as a matter of hazard mitigation as well as prevention of contamination to scarce local water resources. Geographic information systems (GIS) based approaches are being developed as tools to address this concern, focus efforts, and minimize issues associated with stormwater run-off. Understanding influencing watershed factors such as topography, precipitation, soil properties, vegetation, and drainage area, etc. is key to the identification of problem areas for stormwater run-off. These areas, if geographically identified, may be evaluated for stormwater control measure (SCM) development. Structural and nonstructural SCMs are techniques, methods, activities, and practices that have the primary objective of controlling run-off volume and preventing or reducing the quantity of pollutants present in runoff water before reaching natural water systems (Richardson et al. 2018). Implementation of SCMs in strategic locations has the potential to prevent catastrophic run-off events and present a solution for decreasing surface waterway contamination.
Successful implementation of SCM strategies should be based on the unique geographic characteristics, local requirements of each project, and watershed-specific needs. Incorporating all of these potential decision making criteria for planning optimal SCM sites necessitates a system which allows for the comparison and prioritization of diverse criteria on a large scale. An extensively used technique for improving SCM efficiency is incorporating multicriteria decision support (MCDS) systems. MCDS involves the use of scientific and mathematical tools to evaluate complex trade-offs among different environmental and socioeconomic factors. Analytic hierarchy process (AHP) is one such commonly adopted MCDS technique. AHP is a measurement method to derive a scale of relative importance of alternatives or criteria from pairwise comparisons  . AHP input can represent actual quantitative measurements or designations of subjective judgment. A more evolving concept is the integration of MCDS and geographic information system (GIS) processes into SCM suitability studies, which allows for the pinpointing of optimal locations for SCM development.
An evaluation of watershed vulnerability in Bernalillo County, NM was performed by Amankwatia  . The research herein is based on this work, but uses fuzzy AHP (FAHP) to incorporate additional subjective judgment input regarding Bernalillo watershed processes. These subjective inputs were obtained from surveyed experts in the fields of geology and hydrology. Input from local experts was used to obtain a watershed vulnerability assessment of Bernalillo County comparable to the model produced in 2015.
One potential benefit of GIS-based MCDS models is that they can be easily implemented in regions where precise spatial data are not readily available without having to rely on expensive and time-consuming field surveys  . A second potential benefit of these methods is that they use local expert input to customize solutions to regionally specific watershed phenomena. Information from MCDS models can prove useful for overall watershed management and resource planning. This study combines specific soil erosion and infiltration criteria within a GIS-based fuzzy AHP framework to develop a watershed vulnerability map for Bernalillo County, New Mexico. The output could be beneficial for exploring optimization of SCM approaches to watershed management for targeted vulnerable areas within the county.
2. Stormwater Control Measures in the Desert Southwest
Despite widespread support of stormwater control measures (SCM) and projects in the United States, there has been little representation of such approaches in the desert Southwest  . The lack of implementation of SCMs in the Southwest is closely tied to the perception that these regions experience little or no rainfall and, hence, the expense and burden of SCMs seem unjustifiable. However, Gautam et al.  argue that one important reason for establishing SCM in desert environs is that limited, but high-intensity, episodic rainfall does occur. Such events can be a major contributing factor to the elevated levels of harmful pollutants found in natural watercourses. Precipitation in New Mexico is highly variable with some months exhibiting little to no precipitation and summer monsoon storms bringing torrential downpours. These storms typical supply 50 percent of New Mexico’s annual precipitation during the months of July, August, and early September  .
In the semi-arid environment of New Mexico, stormwater is both an invaluable and limited resource. Limited availability of water resources and sustainability problems in this region actually intensify the need for the adoption of SCMs. Gautam et al.  suggested that given the rapid urbanization and severe water availability issues in the desert Southwest, the best strategy for implementing SCMs in such areas should be based on sustainable land use, water conservation, and water reuse techniques. The authors also stated that one critical SCM for long-term sustainability is through low-impact development by reducing impervious areas and limiting directly connected impervious surfaces.
MCDS Approach to Watershed Vulnerability and SCM Using AHP and FAHP
The growing interest in the use of analytical hierarchy process (AHP) in SCM suitability studies is associated with the ability of AHP to deal with imprecision in environmental data. As reported by Young et al.  , use of AHP supported with input from GIS can provide an objective, mathematically-based alternative to the existing subjective SCM approaches. Fuzzy AHP (FAHP) provides the same benefits and incorporates the ability to crowdsource decision prioritization. AHP has been successfully applied in assessing regional scale urban flood hazard   . Research has also been carried out on the effective use of AHP and GIS to manage surface-water and groundwater resources. For example, Agarwal et al.  suggested that AHP can be effective to evaluate groundwater potential risks. A more recent attempt suggested that AHP techniques provide an effective framework for reassessing vulnerability of groundwater resources and aquifer systems  . Additional AHP studies have focused on soil erosion risk  and have shown it to be an effective technique of predicting soil erosion risks at precisions equivalent to field observations  . Fuzzy AHP adds an additional element of flexibility in data integration. The objective of this study was to implement a composite FAHP model based on soil erosion and infiltration with prioritization of decision criteria determined by a fuzzy analysis of expert knowledge of the study region. The purpose of using a weighted decision system is to provide a geographic framework for watershed vulnerability assessment, which is uniquely suited to the study watershed. This technique results in a vulnerability map, which can inform SCM implementation predicated by expert knowledge of the interest area.
3. Factors Influencing Watershed Vulnerability
Watershed vulnerability, defined as high soil erosion and low infiltration potential, is a complex process to model. Several watershed factors, such as slope gradient, soil texture, soil erodibility, soil permeability, land use and land cover, vegetation cover, and drainage density, influence the development of the watershed vulnerability model. The characteristics and effect of these factors on soil erosion and infiltration are outlined in this section.
3.1. Slope Gradient
Slope gradient (S), usually expressed as percent slope or degree slope, is used to represent the degree of variability in elevation and is found to strongly influence the velocity and direction of stormwater runoff. Generally, areas of watersheds with steeper slopes are characterized by high soil erosion susceptibility, whereas areas with mild or gentle slopes are suitable for infiltration and water harvesting practices. The SPA Digital Terrain Model (SPADTM) developed by the United States Army Corps of Engineers (USACE) Geospatial Center (~4.5 m resolution) was obtained with special permission from the New Mexico Department of Geology and Mineral Resources and converted to slope gradient.
3.2. Soil Texture
Soil texture, defined as ratios of sand and clay, is known to have a strong influence on the rate of infiltration and runoff potential of watersheds. The National Resources Conservation Service (NRCS) has identified four distinct soil groups (A, B, C, and D) based on their soil textures. Group A soils are generally comprised of sand and gravel and have higher infiltration potential, whereas Group D soils have a high clay content and lower infiltration potential. Soil Survey Geographic Database (SSURGO) datasets are collected by National Cooperative Soil Survey, distributed by the United States Department of Agriculture (USDA). Hydraulic soil class dataset was downloaded from the USDA database and accessed using Soil Data Viewer.
3.3. Soil Erodibility
The soil erodibility (K) factor is a quantitative soil property that plays a major role in determining the extent to which soil particles are eroded or transported by runoff water. Parameters such as soil moisture, soil permeability, and organic matter content influence the degree of soil erodibility and stability. Soils with low K values are highly stable compared with those with higher K values. K-factor values were obtained through SSURGO using Soil Data Viewer.
3.4. Vegetation Cover
The normalized difference vegetation index (NDVI) is a widely used standardization parameter to represent vegetation cover. NDVI values typically range from 1.0 to −1.0 with non-vegetated areas such as streams and rivers producing small or slightly negative values, and vegetated areas producing values ranging from 0.4 to 1.0  . Landsat 8 data was downloaded for the area covering Bernalillo County using Earth Explorer provided by the United States Geologic Survey (USGS). Data was collected from cloud-free fly-overs in June 2018. June represents the driest month of the year in New Mexico and conditions which would be present at the start of the annual monsoon season  . Landsat 8 data was converted to reflectance, corrected for sun angle, and NDVI was calculated from the red and near infrared data bands.
3.5. Soil Permeability
Soil permeability, or more correctly saturated hydraulic conductivity (Ksat), is an important physical soil property that refers to the ability of soils to conduct water. High Ksat values represent the portions of watersheds that are highly drained and are suitable for infiltration practices. Conversely, lower Ksat values correspond to watershed regions that experience high runoff and are prone to soil erosion. Ksat values were obtained through SSURGO using Soil Data Viewer.
3.6. Land Use
Land use for purposes herein refers to the use of land for human activities. Land use directly translates to the amount of water penetration into soil layers or being produced as runoff water. For example, rangeland or grassland will have different infiltration or erosion potential than land which is used for industrial purposes. Land use data was downloaded from the USDA Geospatial Data Gateway and classified into risk categories.
Area precipitation characteristics (type, duration, intensity, and amount) are highly important in drainage projects because of their direct influence on the volume of generated runoff. Precipitation characteristics play a significant role in soil erosion. In civil and environmental engineering, storm classification is used to determine appropriate stormwater control measures. Storm intensity and duration are important factors in determining the volume of water precipitated by a storm and the erosional or run-off risk that the storm will pose to a watershed. In semiarid watersheds, high-intensity precipitation storms are believed to be the main component for soil erosion. Precipitation frequency (PF) estimates are used as design criteria for a wide variety of civil and environmental engineering applications, such as the design of hydraulic structures such as culverts, roadway drainages, and bridges. PF estimate data were obtained from NOAA’s Hydrometeorological Design Studies Center—Precipitation Frequency Data Server (PFDS). This dataset consists of spatially distributed point PF estimates formatted as a georaster of the Southwestern United States. A high intensity, high volume storm is likely to pose a higher risk for erosion and storm water run-off. For this project, a precipitation depth for a 10-year, 6-hour storm was used to approximate precipitation for a severe monsoon storm, but alternate storm design frequencies could easily be substituted in the FAHP-GIS method.
3.8. Drainage Density
Drainage density is a quantitative indicator that represents the ratio of the total length of streams to the total watershed area. Low drainage density values are favorable for infiltration, whereas larger values indicate high runoff potential. High drainage density values are an indication of high soil erosion vulnerability, whereas lower values signify minor risk. The linear drainage file was downloaded from USDA Geospatial Data Gateway and drainage density was determined using the ArcMap line density function.
4. Site Description
This study was conducted using maps and digital data from Bernalillo County, New Mexico. The county is located in North Central New Mexico, centered at approximately longitude 106˚41'56.40''W and latitude 35˚3'14.4''N. Bernalillo County has a total area of approximately 3023 km2 (1167 mi2). Bernalillo County is the most highly populated county in the State of New Mexico, and is home to the City of Albuquerque, its county seat. The region is covered by shrubs and scrubs (71%), deciduous, mixed, and evergreen forests (18%), associated urban (9%), and open water (1%). The two dominant soil orders are Entisols (47%) and Aridisols (32%)  . Albuquerque’s urban sprawl consisting of residential, commercial, military, and industrial land perpetuates a high risk of surface run-off events. Urban surface conditions result in contamination of surface water run-off, which poses risk nearby waterways and nearby agricultural development. Limited water availability, unique soil characteristics, dominant shrub-land vegetation, and intricate climate, coupled with a forecasted 35% population growth over the next 20 years  , makes Bernalillo County a suitable case study area for the current research.
5. Modeling Approach
The overall modeling approach used for this project is based on Fuzzy AHP multicriteria decision making utilizing expert opinions collected by survey (Figure 1). The principal benefit of this method as an approach for evaluating watershed vulnerability is its ability to objectively and simultaneously consider an unlimited number of relative criteria.
For this project, two submodels are used to evaluate watershed vulnerability.
Figure 1. Fuzzy AHP Model development.
One submodel assesses soil erosion potential and the other assesses infiltration potential; each model incorporates multiple criteria. Coupled, these processes influence watershed vulnerability for stormwater run-off. The premise herein is that areas with high soil erosion potential and low infiltration potential are most vulnerable.
5.1. Soil Erosion Submodel
Soil erosion poses an imminent threat to infrastructure and property. As recently as 2018, heavy monsoon rains in the area flooded homes, eroded dirt roads, blocked access to homes and businesses, and left thick blankets of sediment on agricultural crops. The City of Corrales, New Mexico in Sandoval County adjacent to Bernalillo County saw up to three feet of silt and sand deposited in areas and extreme erosion in areas during the July 2013 floods  . The soil erosion submodel combines precipitation, soil erodibility, land slope, land use, and vegetation cover.
5.2. Infiltration Submodel
Infiltration is the entry of water into the soil surface and its subsequent vertical motion through the soil profile. Areas which experience low soil infiltration will be more highly vulnerable for surface run-off events. The soil infiltration submodel uses drainage density, hydrologic soil classification (HSC), land slope, land use, and saturated hydraulic conductivity as criteria to establish run-off potential.
6. Tools and Data Processing
A variety of software tools and databases are required to facilitate implementation of an AHP approach to assess watershed vulnerability. The following sections summarize the tools and databases used for this research.
ArcGIS is geographic information software system, which provides an extensive toolset to manage multiple layers of spatially-distributed geographic data. Modeling approaches based on GIS have been shown to be successful in assessing regional-scale soil erosion risks  . GIS-based models have also been extensively used in evaluating groundwater vulnerability in order to aid policy formulation. ArcGIS provides an infrastructure for making maps and working with large data sets of geographic and remotely sensed data. In this study, ArcMap was used for a range of geoprocessing tasks including spatial referencing, resampling, mosaicking, clipping, reformation, and classifying geographic data. ArcGIS’s raster calculator function was used to combine thematic data layers, creating erosion and infiltration models. Data analysis and decision support tools were critical in creating custom watershed vulnerability maps.
6.2. Soil Data Viewer
The NRCS Soil Data Viewer is an ArcGIS extension for analyzing soil data and building soil-based thematic maps. The program enables map makers to interpret soil attributes and properties without having to deal with the complexity associated with querying and processing soil databases. Soil Data Viewer was used to create and export thematic layers such as soil texture, hydrologic soil groups, soil depth, and saturated hydraulic conductivity maps for the study area.
6.3. Data Processing
All layers were projected in ArcGIS to North American Datum (NAD) 1983 UTM zone 13 projection prior to assimilation into the model. A grid size of 30 × 30 (m) pixels was chosen for all spatial layers to cover the study area extent. All feature datasets were converted to raster format and all raster datasets were spatially resampled to 30 (m) resolution using the Data Management toolbox of ArcGIS. The Data Management toolbox was then used to clip the input thematic layers to conform to the spatial extent of the study area. In order to produce the soil data maps (including soil texture K-factor and Ksat), SSURGO data were loaded into Soil Data Viewer and individual layers were selected for export.
6.4. Expert Survey
A survey conducted of experts in the fields of geology and hydrology was conducted to generate comparison values for each of the criteria influencing soil erosion and infiltration. Surveys were distributed to experts of three local organizations familiar with watershed processes in Bernalillo County, New Mexico. These organizations included the New Mexico Institute of Mining and Technology, the New Mexico Bureau of Geology and Mineral Resources, and the United States Geologic Survey (USGS) in Albuquerque. Participants were provided with a description of each of the criteria used in developing the watershed vulnerability model. The survey was constructed in two parts. Part 1 asked participants to compare criteria influencing erosion and Part 2 compared criteria influencing infiltration. Participants compared each influencing criterion (A) to each other criterion (B) and evaluated whether (in their opinion) criterion A was equally as, slightly more, moderately more, strongly more, or absolutely more important than criterion B. Or conversely, participants indicated whether criterion B was equally as, slightly more, moderately more, strongly more, or absolutely more important than criterion A. In this distributed survey (soil erosion and infiltration, respectively), each of the five influencing criteria is compared against each other criteria for the same phenomenon. Table 1 shows the relative scale of importance used for the respective pair-wise comparisons.
7. Development of Fuzzy AHP Model
7.1. Initial Comparison Matrix
Unlike “crisp” AHP, each comparison criteria is given as a set of possible values or distribution. A triangular distribution was applied to the selected fuzzy numbers indicating importance. The triangular fuzzy number (TFN) is applied through a triangular membership function (μA) from 0 to 1 defined by three real numbers, expressed as a lower value, mean value, and an upper value (l, m, and u), where
Applying a triangular distribution to the numbers results in a minimum, mean, and maximum value for each importance comparison. The composite initial matrix is transformed into an integrated fuzzified pair-wise comparison matrix (FPCM) using a geometric mean of fuzzy triangular functions  . An example of the integrated FPCM for soil erosion is provided as Table 2, where C1-C5 represents the five criteria: precipitation, K factor, slope, land use, and vegetation, respectively. A similar integrated FPCM was developed for infiltration for its five criteria.
Table 1. Fuzzy judgement for pair-wise comparisons.
Table 2. Integrated fuzzified pair-wise comparison matrix for soil erosion.
Figure 2. Expert survey initial comparison matrix for soil erosion.
7.2. Consistency Evaluation
The principal behind using FAHP as a multicriteria decision making process (MCDS) relies on the consistency of input for the pairwise comparison matrix. In order to ensure that AHP judgments are accurate and reliable, a consistency ratio (CR) for each submodel was calculated according to method proposed by 
Figure 3. Expert survey initial comparison for infiltration.
where n is the number of criteria in the comparison, λ max is the maximum Eigenvalue, and RI is a random consistency index based on the number of criteria. The more consistent a comparison matrix is, the closer the value of computed maximum Eigenvalue (λmax) is to the number of criteria. Ideally, the consistency ratio should be less than 10% for each submodel  . For both erosion and infiltration the consistency ratio was found to be less than 0.02% or 2%. Because there are five comparison criteria for each submodel, the maximum Eigenvalue should approach a value of 5 for both. The calculated consistency ratios and maximum Eigenvalues for the erosion and infiltration sub models are CRerosion = 0.0161; λmax erosion= 5.07 and CRinfiltration = 0.0146; λmax infiltration = 5.06, respectively.
7.3. Defuzzification of Fuzzy Pair-Wise Comparison Matrix
Determining the crisp weight of each criteria from the integrated fuzzy comparison matrix uses the fuzzy extent analysis (FEA) method by Chang  with the Wang et al.  modification of the normalization process applied to calculate the values in the n x n fuzzy synthetic extent matrix (FSEM). The centroid method was then used to calculate crisp criteria weights for each of the influencing criteria. The centroid method calculates the average for row i in the FSEM as (li + mi + ui)/3, which are summed to give a total for n rows. Respective crisp criteria weights (wi) are determined by dividing row i centroid average by the total.
8. Vulnerability Rating of Fuzzy AHP Model
A five-category rating system was used to classify vulnerability for this research. Each category is used to define a condition of vulnerability from “not vulnerable” to “extremely vulnerable”. Vulnerability categorizations are as follows: 1) Not Vulnerable; 2) Slightly Vulnerable; 3) Moderately Vulnerable; 4) Highly Vulnerable; and 5) Extremely Vulnerable.
The next step is to assign vulnerability ratings for each criterion. Ranges were established to classify each criterion into risk categories for vulnerability. Rating categories were modeled after the works of Richardson and Amankwatia  . Categories were defined to reflect the level of vulnerability a criterion poses for stormwater run-off. For example, a slope of 0 - 2 degrees is considered not vulnerable, 2 - 5 degrees is slightly vulnerable, 5 - 10 degrees is moderately vulnerable, 10 - 18 degrees is highly vulnerable. and any slope greater than 18 degrees is extremely vulnerable. Each class is attributed to a value 1 - 5. All of the weight values for each submodel are scaled to add up to 1, such that after all the criteria are weighted and combined, the maximum vulnerability value remains on a scale from 1 to 5. Also note that land cover classifications “herbaceous”, “hay and pasture”, and “cultivated crops” were reclassified for this research as “cultivated”. Criteria vulnerability ratings and the crisp weights derived from fuzzy AHP are presented in Table 3 and Table 4.
9. Results and Discussion
9.1. Evaluation of Watershed Vulnerability
The weighted overlay of thematic map layers created a map of Bernalillo County, which indicates erosion vulnerability and infiltration vulnerability due to a storm event. The map of soil erosion potential for Bernalillo County is shown in Figure 4. Infiltration potential for Bernalillo County is shown in Figure 5. These two maps were combined to produce a map of overall watershed vulnerability for run-off related events. The scaled watershed vulnerability of Bernalillo County is given in Figure 6. Warm or red and orange-colored areas indicate regions of high vulnerability and cool-colored areas or blue and green areas indicate areas of lower vulnerability. For this study, the relationship of erosion and infiltration processes is combined with equal weighting to produce a map of combined vulnerability.
Weighting factors for each influencing parameter included in the study of erosion potential and infiltration potential are shown in Table 3 and Table 4, respectively. NDVI is a dominant criterion influencing soil erosion, followed by precipitation. Saturated hydraulic conductivity and soil texture are equally dominant as influencing infiltration.
Areas of high infiltration vulnerability (or low infiltration) and areas of high erosion vulnerability designate areas are highly or extremely vulnerable for the occurrence of run-off related events. These high vulnerability areas, shown as red in Figure 6, are areas were stormwater control measures (SCM) can be focused for the protection of property and infrastructure as well as prevention of surface waterway contamination.
Table 3. Soil erosion criteria vulnerability rating.
Table 4. Infiltration criteria vulnerability rating.
Figure 4. Erosion vulnerability for Bernalillo County, New Mexico.
Figure 5. Infiltration vulnerability for Bernalillo County, New Mexico.
Figure 6. Combined watershed vulnerability for Bernalillo County, New Mexico.
One benefit of working within a GIS platform is the potential for inclusion of areas of special interest in SCM implementation. As an example, infrastructure such as roadways may be added to a vulnerability map to identify target areas of importance for SCM development. Figure 7 shows watershed vulnerability for a section of eastern Bernalillo County with overlaid roadways, wherein there are a few neighborhoods and sections of highway which are designated as high vulnerability areas. High vulnerability areas occurring at a great distance from target features may be designated as a low priority for SCM. In addition, specialized targets such as roadways, special properties, or contamination hotspots may be introduced to the model in order to assess which vulnerability areas are of the highest priority for development of SCM.
Most of the area of Bernalillo County represents moderate watershed vulnerability. The distribution of vulnerability values for the Bernalillo County area (in model pixels) is given in Figure 8. The mean vulnerability value of the model was determined to be 3.07 (on a scale of 1 to 5) with a standard deviation of 0.28. The majority of the pixels lie within two standard deviations of the mean. The distribution indicates that the “not vulnerable” and “extremely vulnerable” categories do not represent a significant portion of total land area within Bernalillo County.
9.2. Advantages of Fuzzy AHP Model Approach
Multicriteria decision support approaches such as FAHP provide a relatively easy and customizable way of determining watershed vulnerability distributions on a regional scale. These methods can provide several strategic benefits to decision makers including the ability to incorporate a variety of quantitative data, the experience of local experts, and project specific objectives. An advantage of FAHP is the method’s ability to deal with data uncertainty and lack of sampled erosion, infiltration, or other watershed data. However, developing pairwise comparisons for AHP hinge on subjective judgment and can create a great deal of variability in model decision making priorities. FAHP provides the opportunity to combine the judgement of multiple information sources, to provide additional stability to the normalized priority vector values. A sensitivity analysis can be performed wherein the input data are slightly modified in order to observe
Figure 7. Roadway overlay on combined vulnerability map (zoomed-in view).
Figure 8. Scaled combined vulnerability value distribution.
the impact on the weighted priority results. If the rank order does not change, the results are said to be robust; otherwise, they are sensitive.
The FAHP approach presented herein offers several benefits in developing an understanding of regional-scale watershed vulnerability. The integration of FAHP and GIS techniques into watershed vulnerability studies for SCM suitability can be easily implemented in regions where precise spatial data are not readily available without having to rely on expensive and time-consuming field surveys. Information from MCDS models can prove useful for large scale watershed management and resource planning, narrowing down areas which can later be tested for SCM suitability.
This study combines specific soil erosion and infiltration criteria within a GIS-based fuzzy AHP (FAHP) framework to develop a watershed vulnerability map for Bernalillo County, New Mexico. Quantifying a regional-scale risk index map can provide spatially relevant information for engineering management and decision making. In view of this, a 30 × 30 (m) watershed vulnerability risk map was generated for the Bernalillo County, New Mexico using an FAHP methodology based on judgements of experts on local watershed processes. The approach incorporated drainage density, slope, precipitation, vegetative cover, land use, soil erodibility, soil texture, and soil permeability as the primary factors contributing to watershed vulnerability. The final map was produced by merging two submodels representing soil erosion and infiltration potential using the ArcGIS software. Analysis of the vulnerability map identifies that the most problematic areas are present in the Eastern and Northwestern portions of the study area. However, based on the analysis less than 1% of the total area lies within the lowest and highest vulnerability categories with the majority of land area centered around moderate vulnerability.
Results were compared with those using AHP derived weights from the Amankwatia (2015) effort. Both weighting methods showed similar regional vulnerability trends. The percentages associated with each vulnerability category were different, especially with the not vulnerable and extremely vulnerable categories. The previous study showed approximately 4% of the total area was within the lowest and highest vulnerability classifications. This result is expected due to the differing weights (AHP versus FAHP) applied to each criterion. Further analysis is required to determine model sensitivity to each influencing criterion.
Fuzzy AHP is, therefore, an easily applied numerical method for evaluating criteria weights for MDCS systems. It allows for the quantitative evaluation of data from multiple, sometimes incomparable, sources by creating a system of numeric weights for a quantitative comparison and data incorporation. The method used in this process, couples two key hydrologic phenomena to assess watershed vulnerability (erosion and infiltration). The purpose of the application of this process is to identify vulnerable areas in need of implementation of stormwater management measures. FAHP relies on the experience of local experts and allows for the inclusion of additional evaluation datasets (criteria) customizable to an area’s regionally specific watershed phenomena. This method produces site-specific vulnerability estimations that can be used for local discussions about individual risk. Information from the watershed vulnerability model can be useful for managing new structural or nonstructural SCM projects or improving existing ones to mitigate impact to stormwater contamination and run-off volume for identified vulnerable areas. Better prediction of problems allows for better planning and management of hazards, which contribute to watershed vulnerability.
11. Recommendation for Future Work
Future work should include field investigations of regional susceptibility to stormwater run-off and erosion in order to validate the usefulness of the model and for use in stormwater and watershed management. In addition to photographic documentation from various field sites, and quantitative erosion and infiltration data points could be used to adjust the vulnerability model and test for accuracy using a confusion matrix approach  .
Impacts of variation in spatial resolution of datasets should also be explored. While most datasets including SPADTM, SSURGO, Landsat 8, and USDA landuse data had a resolution of 30 (m) or less, the climate precipitation data was provided in an 800 (m) raster. This results in reduced spatial variation in the precipitation data when it is incorporated into the model and may impact interpreted vulnerability by area. Amankwatia  used a 100 yr 24 hr precipitation depth criterion versus a 10 yr 6 hr precipitation depth criterion herein. NDVI data appropriate to the critical premonsoonal period should be extracted for multiple years and averaged to better address soil erosion potential. Amankwatia  based the NDVI criterion on an average annual maximum from 1995-2009, wherein this study used one premonsoonal month (June 2018). Additional expert opinion could also be solicited for inclusion into the FAHP decision criteria weighting process.
Future studies should explore the relative weighting of erosion and infiltration submodels as they pertain to watershed vulnerability. For this study infiltration and erosion submodels were weighted with equal importance. Further research should be conducted to determine the respective influence of each submodel on watershed vulnerability.
In order to improve on the functionality and applicability of the model, other factors such as public awareness and economic indicators could be incorporated into the model. Integration of watershed vulnerability with infrastructure assessment and property maps can help to establish where implementation of SCM would be most critical. An example of this was provided with roadway overlays. Due to the large-scale nature of this data modeling process, it is feasible to extend this study beyond the boundaries of Bernalillo County to model the State of New Mexico.
The high-resolution digital elevation data (SPADTM) was obtained with special permission from the New Mexico Department of Geology and Mineral Resources.
 Blackwell, W. (2018) Mitigation. New Mexico State Hazard Mitigation Plan 2018. New Mexico Department of Homeland Security and Emergency Management (NM DHSEM).
 Richardson, C.P. and Amankwatia, K. (2018) GIS-Based Analytic Hierarchy Process Approach to Watershed Vulnerability in Bernalillo County, New Mexico. Journal of Hydrologic Engineering, 23, Article ID: 04018010.
 Amankwatia, K. (2015) A GIS-Based Multi-Criteria Decision Support Approach to Stormwater Best Management Practices (SCMS): Identifying Potential SCM Vulnerable Sites for Effective Water Conservation and Water Reuse in Bernalillo County, New Mexico. M.S. Thesis, New Mexico Institute of Mining and Technology, Socorro.
 Gautam, M.R., Acharya, K. and Stone, M. (2010) Best Management Practices for Stormwater Management in the Desert Southwest. Journal of Contemporary Water Research & Education, 146, 39-49.
 NOAA (2018) North American Monsoon Highlights.
 Young, K.D., Tamim, Y., Randel, D.L., David, K.F. and David, L.H. (2010) Application of the Analytic Hierarchy Process for Selecting and Modeling Stormwater Best Management Practices. Journal of Contemporary Water Research & Education, 146, 50-63.
 Wang, Y., Li, Z., Tang, Z. and Zeng, G. (2011) A GIS-Based Spatial Multi-Criteria Approach for Flood Risk Assessment in the Dongting Lake Region, Hunan, Central China. Water Resource Management, 25, 3465-3484.
 Agarwal, E., Agarwal, R., Garg, R. and Garg, P. (2013) Delineation of Groundwater Potential Zone: An AHP/ANP Approach. Journal of Earth System Science, 122, 887-898.
 Al-Abadi, A. and Al-Shamma’a, A. (2014) Groundwater Potential Mapping of the Major Aquifer in Northeastern Missan Governorate, South of Iraq by Using Analytical Hierarchy Process and GIS. Journal of Environment and Earth Science, 4, 125-149.
 Kachouri, S., Achour, H., Abida, H. and Bouaziz, S. (2014) Soil Erosion Hazard Mapping Using Analytic Hierarchy Process and Logistic Regression: A Case Study of Haffouz Watershed, Central Tunisia. Arabian Journal of Geosciences, 8, 4257-4268.
 Bulut, G., Cal, M., Richardson, C. and Gallegos, J. (2012) A GIS-Based Soil Erosion Risk Map for New Mexico. World Environmental and Water Resources Congress, Reston, May 2012, 3754-3763.
 Richardson, C.P. and Amankwatia, K. (2019) Assessing Watershed Vulnerability in Bernalillo County, New Mexico Using GIS-Based Fuzzy Inference. Journal of Water Resource and Protection, 11, 99-121.
 Mitra, B., Scott, H., Dixon, J. and McKimmey, J. (1998) Applications of Fuzzy Logic to the Prediction of Soil Erosion in a Large Watershed. Geoderma, 86, 183-209.
 Forman, E. and Peniwati, K. (1998) Aggregating Individual Judgments and Priorities with the Analytic Hierarchy Process. European Journal of Operational Research, 108, 165-169.
 Krois, J. and Schulte, A. (2014) GIS-Based Multi-Criteria Evaluation to Identify Potential Sites for Soil and Water Conservation Techniques in the Ronquillo Watershed, Northern Peru. Applied Geography, 51, 131-142.
 Wang, Y.M., Elhag, T.M.S. and Hua, Z.S. (2006) A Modified Fuzzy Logarithmic Least Squares Method for Fuzzy Analytic Hierarchy Process. Fuzzy Sets and Systems, 157, 3055-3071.