Using GIS to compute a Least-Cost Distance Matrix: A Comparison of Terrestrial and Marine Ecological Applications

Patrick N. Halpin1 and Andrew G. Bunn2

Presented at the Twentieth Annual Esri User Conference, June 2000

1. Corresponding author: Landscape Ecology Laboratory, Nicholas School of the Environment, Duke University, Box 90328, Durham NC 27708-0328

2. Mountain Research Center, Montana State University, Bozeman MT 59718

Abstract

Modeling habitat connectivity on both landscapes and seascapes requires the analysis of potential pathways between habitat patches. In order to assess the importance of individual pathways, a complete set of possible paths must first be developed. In terrestrial situations, least-cost path algorithms can be used in an iterative manner to create a set of all potential paths between patches. In marine examples, directionality due to ocean currents requires a vector network approach to create the relative paths between habitat patches.
An example of the terrestrial problem is presented using wetland forest patches in Eastern North Carolina. In order to analyze a full set of potential paths between all patches we developed an ArcInfo AML macro to compute a distance matrix D whose elements dij are the functional distances between points i and j. Because many ecosystem elements do not move through the landscape in a straight line, we measured dij using least-cost paths. We used grid functions inside a macro to iteratively loop through a spatial array and compute d ij for each unique pair of nodes in the array. These functions are similar to Euclidean distance functions but instead of working in geographical units they work in cost units. We demonstrate the usefulness of this method in for two wetland habitat species (mink and prothonotary warbler) in Eastern North Carolina.
An example of the marine connectivity problem is presented using a network of connecting flow paths between marine reef habitats. This example is portrayed using a set of constructed artificial reefs and ship wrecks off the coast of North Carolina. In this example we use an easterly wind driven current representing general conditions in early spring (March-April). Changes in current directions and velocities require the development of new network solutions for each current regime.
In addition to requiring two different types of path analysis approaches, marine and terrestrial analyses also differ in the development of the distance matrix. In the terrestrial example, species traveling between patches are expected to move equally well in either direction, so the distance between patches i and j are equivalent to j and i. This results in needing to calculate half of the possible distance matrix. In the marine example, the entire distance matrix must be calculated because the cost-distance from i to j may be different than j to i due to ocean current impedance. In addition, the terrestrial network can be assumed to be entirely connected into a single patch at some travel distance, while the marine network is never considered to be a single patch at any distance.
The methods contrasted here provide an objective framework for developing and assessing spatial connectivity for patchy environments in both nterrestrial and marine ecosystems.

Introduction

Modeling habitat connectivity on both landscapes and seascapes requires the analysis of potential pathways between habitat patches. In order to assess the importance of individual pathways, a complete set of possible paths must first be developed. In terrestrial situations, least-cost path algorithms can be used in an iterative manner to create a set of all optimal paths between patches. In marine examples, directionality due to currents requires a vector network approach to create the relative paths between habitat patches.
 
  1. Figure 1. Location of (a) the terrestrial study area wetland forest patches in the vicinity of the Alligator River NWF; and (b) the marine study area artificial reefs and ship wrecks in the Cape Lookout region of the South Atlantic Bight.
In order to test the development of connectivity models for both terrestrial and marine environments, two study sites were selected in close proximity to one another. The terrestrial study area represented a patchy environment of wetland forests in the Alligator River National Wildlife Refuge in Eastern North Carolina. The marine study area was represented by patchy artificial reefs and wrecks in the Cape Lookout region of the South Atlantic Bight off the coast of North Carolina.
Until now, approaches to building models of patchy environments have focused primarily on two types of spatial data. The two data structures most familiar to conservation biologists are lattices either as coverages of vectors (polygons) or as raster grids. In the case of vector coverages landscapes are commonly represented as island models consisting of habitat patches in a non-habitat matrix. With raster grids landscapes are modeled as habitat mosaics, where each cell is assigned a value indicating habitat suitability.
We demonstrate the utility of a unfamiliar type of lattice, the graph (Harary 1969), in determining landscape connectedness using focal species analysis in an island model. A graph represents a binary landscape of habitat and nonhabitat where patches are described as nodes and the functional connections between them as edges. Graph theory is a widely applied framework in geography, information technology, and computer science. It is primarily concerned with maximally efficient flow or connectedness in networks (Gross and Yellen 1999). To this end, graph-theoretic approaches can provide powerful leverage on ecological processes concerned with connectedness as defined by dispersal (Urban and Keitt, in review). In particular, graph theory has great potential for use in applications in a metapopulation context.

