Using landscape metrics to characterize towns along an urban-rural gradient

Urban-rural gradients are useful tools when examining the influence of human disturbances on ecological, social and coupled systems, yet the most commonly used gradient definitions are based on single broad measures such as housing density or percent forest cover that fail to capture landscape patterns important for conservation. We present an approach to defining urban–rural gradients that integrates multiple landscape pattern metrics related to ecosystem processes important for natural resources and wildlife sustainability. We develop a set of land cover composition and configuration metrics and then use them as inputs to a cluster analysis process that, in addition to grouping towns with similar attributes, identifies exemplar towns for each group. We compare the outcome of the cluster-based urban-rural gradient typology to outcomes for four commonly-used rule-based typologies and discuss implications for resource management and conservation. The resulting cluster-based typology defines five town types (urban, suburban, exurban, rural, and agricultural) and notably identifies a bifurcation along the gradient distinguishing among rural forested and agricultural towns. Landscape patterns (e.g., core and islet forests) influence where individual towns fall along the gradient. Designations of town type differ substantially among the five different typologies, particularly along the middle of the gradient. Understanding where a town occurs along the urban-rural gradient could aid local decision-makers in prioritizing and balancing between development and conservation scenarios. Variations in outcomes among the different urban-rural gradient typologies raise concerns that broad-measure classifications do not adequately account for important landscape patterns. We suggest future urban-rural gradient studies utilize more robust classification approaches.

Objectives We present an approach to defining urbanrural gradients that integrates multiple landscape pattern metrics related to ecosystem processes important for natural resources and wildlife sustainability. Methods We develop a set of land cover composition and configuration metrics and then use them as inputs to a cluster analysis process that, in addition to grouping towns with similar attributes, identifies exemplar towns for each group. We compare the outcome of the clusterbased urban-rural gradient typology to outcomes for four commonly-used rule-based typologies and discuss implications for resource management and conservation. Results The resulting cluster-based typology defines five town types (urban, suburban, exurban, rural, and agricultural) and notably identifies a bifurcation along the gradient distinguishing among rural forested and agricultural towns. Landscape patterns (e.g., core and islet forests) influence where individual towns fall along the gradient. Designations of town type differ substantially among the five different typologies, particularly along the middle of the gradient. Conclusions Understanding where a town occurs along the urban-rural gradient could aid local decision-makers in prioritizing and balancing between development and conservation scenarios. Variations in outcomes among the different urban-rural gradient typologies raise concerns that broad-measure classifications do not adequately account for important landscape patterns. We suggest future urban-rural gradient studies utilize more robust classification approaches.

Introduction
Urban-rural gradients are useful tools when examining the influence of human disturbances on ecological (Nagy and Lockaby 2011), social (Timm et al. 2015), and coupled natural-human systems (Liu et al. 2007;Ostrom 2009). Urban areas typically have dense populations, few natural lands, and abundant impervious surface, whereas rural areas contain low population density, much natural or cultivated land, and little impervious surface. Often rural areas are defined as those that are 'not urban' (Ratcliffe et al. 2016). However, social and ecological scientists have long acknowledged the merit of looking beyond the simple urban-rural dichotomy, recognizing the advantages of more nuanced characterizations of urbanization (McDonnell and Pickett 1990;Isserman 2005). For example, research sites are frequently selected along an urban-rural gradient to assess the effects of urbanization on ecological processes or outcomes (Ahrné et al. 2009), and a recent review of gradient studies found wildlife species' responses to urbanization to vary among: (i) positive, (ii) negative, (iii) intermediate, (iv) punctuated, (v) bimodal, and (vi) no response (McDonnell and Hahs 2008). This complex response of wildlife to urbanization highlights several important lessons for conservation policy-makers and practitioners. For one, a suite of conservation efforts may be needed to achieve widespread sustainable outcomes across the gradient (Norton et al. 2016). More fundamentally, reflecting the way that wildlife experience and use the landscapes along the urbanrural gradient is important for constructing effective land management and conservation policies (DeStefano and DeGraff 2003; Thornhill et al. 2017;Xun et al. 2017). A simple urban-rural dichotomy does not fully capture the landscape patterns influencing ecosystem processes and wildlife responses.
However, finding an appropriately nuanced gradient definition for a specific research purpose can be challenging. A standard system for characterizing urban-rural gradients would allow for more consistent public policy across regions and clarify comparisons of outcomes among research studies (MacGregor-Fors 2011); yet, differences in geographic scale, climate, ecosystem services, culture, and study objectives (i.e., the intended application of the gradient) make universal application of a single definition unlikely (McDonnell and Hahs 2008). Instead, McDonnell and Hahs (2008) describe three categories of measures used to define urban-rural gradients: (1) transects, consisting of one or more straight lines between an urban site and a rural site, with distance as the sole underlying metric; (2) broad measures (e.g., human population density or percent forest cover), which attempt to encapsulate the urbanization process into a single metric; and, (3) specific measures, which attempt to elucidate a more direct relationship between urbanization processes and outcomes of interest, and that may involve complex landscape pattern metrics (e.g., patch size and connectivity) or the combination of multiple metrics.
These gradient measures can be implemented in continuous or categorical form. Continuous metrics such as percent urban land, percent impervious surface, percent forest cover, or distance from city center are commonly used to select and characterize study sites along an urban-rural gradient (Supplementary Materials Table S1). In comparison, a number of studies have used categorical measures to define urban-rural gradients with typically three or four classes (e.g., urban, suburban, exurban, and rural), although as many as eight or nine (Table S2). Two approaches, segmentation rules and cluster analysis, have been used to generate categorical urban-rural gradient typologies. Rule-based typologies, in which the researcher chooses specific metric thresholds to define each class, are most common (Table S2). Often based on broad measures such as population or housing density, rule-based typologies are easily reproducible and allow for direct comparisons of study outcomes over space and time (McDonnell and Hahs 2008;Padilla and Sutherland 2019). However, applying specific metric thresholds unilaterally can mask regional diversity. In addition, beyond two or three metrics, rules can become prohibitively complex or generate numerous classes too difficult for stakeholders to differentiate or interpret. As an alternative, a small number of recent studies have used cluster analysis to generate urban-rural gradient typologies (e.g., Owen et al. 2006;. Cluster algorithms readily accommodate multiple complex metrics and algorithmically identify class thresholds. Unlike rule-based approaches, typologies generated using cluster algorithms adjust to the set of underlying data (e.g., towns in New England versus towns across the continental United States), potentially offering flexible solutions that adapt to the specific management or policy application.
The spatial unit of analysis also varies among urban-rural gradient studies (Tables S1 and S2). Some studies use map pixels or well-defined buffers around study sites to divide their study region into components of the same size and shape, fostering easy comparisons among units and smooth transitions over a region. The key decision is what size pixel or buffer to use. Other studies use administrative units (e.g., census blocks, municipalities, or counties) or environmental units (e.g., watersheds) because they better facilitate jurisdictional resource management or conservation policy analysis and allow for easier interpretation by the general public. Pixels rarely match administrative or environmental units and may require complex aggregation (if pixels are too small) or resampling (if pixels are too large). When defining an urban-rural gradient, it is important to match both the type of metrics and the spatial unit of analysis to the objectives of the study (McDonnell and Hahs 2008).
In a recent review of 250 urban-rural gradient studies, Padilla and Sutherland (2019) revealed that studies that define gradients using multiple metrics were more likely to see a significant relationship between urbanization and the ecological processes being studied. While many studies define and use urban-rural gradients, the majority use a rules-based approach typically based on 1-3 broad measures (e.g., housing density). Only a few studies define the gradient using a clustering process (Table S2). Incorporating large numbers of metrics using a rules-based approach is onerous, at best. Alternatively, some researchers may find the use of cluster analysis daunting. Our goal is to show that this latter concern need not be the case.
In this paper, we present a cluster-based approach to defining urban-rural gradients that integrates multiple landscape pattern metrics that drive ecosystem processes important for sustainability of local natural resources and wildlife. We use the resulting typology to characterize administrative units along the gradient. The use of several landscape metrics means that our approach captures more landscape complexity and better reflects the unique combination of attributes found in focal landscapes than simple rule-based classifications. The use of administrative spatial units facilitates analysis and decision making in resource management and conservation. While our urban-rural gradient approach can be applied to many landscapes and natural resource contexts, here we illustrate the process for small local resources in New England, United States (U.S.). The resulting typology will ultimately be used to facilitate the assessment of policies and programs targeting management and conservation of small natural features (e.g., seasonal pools and rocky outcrops) in urbanizing landscapes (Hunter et al. 2017). We define our urban-rural gradient by first developing a set of land cover composition and configuration metrics that typify New England landscapes. Unlike most previous studies (Table S2), we intentionally avoid mixing socio-economic and ecological metrics in our classification because it can inhibit the ability to examine the links between drivers of anthropocentric disturbance and ecological response outcomes (Cadenasso et al. 2007). We then use these metrics as inputs to a clustering process that produces an urban-rural gradient typology consisting of five town types: urban, suburban, exurban, rural, and agricultural. We compare our resulting cluster-based typology to four commonly-used, rule-based typologies using broad measures (population density, housing density, proportion impervious surface, and percent forest cover) and discuss implications for natural resource management and conservation.

Study area
Our study area is the six-state New England region in northeastern United States ( Fig. 1; states include Connecticut, Maine, Massachusetts, New Hampshire, Rhode Island, and Vermont). The region varies considerably in land cover and demographics. Forest cover ranges from 10% in densely populated eastern Massachusetts to 80% in northern New England (Ducey et al. 2016). Wetland area encompasses nearly 25% of Maine, however, only approximately 5% of Connecticut and Vermont (Fretwell et al. 1996). Likewise, impervious surfaces account for nearly 11% of Rhode Island and only 2% of Vermont (Nowak and Greenfield 2012). Demographically, population density varies from 0.0014 to 7123 persons per squarekilometer, and housing density varies from 0.015 to 3171 units per square-kilometer, with greater densities occurring in major metro areas, primarily in southern New England, and lesser densities occurring in northern New England. From 2010 to 2017, 78% of New England counties saw housing growth above the national median (Foster 2017; Census Bureau 2018).
We conduct our analysis using town-level administrative units because in New England, local town governments, in the form of democratic town meetings or town officials (e.g., mayors, town councils, planning boards, and conservation commissions), are responsible for making decisions about land use regulations and infrastructure investments. Although state regulations apply to all towns, local regulations, in combination with various social and ecological processes, can create distinctive land cover composition and configuration across towns. Towns, as referred to here, encompass minor civil divisions as defined by the U.S. Census and include: cities, towns, townships, gores (irregular, non-surveyed parcels), grants, Federally-recognized tribal lands, locations, plantations, purchases, survey townships, and unincorporated territories. We avoided computational complications from outliers by excluding 16 very small towns (\ 5.0 km 2 ) from our analysis, leaving 1590 study towns.

Landscape metric definition
We characterized our urban-rural gradient with seven land cover composition and configuration metrics (Table 1 A) that vary spatially across our study area (Fig. S2). We generated these metrics by first creating an aggregated land cover map from the 2011 U.S. Geological Survey National Land Cover Database (NLCD) 30-meter resolution land cover and percent developed imperviousness data layers (Xian et al. 2011;Homer et al. 2015). We included multiple land cover classes to capture the mix of types in our study area that are relevant for conservation. Our aggregated NLCD land cover classes include: (1) forested lands: deciduous forest, evergreen forest, mixed forest, woody wetland, and emergent herbaceous wetlands, (2) developed areas: lands where imperviousness is [ 20%, (3) agricultural lands: pasture/hay and cultivated crops, (4) open water, and (5) other land covers (Fig. S1). We used forest and wetlands as proxies for natural lands, understanding that not all forested lands are natural. We chose a threshold of 20% imperviousness to capture areas of residential, commercial and industrial development, because imperviousness [ 20% has been shown to have measurable effects on surrounding ecosystems (Arnold Jr and Gibbons 1996;Brabec 2009). Agricultural lands also affect wildlife, water quality, and other ecosystem services (Christin et al. 2013;Chaplin-Kramer et al. 2015;Lee and Carroll 2015), although not necessarily in the same manner as impervious surface, so we created a separate agricultural land cover class.
We calculated three land cover composition and four landscape pattern metrics for each our 1590 towns (Table 1 A). Our land cover composition metrics are calculated as the proportion of town land area that falls into one of three aggregated land cover classes: proportion forested lands, proportion agricultural lands, and proportion developed areas. We created our land cover configuration (i.e., pattern) metrics using Morphological Spatial Pattern Analysis (MSPA), an image analysis technique that assigns pixels to a pattern type based on their spatial relationship with other pixels of the same cover type (Vogt et al. 2006;Soille and Vogt 2009). MSPA has been applied in landscape spatial analyses to rank riparian corridors (Burton et al. 2005), evaluate green infrastructure networks (Wickham et al. 2010), illustrate beetle outbreak patterns (Chen 2014), and characterize forest fragmentation (Rogan et al. 2016). We used the open-source software GuidosToolbox (Vogt and Ritters 2017) that uses a single size threshold parameter to segment binary patterns into seven mutually-exclusive categories: core, islet, bridge, loop, edge, perforation, and branch (see Soille and Vogt 2009 for details). After investigating different values for the size parameter and all seven pattern types across our study region, we selected a size parameter of one pixel (30-m) and core and islet patterns for forested lands and developed areas to represent key spatial patterns in our towns. The toolbox designates contiguous land cover that includes a complete edge and additional pixels as core; if there are no additional pixels, the patch is designated islet. An abundance of developed core suggests densely built areas such as commercial developments or industrial parks, while an abundance of developed islets suggests isolated buildings and sparse development scattered across the landscape. Large swaths of core forested lands may provide relatively undisturbed natural resources and wildlife habitat, while predominantly islet forested lands represent more fragmented habitat or urban parks. These landscape pattern metrics are representative of an urban-rural gradient that reflects human disturbances (e.g., forest fragmentation) and wildlife responses (e.g., altered biodiversity; Andren 1994;Galster et al. 2001;Cushman 2006;Hanski 2015). Similar to our composition metrics, we calculate our land cover pattern metrics as the proportion of town land area.

Principal component analysis
Several of our landscape metrics were highly correlated (Table S3). In particular, forest and core forest were highly positively correlated, as were developed and core developed (Pearson correlation coefficients [ 0.9). Developed was moderately negatively correlated with forest and core forest, as was forest with core developed (Pearson correlation coefficients between -0.7 and -0.9). Notably, islet forest was modestly positively correlated with developed (Pearson correlation coefficient: 0.69). Unsurprisingly, agriculture was not correlated with any of the other landscape metrics (Pearson correlation coefficients all \ |0.25|). Interestingly, islet developed was also not correlated with any of the other landscape metrics (Pearson correlation coefficients all \ |0.44|). Developable land is all land in the town that is: (1) not part of the protected area database AND (2) not water Rather than arbitrarily eliminate one or two of the most correlated metrics from further analysis, we reduced dimensionality in our dataset with a principal component analysis (PCA). PCA is a multivariate statistical technique that extracts the most important information from a set of explanatory variables into a smaller set of principal components that explain the majority of the variation in the original variables (Abdi and Williams 2010). The principal components are ordered in terms of explanatory power; the first principal component, PC1, explains the highest amount of variation, the second component, PC2, the second highest amount of variation, and so on. The principal components can be interpreted directly or used as replacements for the original set of explanatory variables in subsequent analyses. PCA can also be used to represent the pattern of similarity among the original observations by displaying them as points in ''maps'' typically in the form of a scatter plot of the first two principal components. PCA is used in many fields of study including ecology and land use planning (Hahs and McDonnell 2006;Owen et al. 2006;. We used principal components as inputs to the town-level clustering process described below.
Cluster-based urban-rural gradient typology We used an affinity propagation clustering technique (Frey and Dueck 2007) to generate our town-level urban-rural gradient typology. Previously applied to landscape metric datasets (Cardille and Lambois 2010;Cardille et al. 2012;Partington and Cardille 2013), affinity propagation has been shown to identify clusters faster and with less error and sensitivity to outliers than other clustering methods (e.g., k-centers clustering). Affinity propagation clustering also identifies representative exemplars, which we use to illustrate our urban-rural gradient. We evaluated a range of solutions from two clusters through ten clusters, the range of distinct classes used in previous urban-rural gradient studies. We guided our decision in the final number of clusters by: (1) visually inspecting scatter plots of the first two principal components with individual town data points colorcoded by assigned cluster, (2) examining frequency histograms and box plots of the underlying landscape metrics by cluster, and (3) testing for significant differences between typology classes using one-way ANOVA tests for each of the landscape metrics and performing pair-wise comparisons of class means using Tukey HSD tests.

Rule-based urban-rural gradient typologies
We compared our cluster-based typology to four alternative rule-based typologies (Table 1B). First, following Isserman (2005), we defined a typology based on population density and ''urban areas'' and ''rural areas'' as designated by the U.S. Census Bureau at the town level (Ratcliffe et al. 2016): (1) urban: population density is [ 500 people per square mile AND 90% of the population is urban; (2) mixed urban: not rural AND not urban AND population density is [ 320 people per square mile; (3) mixed rural: not rural AND not urban AND population density is\ 320 people per square mile; and (4) rural: population density \ 500 people per square mile AND 90% of the population is rural.
Second, following Theobald (2005), we defined a typology based on lot area per housing unit (i.e., housing density): (1) urban: \0.10 hectares; (2) suburban: 0.10-0.68 hectares; (3) exurban: 0.68-16.18 hectares; and (4) rural: [16.18 hectares. Housing density is calculated within 'developable land,' which is defined as land that is not water and not part of the USGS Protected Areas Database (https:// www.usgs.gov/core-science-systems/scienceanalytics-and-synthesis/gap/science/protected-areas). Theobald (2005) presents two different thresholds for distinguishing between rural and exurban. We chose the recommended threshold of 16.18 hectares for the lower bound on rural home parcels. However, because crop and timber production can be productive on small land holdings in New England, we also assessed the smaller threshold of 8.09 hectares.

Principal component analysis
Results of the PCA (Table 2) show the importance of each principal component (i.e., the proportion of variation among the PCA input variables explained by that component) and the loading scores for each input variable, which indicate the strength of the relationships between the input variables and the individual components. We consider loading scores greater than 0.4 (or less than -0.4) to indicate a strong relationship. The first principal component (PC1) explains 61% of the variation among the original landscape metrics. Forest and core forest are positively associated, while developed and islet forest are negatively associated with PC1. As such, PC1 captures the common urban-rural gradient dichotomy with forest lands at one end and developed areas at the other, with the additional contribution of small urban parks (i.e., islet forests). In comparison, core developed is positively associated, while islet developed and agricultural lands are negatively associated with PC2. That is, PC2, which explains an additional 21% of the variation among the landscape metrics, captures an alternative representation of the urban-rural gradient dichotomy with developed areas at one end and agricultural lands at the other, with the additional contribution of small islets of developed lands within the farming landscape. Recalling that both agriculture and islet developed were the only landscape metrics not correlated with any others, PC3 indicates that these two metrics are diametrically opposed along a third component. We elected to use the first three principal components from our PCA as inputs to the cluster analysis process, because they collectively explain 94% of the variation among our landscape metrics (Table 2).
Cluster-based urban-rural gradient typology We chose the five-cluster solution to represent the urbanrural gradient (Fig. 1), because it captures key land use disturbances (development and agriculture) in our study area, distinguishes among land cover composition and patterns relevant for management and conservation of small natural features, and can be meaningfully compared to other urban-rural gradient typologies. Visual inspection of a scatter plot of the first two principal components with data points color-coded by town type (Fig. 2), as well as frequency histograms and box plots of the underlying landscape metrics (Fig. 3), indicated distinct clusters of towns. Statistically, means for all seven underlying landscape metrics are significantly different among the town types, and pairwise comparisons of means of each metric differed significantly between nearly all town types (Tables S4, S5).
We named our five town types in accordance with commonly used urban-rural gradient labels (Table S2), informed by values of the seven underlying landscape metrics (Fig. 3). Urban towns (6% of towns in our study area) are concentrated in southern New England and have the greatest proportion of developed areas and core developed areas and the smallest proportion of forested lands of all town types. Suburban towns (14%), primarily fringing urban areas and interstate highways, are less developed than urban towns, with smaller proportions of developed areas and larger proportions of forested lands. Suburban towns also appear to be more fragmented than exurban or rural towns, having more of both islet developed areas and islet forested lands. Exurban towns (31%) neighbor suburban areas, border other main roadways, and are less developed than suburban towns, having on average smaller proportions of developed areas, both core and islet, and greater proportions of forested lands and core forested lands. Exurban towns also contain more agricultural lands than either suburban or rural towns. Rural towns (46%), characterized by large proportions of forested lands and core forested lands, small proportions of developed areas, and little to no agricultural lands, dominate northern New England, but are also abundant in western Massachusetts.
Agricultural towns (3%) are uniquely characterized by abundant agricultural lands, amounts of developed areas similar to exurban towns, and amounts of forested lands similar to suburban towns. There are two main regions of agricultural towns, one in northeastern Maine and a second in northwestern Vermont, and smaller agricultural regions in southeastern Massachusetts and the Connecticut River valley. The four land cover pattern metrics vary as expected along our urban-rural gradient (Fig. 3). The proportion of core developed areas decreases, while the proportion of core forested lands increases along the urban-suburban-exurban-rural gradient (Fig. 3c,  d). Further, the proportion of islet developed areas is greater in suburban, exurban, and agricultural towns and less in urban and rural towns (Fig. 3e). In urban towns, where impervious surface covers comparatively more of the land area, isolated development is not as common. Counterintuitively, rural towns also have small proportions of developed islets because they have very few developed areas of any type. It is perhaps more intuitive to note that the proportion of developed areas that are islets (rather than the proportion of the entire town land area) is greatest in Fig. 2 Results of fivecluster solution presented as a scatter plot of the first two principal components used in the cluster analysis. PC1 is driven by forest and core forest landscape metrics in the positive direction and developed, core developed, and islet forest landscape metrics in the negative direction. PC2 is driven by the core developed landscape metric in the positive direction and islet developed and agricultural landscape metrics in the negative direction. Loading scores are provided in Table 2 rural towns, followed by agricultural and exurban towns (Fig. S3). While the proportion of islet forested lands is relatively small across all town types, the proportion tends to decrease from urban to rural towns along the gradient (Fig. 3f). That is, natural land cover in urban areas is more likely to be provided by small, isolated parks than in other town types.
The set of exemplar towns distinctly illustrate the urban-rural gradient across our study area ( Fig. 4; Table S6). The urban exemplar is dominated by core developed areas, while the suburban exemplar features relatively equal amounts of developed areas and forested lands in a highly fragmented pattern. The exurban exemplar is primarily forested lands but also features a small core developed area and several small patches of agricultural lands. The rural exemplar is dominated by forested lands, while developed areas and a few patches of agricultural lands border roads. The agricultural exemplar is dominated by agricultural land and has a much smaller amount of forested lands than both the rural and exurban exemplars. Although the agricultural exemplar has a large amount of water, recall that water is not included in our metric calculations. Note also that the ''donut hole'' in the agricultural exemplar town is a feature of Vermont town governance, whereby highly developed village centers are considered separate towns-in this case, the 'missing piece' of St. Albans Town is a small urban town called St. Albans City (Fig. 1, agricultural  exemplar).

Alternative rule-based urban-rural gradient typologies
Applications of the four rule-based urban-rural gradient typologies to the New England region reveal interesting comparisons among themselves and with our multi-metric, cluster-based typology (Figs. 1, 5  and 6). The population-based typology (Figs. 5a and 6) is characterized by the greatest percentage of urban towns compared to the other typologies and very few towns designated as either mixed urban or mixed rural. Towns designated as urban or rural in the populationbased typology generally correspond to towns labeled urban or rural in our cluster-based typology. However, because there are many more urban and rural towns generated by the population-based typology, many population-based urban towns are designated as suburban and many population-based rural towns are designated exurban in our cluster-based typology. Our cluster-based agricultural towns are most often designated as rural in the population-based typology (Figs. 1 and 5a).
The housing-based typology (Fig. 5b) is dominated by exurban towns, which account for 67% of all towns in the study area, compared with 31% of towns in the cluster-based typology and even less in the other rulebased typologies (Fig. 6). Exurban towns in the housing-based typology fall relatively equally into the rural or exurban category in our cluster-based typology. Beyond that, there is substantial overlap between rural and urban town types in the two typologies and, to a lesser extent, there is also overlap in suburban town types. Surprisingly, two towns designated as urban in the housing-based typology are designated rural in our cluster-based typology; however, this is purely a function of the use of developable area rather than total land area in the calculation of housing density. These two towns are small and nearly all of the land area is protected, which resulted in a very small amount of developable area per housing unit. Our cluster-based agricultural towns are most often designated as exurban in the housingbased typology (Figs. 1 and 5b).
The impervious surface-based typology (Fig. 5c) is dominated by rural towns (76%), illustrating the relatively sparse impervious surface across the study area. The remaining three town types are concentrated along major roadways and the major metro areas of southern New England. Our cluster-based agricultural towns are all designated as rural in the impervious surface-based typology (Figs. 1 and 5c).
The forest-based typology (Fig. 5d) is also dominated by rural towns (63%), but with a substantial portion of exurban towns (26%) particularly in the southern region of the study area, illustrating the broad distribution of forested lands. Only 1% of towns is classed as urban, coincident with the region's major cities. Our cluster-based agricultural towns are most often designated as suburban in the forest-based typology (Figs. 1 and 5d). b Fig. 3 Land cover composition and pattern metrics by town type for cluster-based urban-rural gradient typology. Boxes show middle quartiles, lines show outer quartiles, Xs indicate medians, and dots show outliers

Discussion
Our cluster-based urban-rural gradient typology discriminates among five town types: urban, suburban, exurban, rural, and agricultural (Fig. 1). Distinguishing among developed, forest, and agricultural land cover composition and, more importantly, configuration (i.e., pattern) was essential for achieving our research objective of defining a typology that specifically relates human disturbances to ecosystem processes (e.g., nutrient flows, habitat provisioning) important for sustainability of local natural resources and wildlife in order to facilitate future land management and conservation policy analysis. Below, we discuss the implications of incorporating pattern metrics and agricultural land covers in our methodological approach and the challenges of choosing among available approaches.

Pattern metrics and urban-rural gradient typologies
In general, all five urban-rural gradient typologies we examined illustrate an urban-rural gradient with urban towns surrounded by suburban towns, then exurban  towns, and rural towns filling out the remainder of the region (Figs. 1 and 5). However, the typologies differ in two important ways. First, the broad measures used in the rule-based typologies do not capture variation in landscape patterns that are explicitly incorporated into the cluster-based typology. The effect of including core and islet pattern metrics is clearly illustrated (Fig. 7) with major overlaps of town types along percent forest and percent developed dimensions. The incorporation of islet developed areas and islet forest lands, in particular, is instrumental in distinguishing among urban and suburban, suburban and exurban, and exurban and rural town types in the cluster-based typology (Figs. 2, 3 and S3). Similar substantial overlaps of cluster-based town types exist along population density and housing density dimensions, with substantial overlap among rural and exurban towns at the low-density end of the spectrum and among urban and suburban at the high-density end. Note, for example, that the housing-based typology designates the cities of Manchester (NH), Portland (ME), Springfield (MA), and Worcester (MA) as suburban rather than urban because it fails to capture non-residential development (Fig. 5b). These four towns are characterized by large proportions of core developed areas (including commercial, industrial, and transport land uses in addition to residential uses) and small proportions of core forested lands. Thus, the rule-based typologies based on broad measures are not capturing landscape patterns that are important considerations for resource management and conservation Bode and Maciejewski 2014;Crooks et al. 2017;Rosetti et al. 2017).
Second, the quantity of towns falling within each class is dramatically different (Fig. 6). Excluding the relatively small class of agricultural towns, the clusterbased (land cover pattern) typology results in a relatively linear increase in numbers of towns in each town type as might be expected from common conceptual models illustrated by linear transects (McDonnell and Pickett 1990) or concentric circles of gradually increasing size emanating from the urban core (von Thünen 1966). In contrast, each rule-based typology results in a relatively skewed distribution with a single town type dominating the landscape, which may identify too few or too many towns for subsequent policy assessment, particularly if the research objectives are targeted at the middle of the urban-rural gradient, for example because forest  (Table 1A). Population-based typology is based on population density rules following Isserman (2005). Housing-based typology is based on housing density rules following Theobald (2005). Impervious-based typology is based on percent developed rules loosely following Ahrné et al. (2009). Forest-based typology is based on percent forested lands using simple rules. Specific rules are described in the main text and Table 1B. Agricultural towns are only applicable in the cluster-based (land cover pattern) typology fragmentation is higher (Irwin and Bockstael 2007). In our cluster-based typology, suburban and exurban towns comprise 45% of towns in the study area, while the population and impervious surface-based typologies identify less than 25% of towns, and the housingbased typology identifies more than 85% of towns in the middle of the gradient. It may be argued that the thresholds for each of the rule-based typologies could be set such that the distribution of towns among town types resembles that of the cluster-based typology. In fact, the high correlation among the four broad measures (0.85 to 0.99) would easily enable close alignment of town designations among the four rulebased typologies. However, we remained faithful to the original authors' thresholds because they chose them deliberately (Isserman 2005;Theobald 2005, Ahrné et al. 2009). Furthermore, despite any reassignment of thresholds, the important observation above remains valid. That is, because of the substantial overlap of the cluster-based town types along all four broad measures, there is no alternative rule-based threshold selection that would align any of the rulebased town types with those of the cluster-based typology, indicating that the rule-based typologies do not capture the landscape patterns found across the New England region.
Rural towns, agricultural towns, and land use intensity While agricultural lands are used to qualitatively describe areas along urban-rural gradients (e.g., Saito and Koike 2013;Wadduwage et al. 2017), they are rarely explicitly incorporated into quantitative gradient definitions (although, see Nagy and Lockaby 2011, Lee and Carroll 2015. In most classification systems, towns with abundant agricultural lands typically fall into either rural or exurban categories (Table S2). Depending on the broad-measure used, the rule-based typologies bounds around each town type, illustrating the overlap along these axes due to the inclusion of pattern metrics and agricultural land covers in the cluster analysis. Exemplar towns for each town type are labelled and identified by a black dot designated agricultural lands in our study area as rural, exurban, or suburban (Figs. 1 and 5). Our clusterbased urban-rural gradient typology placed agriculturally dominated towns in a class of their own rather than embedding them within another class, a result of our including agricultural land cover types as inputs to the classification process. Towns in our New England study area that contain large proportions of agricultural lands were consistently grouped together by the clustering algorithm, regardless of the number of clusters requested, with the important consequence that our typology identifies a distinct bifurcation in the urban-rural gradient (Fig. 7). While the urban-suburban-exurban-rural gradient is characterized linearly by decreasing developed areas and increasing forested lands, albeit with overlapping classes, our agricultural towns do not fall along this gradient. In fact, they fail to resemble rural forested towns using any of our landscape metrics, having proportions and patterns of developed areas similar to exurban towns and distinctly fewer and more isolated forested lands similar to suburban towns (Figs. 3 and S3). These findings are similar to that of Warren et al. (2018), who demonstrate empirically that population density and impervious land cover are poor predictors of forest cover versus cropland at the 'rural' end of the gradient. The distinction between agricultural towns and other town types is important if the application of the urban-rural gradient needs to capture patterns of mixed land uses across a region or understand how certain public policies might differentially affect or be affected by agricultural land use and conservation management. For example, agriculturally dominated towns often contain fewer natural lands that provide wildlife habitat. In addition, because of heavy chemical use, agricultural lands themselves affect ecosystem processes and are generally unsuitable spaces for wildlife (Henry et al. 2012;Christin et al. 2013;Knutson et al. 2018), although field margin habitat may be managed to support some species (Heath et al. 2017). Beyond this, agricultural land uses have the potential to be affected by federal and state policies in ways that nonagricultural rural areas may not. For example, many Farm Bill programs provide incentives for conservation on working farms, however, they do not provide financial assistance to large-lot residential land owners (Hellerstein 2017). It is noteworthy that our aggregated agricultural land metric encompasses a range of agricultural intensity, from lightly managed pasture to organic farms to chemically-intensive operations. It follows that, while our agricultural towns are similar in terms of a broad category, they likely vary with respect to intensity of practices. If our analysis took place in an intensive agricultural region like the Midwestern U.S., it might be of interest to incorporate multiple agricultural metrics. Given the potential difference in the environmental and social conditions between agricultural and forested rural towns, a clear distinction of the ''rural'' end of an urban-rural gradient may be warranted, as it is for our broader research objectives.

Challenges with urban-rural gradient classification approaches
Generating an urban-rural gradient typology with a cluster analysis approach, rather than simpler rulebased calculations, has advantages and disadvantages. Cluster-based typologies can easily be informed by several land cover and pattern metrics that capture greater complexity of the system and more readily address local resource management and wildlife conservation goals where habitat fragmentation, patch size, and isolation are important considerations. In comparison, using rules to define thresholds based on multiple metrics can be quite cumbersome or lead to numerous types that are difficult to distinguish consistently. Additionally, using a cluster-based approach produces classes that rescale as the extent of the focal landscape changes. For example, applying our cluster analysis process using the same metrics calculated for only the state of Maine distributed the towns among the five types differently than when we defined the typology with data from across New England (Figs. 1 and S4). Thus, cluster-based typologies produce relative definitions of the urban-rural gradient at different spatial scales. While some researchers view this feature of cluster analysis as a disadvantage (e.g., Nowosad and Stepinski 2019), we see its potential for facilitating local and regional conservation policy objectives where perceptions of urbanization matter, such as in developing landscapes (Short- Gianotti et al. 2016). In comparison, the rule-based approaches produce the same exact town grouping for the state of Maine alone and for Maine in the New England extent, i.e., nearly all towns in Maine are designated rural based on population density or no towns are designated urban based on housing density, illustrating the scale-invariant nature of rule-based approaches. Even Theobald (2005) presented two different housing density rules for defining rural town types, with the second definition identifying fewer towns as exurban and more as rural, although still failing to capture some smaller cities in our study region as described earlier (Figs. 5b and S5).
Cluster analysis does, however, have limitations. The different classes resulting from a cluster-based analysis overlap each other when illustrated on any one or pair of underlying landscape metrics (Fig. 7), while rule-based typologies, particularly those using single metrics, can easily illustrate the gradient along a single axis with clear breakpoints between classes. The more metrics considered, the more difficult it becomes to visually illustrate the gradient in a single graph. With many underlying landscape metrics and final clusters, even comparative visuals such as Fig. 3 become difficult to assimilate.
Further, cluster-based analysis requires thoughtful decisions regarding metric selection, unit of analysis, and number of clusters to generate. Numerous landscape pattern metrics exist for capturing the relationship between urbanization and ecosystem processes (Hahs and McDonnell 2006;Padilla and Sutherland 2019) and several tools exist to facilitate creation of these metrics (e.g., McGarigal et al. 2002). We used the open source GUIDOS Toolbox to develop our spatial pattern metrics because of its relative ease of use in generating multiple pattern metrics that are relevant for management and conservation of small natural features and associated wildlife. However, we ultimately did not incorporate all of the pattern metrics available in the toolbox because the resolution of our data did not discriminate among the linear pattern metrics (bridges, loops, and branches). We used webavailable land cover data at 30-m spatial resolution to define our landscape pattern metrics because it was the smallest scale consistently available across the study region. A test of the toolbox on 5-m resolution data demonstrated that the toolbox was able to more accurately generate all the pattern metrics. While the spatial resolution of our data did restrict the number and type of pattern metrics that we were able to include in our cluster process, we contend that incorporating some pattern metrics is better than none and that large and small patches were the most relevant patterns for our investigation of the management and conservation of small natural features.
We conducted our analysis at the town scale, the relevant administrative unit for land management and conservation policies our study area (New England, U.S.); however, if municipal-level detail is not relevant or if decision making takes place at a different level of governance (e.g., county), the cluster analysis could be performed at an alternative level. For example, Stokes and Seto (2019) characterized the sustainability of small ''neighborhood'' units within cities based on natural and built structures, while Boarnet et al. (2011) assessed growth management legislation at the county level.
Most cluster algorithms require the user to specify the number of clusters, which, depending on the ultimate application, may require investigation among multiple clustering solutions. A number of statistical procedures exist for identifying the ''optimal'' number of clusters in a given data set (Milligan and Cooper 1985), some of which have recently been used in landscape pattern cluster analysis (e.g., Long et al. 2010). However, in their comparison of 30 procedures for determining the correct number of clusters on an artificial data set containing distinct and non-overlapping data, Milligan and Cooper (1985) found none of the procedures were able to predict the true number of clusters 100% of the time; the best procedure was 90% accurate. Therefore, we suggest that the application of the cluster solution-in our case, the characterization of towns in New England along an urban-rural gradient for use in resource management and conservation policy making-can guide selection of the number of clusters based on somewhat subjective ''usability'' criteria. Current research defining urbanrural gradient typologies, our particular application of cluster analysis, typically use three or four classes, some subset of urban, suburban, exurban, and rural; typologies with more than four classes typically include one or more agricultural (e.g., crop or pasture) classes (Table S2). We used a combination of visual inspection of graphs and maps (Figs. 1, 2 and 3) and statistics on the underlying landscape metrics (Tables S4 and S5) to guide our final selection of the number of clusters (i.e., town types). We chose five rather than four town types because our inclusion of agricultural land cover resulted in a separate agricultural cluster with both the four-cluster and five-cluster solutions. By selecting the five-cluster solution, we were able to more fully capture the heterogeneity of land cover patterns between highly built urban centers and sparsely developed rural towns in our study area and make comparisons with other four-class urbanrural typologies that do not incorporate agricultural lands. That is, we contend that our five-cluster typology is more aligned with the four-class rulebased typologies than our four-cluster typology would have been because we can still make comparisons among the urban, suburban, exurban, and rural classes across the different typologies, as well as include agricultural towns. There are distinct and relevant differences among suburban, exurban, and rural towns in our study area (Figs. 2 and 3). Ultimately, thoughtful methodological decisions will lead to urban-rural gradient typologies that are more useful for their final application.

Conclusions
We demonstrated a multi-step cluster analysis approach that combines multiple land cover composition and pattern metrics to generate a typology that classifies towns along an urban-rural gradient. The resulting typology consists of five town types, four that coincide with other common urban-rural gradient typologies-urban, suburban, exurban, and ruraland a fifth agricultural town type seldom used in gradient typologies. Our cluster-based typology captured differences in the quantity and pattern (e.g., core versus islet patches) of developed and forested land covers across our study region (New England, USA) and identified exemplar towns for each town type that succinctly illustrate the urban-rural gradient. The pattern metrics were instrumental in differentiating among the town types, particularly those along the middle of the gradient. The distribution of towns among town types using the cluster-based approach was substantially different than the distributions resulting from four rule-based approaches based on broad measures (population density, housing density, percent forest, and percent impervious surface), suggesting that rule-based typologies are not capturing landscape patterns relevant for sustainable ecosystems. Thus, we suggest that future urbanization research use gradients developed using methods that explicitly incorporate landscape pattern metrics. While our illustration has a particular regional focus, the same cluster analysis approach can be used to define urban-rural gradients across a broad range of research objectives and geographic locations. Understanding where a town occurs along the urban-rural gradient could aid local decision-makers in prioritizing and balancing between development and conservation scenarios across diverse towns and landscapes.