Terrestrial Distance Matrix Example

Here we explore landscape connectedness for two focal species in the Southeastern Coastal Plain, American mink (Mustela vison) and prothonotary warbler (Protonotaria citrea). Both depend on wetlands and bottomland forests and both are species of concern to conservation managers. The Alligator River and Pocosin National Wildlife Refuges (NWR) in North Carolina and the surrounding counties are an ideal study site for this approach. It is a large area with significant human-induced changes to the pattern-forming agents of the landscape. The alteration of the hydrologic regime by diversion and conversion of wetlands has changed the physical template of the landscape. Biotic processes have been greatly altered by extirpation of native species and the introduction of exotic species. The area contains some of the last wildlands in North Carolina's coastal plain and is under increasing development pressure. A well-developed system of reserves is in place, but we are not aware of any other study that examines functional connections within and among those areas.

Terrestrial Study Area and Methods

Study Area
The terrestrial study focuses on the Alligator River NWR and surrounding counties in the Coastal Plain of North Carolina (35°50'N; 75°55'W). It is a riverine and estuarine ecosystem with an area of almost 580,000 ha and over 1400 km of shoreline. The area is rich in wildlife habitat, dominated by the Alligator River NWR, the Pocosin Lakes NWR, Lake Mattamuskeet NWR, Swanquarter NWR, and a variety of other federal, state, and private wildlands. The western portion of the study area is largely in small (<1000 ha) private holdings. The area is bounded on the east by Croatan Sound. The Pamlico River and Pamlico Sound form the southern boundary and the Chowan River delta and Albemarle River form the northern boundary. The western boundary of the study area is formed by the Beaufort County line. The area encompasses much of Beaufort and Dare Counties, and all of Hyde, Washington, and Tyrell Counties, North Carolina (figure 2).
 
 
  1. Figure 2. Terrestrial Study area with hydrography, major roads, and suitable habitat.
The vegetation is characterized by the Southern Mixed Hardwoods forest community. The area has many diverse vegetation types, including fresh water swamps, pine woods, and coastal vegetation. In the upland community dominant species include many types of oak (Quercus spp.), American Beech (Fagus grandifolia) and evergreen magnolia (Magnolia grandiflora). Mature stands may have 5 to 9 codominants. The wet lowlands are dominated by bald cypress (Taxomodium distichum). The pine woods are dominated by longleaf pine (Pinus palustris), but loblolly pine (P. taeda) and slash pine (P. elliottii) are also important (Vankat 1979).
Focal species
Because connectivity occurs at multiple scales and multiple functional levels (Noss 1991), We have chosen two focal species to apply a graph-theoretic approach to connectedness. Focal species analysis is an essential tool for examining connectivity in a real landscape as individual species have different spatial perceptions (O'Neill et al. 1988). The American mink and the prothonotary warbler are appropriate candidates for focal species analysis as they share very similar habitat but have different ecological requirements and fall into different categories as focal species. Both species are wetland dependent and indicators of wetland quality and abundance in a landscape, both are charismatic and enjoy flagship status. Furthermore, as meso-predators mink have small but important roles as umbrella species and keystone species (Miller et al. 1999).
American mink are meso-level, semiaquatic carnivores that occur in riverine, lacustrine and palustrine environments (Gerell 1970). They are chiefly nocturnal and their behavior largely depends on prey availability. They have a great deal of variation in their diet according to habitat type, season and prey availability (Dunstone and Birks 1987). Muskrats (Odantra zibethicus) are a preferred prey item (Hamilton 1940, Wilson 1954), but mink diets in North Carolina are composed of aquatic and terrestrial animals as well as semiaquatic elements (e.g. waterfowl) (Wison 1954). In the southeast they have home ranges on the order of 1 ha and a dispersal range of roughly 25 km (Nowak 1999).
Prothonotary warblers are neotropical migrants that breed in flooded or swampy mature woodlands. They have two very unusual traits in wood warblers in that they are cavity nesters and prefer nest sites over water. They are forest interior birds that experience heavy to severe parasitism by brown-headed cowbirds (Molothrus ater) (Petit 1999). They are primarily insectivorous, occasionally feeding on fruits or seed (Curson et al. 1994). Preliminary data indicate that natal dispersal ranges from less than 1 km to greater than 12 km (n=15) (Petit 1999). Although this is formulated from a small sample, it is on the same order as other song bird dispersal (e.g. Nice 1933). Here we posit warbler dispersal to be 5 km and return to the uncertainty of this statement later.
Geospatial Data
To my knowledge there are no current data on the spatial distribution of the focal species. The decline in trapping of the mink has perversely led to a decline in good biological information on the species in the study area. We am unaware of any work done with mink in the study area since Wilson's (1954) study. The Breeding Bird Survey indicates that this study area contains the one of the highest concentrations of prothonotary warblers in the Southeast (Price et al. 1995). Finer scale spatial information is not readily available.
Mink and prothonotary warblers are habitat specialists and use the same habitat. To identify habitat patches in the landscape we combined data from the US Fish and Wildlife Service's National Wetlands Inventory and land use coverage from the North Carolina Center for Geographic Information and Analysis (1996). Cells that were defined as being bottomland hardwood swamp or oak gum cypress swamp, and riverine, lacustrine, or palustrine forested wetland were selected as habitat patches. These cells were then aggregated into regions using an eight-neighbor rule and the intervening matrix described as nonhabitat. A unique patch identifier, area, core area, and edge area described patches. Core area was defined as the sum of the cells in the region with no neighbors touching nonhabitat cells using an eight-neighbor rule. Edge area equaled area minus core area. Cell size was 30 m. Transportation and hydrography digital line graphs were obtained for the study area from the US Geological Service at 1:100,000 resolution.
Because it is ecologically unlikely that either focal species will move from one patch to another in a straight line (e.g. over large expanses of open water or high development), we measured dij as centroid to centroid distances computed to reflect the navigability of the intervening non-habitat matrix. We used least-cost path algorithms on a simple cost-surface with a resolution of 90 m. Cost was defined by a surface comprised of x, y, and z where z was a uniform impedance that represented the cost of moving through that cell. Cells corresponding to areas of habitat were given a weight of 0.1, all other forest types were given a weight of 1. Cells classified as riverine/estuarine herbaceous were given a weight of 2. Shrubland was given a weight of 3. Sparsely vegetated cells (cultivated, managed herbaceous) were given a weight of 4. Areas of development and large water bodies were given a weight of 5. Streams were defined with a weight of 1 for both species while roads were assigned the high cost of 5 for mink only. We used global grid functions inside a macro in ArcInfo 7.2.1 (Esri 1998), with embedded Fortran programs, to iteratively loop through the array of patches and compute d ij for each unique pair of nodes in the array. Figure 3 shows a single patch centroid and the least-cost paths drawn to four other patch centroids.

  1. Figure 3. A source patch with least-cost paths drawn to four other patches.
Computing the lower triangle of the distance matrix using least-cost path algorithms is computationally challenging. For n patches there are n(n-1)/2 distances to compute. For this reason, and the fact that prothonotary warblers are not likely to persist in forest patches less than 100 ha (Petit 1999), we chose to only explicitly incorporate patches greater than 100 ha in my analyses. Using habitat patches greater than 100 ha results in 83 patches, roughly 82.5% of the 53,392 ha of possible habitat. Because all habitat patches, regardless of size, are given the lowest value on the cost surface, they are implicitly included in all analyses in that the species can traverse them easily as stepping stones accruing minimal cost. Although, this cost path analysis used slightly different costs for the two species the distance matrices for mink and warbler were virtually identical (r 2 = 0.99). Therefore, to simplify this analysis we chose to use the more conservative mink matrix for both species.

Terrestrial Example Results

The mean distance between patches in the 3403 x 3403 matrix is 62.7 km. The study area and habitat patches are illustrated in figure 2.
Edge Removal
The graph begins to disconnect and fragment into subgraphs at a 19 km edge distance and quickly fragments into numerous components containing only a few nodes (figure 3). The edges are drawn as straight lines between patch centroids with 5 km, 10 km, 15 km and 20 km thresholding distances in figure 4, even though the actual paths are computed by least-cost and are circuitous.
  1. Figure 4. Number of graph components and graph order as a function of effective edge distance.
  2. Figure 5. All graphs edges with increasing thresholded distances from 5 km to 20 km.

Terrestrial Network Discussion

Distance between patches can be measured in several different ways: edge to edge, centroid to centroid, centroid to edge, etc. However, measuring these as Euclidean distance makes little sense when the variance in mortality cost associated with traversal of the intervening habitats is large and cost associated with traversal of the intervening habitats is large and spatially heterogeneous. Few organisms or even ecosystem processes, such as groundwater movement or wildfire spread, move in this way. To differing extents they are all constrained by the landscape. Good multi-dimensional models exist to predict some ecosystem processes (e.g. pollution plumes, Bear and Verruit 1987) but not others. Spatially explicit models that simulate the movement of animals have been explored in some depth. However, most are complex parametric models which are data-hungry. They require specific information and are very hard to parameterize.
For this reason, we have computed D not as Euclidean distances but as a series of least-cost paths on a cost surface appropriate to the organism in question. These paths are designed to simulate the actual distance the focal species (or any other landscape agent, e.g. fire) covers from one patch to the next. For instance, in this riverine ecosystem the path a mink might take from one side of a river delta to the other would likely involve traversing the shore for 10 km under cover, rather than a 5 km swim across open water. This allows the animal to use stepping-stones of other habitat (low cost) along the way rather than set off into an unknown habitat matrix (high cost). The least-cost modeling combines habitat quality and Euclidean distance in determining d ij .
We explored alternative methods for constructing D including Euclidean distance and resistance-weighted distance between nodes. We found that the topology of the graph is robust and not sensitive to the difference between least-cost path distance and Euclidean distance except at the scale of large obstacles in the landscape. For instance, least-cost paths in this model did not cross the 5 km mouth of the Alligator river when moving from Eastern side of the study area but chose a route through habitat instead. In this case Euclidean distance and least-cost path distance were disparate. However, this particular landscape is rich in forest and focal species habitat and the differences in graph behavior between Euclidean distance D and a least-cost path distance D are small and ecologically appealing. There are so many low-cost cells in this particular landscape that resistance-weighted distance, whereby the least-cost path distance is expressed in cost-units, is largely meaningless.
The construction of these cost surfaces for the mink and prothonotary warbler was general enough so as to avoid committing the animals to movement patterns that are not readily possible to predict at 90 m resolution. The high correlation coefficient between the two matrices (r 2 = 0.99) indicates a paucity of difference between them and allows us to reduce data and present results more clearly. Cost-surface analysis has been only occasionally used by ecologists (e.g. Krist and Brown 1994, Walker and Craighead 1997) but widely used in computer science which is concerned with optimal route planning (e.g.McGeoch 1995, Bander and White 1998). This type of analysis is also common in applications of artificial intelligence (Xia et al. 1997).
Another benefit of a graph-theoretic approach is that dispersal biology does not need to be fully understood for the graphs to be interpretable. To give context to the graph framework we have postulated that dispersal for mink is 25 km and 5km for prothonotary warblers. Dispersal biology is incredibly complex and precise distances are usually unknown. Here these results and their interpretation are largely interpretable despite these uncertainness. Edge and node removal, as well as node sensitivity, are graph descriptors that are useful macroscopic metrics when dispersal can only be estimated (see Keitt et al. 1997 for an additional example).

Marine Connectivity Example

Analysis of the potential connectivity of patchy marine habitats has become an important topic in marine conservation. Better understanding the transport of planktonic larvae from known habitat sites to other suitable habitat sites is critical to proactive marine conservation planning. Generalized regional analyses have been conducted to identify the amount of "upstream' and "downstream" reef area and approximate travel time for large areas (Roberts 1997), but little work has been done on developing spatial analysis tools for assessing connectivity within a reef system.
In order to better assess this problem, a study area within the South Atlantic Bight was chosen and the locations of 178 artificial reefs, ship wrecks and rock bottom formations were digitized into a GIS system.

  1. Figure 6. Locations of artificial reefs, ship wrecks and bottom habita within the South Atlantic Bight study region off the coast of North Carolina.
Because currents contol the travel of free floating larvae, a model of current patterns depicting the direction and velocity must be employed. A physical oceanography model for the Mid-Atlantic and South Atlantic Bights was developed by Werner et al. (1999). A system of finite-element equations were used to create a three dimensional numerical model of ocean circulation in the region. The model is described in Lynch and Werner (1991), Lynch et al. (1996) and Werner et al. (1999). The model is a free-surface, tide-resolving model based on three dimsensional shallow water equations. The model mesh was constructed from irregularly spaced nodes spaced approximately 1km apart near shore and increasing in distance offshore. The original model covered the Mid-atlantic and South Atlantic Bight area using 3335 nodes. The study area used in this analysis represents a subset of 1100 nodes (see figure 7). Flow fields from the model representing the ocean circulation response to an easterly wind (90o direction) were used to represent an early spring (March-April) current regime for this analysis.

  1. Figure 7. Variable node distribution for the three-dimensional numerical physical oceanography model (After Cisco et al. 1999)
The oceanogrphic model provided North and East components of velocity that were added as attributes to a point coverage in a GIS. The 178 points representing the artificial reefs and wreck locations were added to provide to and from points to the model. A grid interpolation was performed to provide estimated North and East components of current velocity to the reef sites from their neighboring nodes. Vector calculations using the North and East components of velocity were performed to derive the direction in degrees of currents away from each node point. Using polar geometry, arcs were generated at standard distances travelling away from each node (see figure 8).

  1. Figure 8. A flow-field of 20km arcs created from North and East components of current velocity
In the next step, all arcs that travelled within a 500m fuzzy distance of more than one artificial reef were retained in the GIS model. This tolerance was set to allow flow paths that travelled near to reef habitats to be included in the network due to the fact that exact vector paths rarely encounter exact point locations. This tolerance simulates potential diffusion that may occur around site locations. Most of the artificial reefs and wreck locations represent small patches of habitat smaller than 500m in diameter. Out of 31684 possible reef connecting paths (number of reefs2 = 1782 = 31684) only 885 paths existed under the tested current regime (see figure 9).

  1. Figure 9. Connected reef paths under the easterly wind current regime (n = 884 paths).
Because marine ecologists are generally interested in the potential survivability of different current paths by marine life, a unit traverse time was calculated for each network segment based on the velocity and length of each segment. Figure 10 depicts the traverse time in hours for 20km long flow path segments in the study region. Notice the shorter traverse times around Cape Lookout and Cape Fear versus longer traverse times offshore.

  1. Figure 10. Segment traverse time in hours for 20km flow paths
In order to assess the pattern of connectivity for patchy marine environments, it is interesting to assess reef connectivity by different travel times. Many marine life forms must reach new suitable habitat within short periods of time to survive and establish. Viewing marine networks within temporal constraints allows marine ecologists to evaluate the possible importance of individual reef sites as they contribute to the connections between distant sites. Figure 11a-d depicts four different sets of reef connections differentiated by required traverse times. The number of individual conections between single reefs greater than 3 days were 767 out of 885 total connections. When we isolate long traverses, greater than 21 days, we have only 22 out of 885 paths.

  1. Figure 11. Connecting reef flow-paths by time: (a) Flow paths > 3 days (n = 767) ; (b) flow paths > 7 days (n = 526); (c) flow paths > 14 days (n = 140); (d) flow paths > 21 days (n = 22).
In the terrestrial example, we demonstrated that edge (path) removal could be used to test the response of the entire connected network to changes in the distance of traversability. If species can only travel a limited distance the components or number of separate units within the network increases. In the terrestrial example, travel distances were used to show how the network graph broke up into an increasing number of components as travel (edge) distances were decreased. In the marine example, travel time will used instead of distance for a similar result. Figure 12 depicts the relationship of increasing traverse time with decreasing number of separate graph components. A siginficant difference between the marine and the terrestrial examples is that the marine network never reaches a single, connected component at any travel time/distance. This is because in a directional, current driven environment, all patch connections ij will never be connected. Downstream patches will not transport materials to upstream patches. The reef connectivity example under the specified current regime shows a range of ~35 to 178 different components or potential management units depending on the amount of time species could tolerate traveling between patches.
  1. Figure 12. The change in the number of graph components or separated patches as a function of travel time in the marine reef example. Note: The number of components never reaches one as in the terrestrial example due to directionality in the current surface.
A final example of marine connectivity analysis involves path tracing. Using the path tracing function in ArcInfo, potential source paths can be calculated to a selected destination point. In the single example below (figure 13), A destination reef was selected and the path tracing functions identify all potential origin paths within the network model. This type of path tracing exercise would allow marine ecologists to calculate potential source locations for individuals or populations of species under specific current regimes.

  1. Figure 13. Marine path tracing example. A single destination reef was selected and all possible source paths flowing to that site are identified.
Other examples of path tracing in the marine environment include tracking the origins of marine mammal or turtle strandings based on the location of beaching and estimated time since death. This type of network analysis could prove to be useful and practical in marine management applications.

Conclusions

In this paper we contrast applications of network models applied to a terrestrial and a marine problem. In terrestrial situations, least-cost path algorithms can be used in an iterative manner to create a set of all potential paths between patches. In marine examples, directionality due to ocean currents requires a vector network approach to create the relative paths between habitat patches.
In addition to requiring two different types of path analysis approaches, marine and terrestrial analyses also differ in the development of the distance matrix. In the terrestrial example, species traveling between patches are expected to move equally well in either direction, so the distance between patches i and j are equivalent to j and i. This results in needing to calculate half of the possible distance matrix. In the marine example, the entire distance matrix must be calculated because the cost-distance from i to j may be different than j to i due to ocean current impedance. In addition, the terrestrial network can be assumed to be entirely connected into a single patch at some travel distance, while the marine network is never considered to be a single patch at any distance.
The methods contrasted here provide an objective framework for developing and assessing spatial connectivity for patchy environments in both nterrestrial and marine ecosystems.

Acknowledgements

We would like to acknowledge the contributions of Dean Urban and Tim Keitt for collaboration on the terrestrial network analysis and Cisco Werner for contribution of ocean current data for the marine example.

Literature Cited

Bander J.L. and C.C. White. 1998. A heuristic search algorithm for path determination with learning. IEEE Transactions On Systems Man and Cybernetics Part A-Systems and Humans, 28: 131-134.
Bear, J. and A Verruit. 1987. Modeling groundwater flow and pollution: with computer programs for sample cases. Kluwer Academic Publishers, Massachusetts.
Curson, J. Quinn, D. and D. Beadle. 1994. Warblers of the Americas: an identification guide. Houghton Mifflin, Massachusetts.
Dunstone, N. and J.D.S. Birks. 1987. Feeding ecology of mink (Mustela vison) in coastal habitat. Journal of Zoology 212: 69-84.
Esri (Environmental Services Research Incorporated). 1998. ArcInfo v. 7.2.1. Redlands, California

Gerell, R. 1970. Home ranges and movements of the mink Mustela vison in southern Sweden. Oikos, 21: 160-173.

Gross, J. and J. Yellen. 1999. Graph Theory and its applications. CRC Press, Boca Raton Florida.

Hamilton, W.J. 1940. The summer food of minks and raccoons on the Montezuma Marsh, New York. Journal of Wildlife Management 4: 80-84.

Hansen, A.J. and F. di Castri. 1992. Preface. Pages v-viii in A. Hansen and F di Castri (eds.), Landscape boundaries: consequences for ecological flows. Springer-Verlag. New York.

Hanski, I. 1998. Metapopulation dynamics. Nature, 396: 41-49.

Harary, F. 1969. Graph Theory. Addison-Wesley, Massachusetts.

Harris, L.D. The Fragmented Forest: Island Biogeography Theory and the Preservation of Biotic Diversity. University of Chicago Press, Illinois.

Harrison, S. 1994. Metapopulations and Conservation. Pages 111-128 in P.J. Edwards, N.R. Webb, and R.M. May (eds), Large-scale Ecology and conservation biology. Blackwell, Oxford.

Keitt, T.H., Urban D.L., and B.T. Milne. 1997. Detecting critical scales in fragmented landscapes. Conservation Ecology, 1: 4 (http://www.consecol.org/Journal/vol1/iss1/art4).

Krist F.J. and D.G. Brown. 1994. GIS modeling of paleo-indian period caribou migrations and viewsheds in northeastern lower Michigan. Photogrammetric Engineering And Remote Sensing, 60: 1129-1137.

Lynch, D.R. and F.E. Werner. 1991. Three dimensional hydrodynamics on finite elements. Part II: nonlinear time-stepping model. Int. J. For Num. Meths. Fluids 12:507-533.

Lynch, D.R., J.T.C. Ip, C.E. Naime, and F.E. Werner. 1996. Comprehensive coastal circulation model with application to the Gulfof Maine. Cont. Shelf Res. 16:875-906.

McGeoch C.C. 1995. All-pairs shortest paths and the essential subgraph. Algorithmica, 13: 426-441.

Miller B., Reading R., Strittholt J., Carroll C., Noss R, Soulé M., Sánchez O., Terborgh J., Brightsmith D., Cheeseman T., and D. Foreman. 1999. Using focal species in the design of nature reserve networks. Wild Earth, 8: 81-92.

Nice, M.M. 1933. Studies in the life history of the song sparrow (1964 reprint). Dover, New York.

Noss, R.F. 1991. Landscape connectivity: different functions and different scales. Pages 27-38 in W.E. Hudson (ed.), Landscape linkages and biodiversity. Island Press, Washington DC

Nowak, R.M. 1999. Walker's mammals of the world. 6th ed. Johns Hopkins University Press, Maryland.

O'Neil, R.V., B.T. Milne, M.G. Turner, and R.H. Gardner. 1988. Resource utilization scale and landscape pattern. Landscape Ecology, 2: 63-69.

Petit, L.J. 1999. Prothonotary Warbler (Protonotaria citrea) in The Birds of North America number 408 (A. Poole and F. Gill, eds). The Birds of North America Inc., Pennsylvania.

Price, J.P., Droege S., and A. Price. 1995. The summer atlas of North American Birds. Academic Press, New York.

Roberts, C.M. 1997. Connectivity and mangement of Caribbean Coral Reefs. Science 278:1454-1457.

Urban, D.L. and T. Keitt. Landscape connectedness: a graph-theoretic approach (in review).

Vankat, J.L. 1979. The natural vegetation of North America. John Wiley and Sons. New York.

Walker, W. and F.L. Craighead. 1997. Analyzing wildlife movement corridors in Montana using GIS. Proceedings of the 1997 Esri User Conference.

Werner, F.E., B.O. Blanton, J.A. Quinlan and R.A. Luettich. 1999. Physical oceanography of the North Carolina continental shelf during the fall and winter seasons: implications for the transport of larval menhaden. Fisheries Oceanography 8:7-21.

Wilson, K.A. 1954. Mink and otter as muskrat predators. Journal of Wildlife Management, 18: 199-207.

Xia Y., Iyengar S.S., and N.E. Brener. 1997. An event driven integration reasoning scheme for handling dynamic threats in an unstructured environment. Artificial Intelligence, 95: 169-186.

Appendix I

The following describes the macro PATH.AML and the process by which it computes the lower triangle of a distance matrix using area-weighted modeling in ArcInfo v. 7.2.1 (Esri 1998).
Least-cost path modeling
PATH uses area weighted distance functions in GRID to calculate least-cost paths. These functions are similar to Euclidean distance functions (e.g. EUCDIST) but instead of working in geographical units they work in cost units. Geographical units are constant across space while cost units are not necessarily linked to space. Cost is defined by a surface comprised of x, y, and z where z can be any uniform impedance that represents the cost of moving through that cell (travel time, dollars, preference, etc.).
The first step in weighted distance modeling is to determine the accumulated cost of moving away from a source cell. The algorithm for the COSTDISTANCE function uses graph theory with the centers of cells acting as nodes and the distance between them as edges. There are eight edges connecting a node to adjacent nodes. The edge distance is in cost units and defined by the z value of the nodes themselves. The edge distance in a cardinal direction from one cell to another is simply calculated as  . For diagonal movement the algorithm is similar but includes the square root of 2, . Defining the values for the COSTDISTANCE output grid is a graph operation where by each cell is assigned a cost to return to the source cell. This operation is designed so that each cell's value is the sum of the edge values for the shortest possible walk connecting it to the source node. The process is iterative and anneals from the source cell. The output grid identifies the accumulative cost for each cell to return to the source cell. However, it does not show how to get there. The cost back-link returns a grid coded for the eight neighbors that can be used to reconstruct the route to the source. The COSTPATH function determines an output grid that spatially defines the least-cost path from a destination cell to the source cell. This is done using the back-link grids to retrace the destination cell to the source.
Logic
In ArcInfo, least-cost paths are most easily computed by interactively choosing the "from" and "to" point on the screen with the mouse and cross-hair cursor. Although useful for some least-cost applications, this is plainly impractical to compute for more than a few paths. The following macro requires a coordinate datafile of the nodes formatted in id , x , y . This file can be produced using the UNGENERATE command (with the last line `end' removed) in ArcInfo. PATH produces least-cost paths for every pair of nodes in the array. The macro is initialized using a simple front end. The user is prompted for the coordinate data file, the cost surface, and a name for the output file. The number of cost paths to compute is written to the screen and the user is asked if they wish to keep them as line coverages.
The program begins by running the Fortran program MATRIX as system call from Arc. The front end is distributed between AML and Fortran. When MATRIX is run the user is prompted for the file containing the coordinate data. That file is reformatted into a column vector of x and y coordinates so that the input file:

is reformatted in columnar format as x i y i x j y j , from i = 2 to n and j = 1 to n-1. This is the lower triangular matrix without the diagonal:

The number of iterations the model will run in ArcInfo (n(n-1)/2) is written to the screen. The spatial array (col_space.dat) and the matrix notation (ij.dat) that goes along with it are written before MATRIX ends.

When MATRIX finishes PATH continues to initialize. The spatial array (col_space.dat) is reformatted slightly using a one line Awk script to a new file, column.dat. In Fortran numeric data is right justified and ArcInfo needs numeric data to be justified left before reading. If data is not justified left it is read in as a character string in Arc.

From this point on PATH begins a conditional loop that begins with the number of the count variable. The `do while' structure is convenient here so that the loop will continue to run until the end of the spatial array (column.dat) is reached. A workspace is created called `cost_space' in which all of the Arc commands are run. To facilitate data management and keep memory usage to a reasonable level this workspace is initialized each iteration. The coordinate data is read into Arc as variables for the source and destination grids using SELECTPOINT. Once again, the vagaries of ArcInfo require an Awk script to manage the data array. There are four variables to read in with each iteration the x and y coordinate for the "from" and "to" points on the cost grid. An easy was to have Arc read these is to delete each line in the array after the variable is written using Awk. By doing this the spatial array is consumed by the program while it runs.

After the source and destination grids are defined as integers the COSTDISTANCE and COSTPATH commands are executed using the destination and source grids respectively. The distance grid is calculated as distance from the source using the user-defined cost surface. The distance grid is used in conjunction with the source point to compute the cost path as a grid. The cost path is converted to a vector using the GRID function GRIDLINE. If the user specified that the line coverages were to be kept then the line file is saved as `path_1', `path_2',...'path_n.' The distance of the path is written to a text file using TABLES. Then the loop begins anew with fresh coordinates.

The loop ends when the entire column of coordinates has been used Then the distance data is pasted together with the appropriate matrix notation (ij.dat) as i, j, d ij using the Fortran program OUTPUT and saved as the user defined outfile. The temporary files are cleaned before PATH ends.