Urban coyotes are genetically distinct from coyotes in natural habitats

Abstract

Urbanization is increasing throughout the world, transforming natural habitats. Coyotes (Canis latrans) are found in highly urban, suburban, rural and undeveloped mountainous habitats, making them an exemplary model organism to investigate the effects of urbanization on animals. We hypothesized that coyotes in natural habitats are more genetically related to distant coyotes in similar natural habitats and less related to coyotes in urban areas due to natal habitat-biased dispersal. We also hypothesized that increasing urbanization would result in decreased genetic diversity due to habitat fragmentation, dispersal barriers and genetic drift. We analyzed 10 microsatellite genetic markers from 125 individual coyotes sampled across a spectrum of highly urban to highly natural areas in southern California. Most coyotes clustered into four distinct genetic populations, whereas others appeared to have admixed ancestry. Three genetic populations were associated primarily with urban habitats in Los Angeles and Orange Counties. In contrast, the remaining population was associated with more naturally vegetated land near the surrounding mountains. Coyotes living in natural areas formed a genetically distinct cluster despite long geographic distances separating them. Genetic diversity was negatively associated with urban/suburban land cover and local road density, and positively associated with the relative amount of natural vegetation. These results indicate that genetic differentiation and loss of genetic diversity coincided with the extremely rapid expansion of Greater Los Angeles throughout the 1900s. Thus, urbanization reduces gene flow and erodes genetic diversity even in a habitat generalist thought to be minimally impacted by land development.

Introduction

Human activities can rapidly alter the ecology and evolution of many species. There is a growing awareness of the need to examine the ecological and evolutionary effects of one major, global, landscape-modifying activity: urbanization (Magle et al. 2012; Thompson, Rieseberg, and Schluter 2018; Rivkin et al. 2019). Cities, representing environments with the highest human footprint, impact species of wildlife by restricting movement patterns (Tucker et al. 2018), increasing exposure to infectious diseases and contaminants (Riley, Serieys, and Moriarty 2014) and rapidly changing the course of evolution (Johnson and Munshi-South 2017).

One way that urbanization influences the evolution of animals is by depleting genetic diversity, which is crucial for the long-term survival of species exposed to a changing environment. Urban populations of wildlife result from one of two processes: some individuals may colonize or recolonize an area following urban development, resulting in a founder effect; alternatively, a pre-urban population may become fragmented into smaller units during urban development. Either way, once formed, urban populations tend to be small, fragmented and isolated—conditions that intensify the effects of genetic drift. These scenarios result in a loss of genetic diversity, which can be harmful to urban populations, putting them at greater risk of local extinction. For example, pumas (Puma concolor) near Los Angeles have been subdivided into fragmented populations due to highways and other human infrastructure. A population viability analysis of pumas predicts a 99.7% probability of local extinction in 50 years if there is inbreeding depression and no dispersal among populations (Benson et al. 2016). Another urban felid in the Los Angeles area, the bobcat (Lynx rufus), experienced loss of genetic diversity following a population bottleneck caused by a disease outbreak associated with exposure to anticoagulant rodenticides (Serieys et al. 2015). Thus, urbanization may be an important driver of rapid loss of neutral genetic variation.

A second way that urbanization influences the evolution of animals is that it reduces dispersal and gene flow among populations. For example, urban highways and high traffic local roads can create physical barriers to dispersal, even in highly mobile animals (Riley et al. 2006; Thomassen et al. 2018). Similarly, features of the built environment disrupt gene flow of mice and rats, affecting their population genetic structure (Munshi-South and Kharchenko 2010; Combs et al. 2018a,b). These studies suggest that urban infrastructure imposes movement barriers on animal populations, often fragmenting them and limiting gene flow among them. Over time, with little to no gene flow, genetic drift acts differently in the isolated populations, and they become genetically distinct from one another, not requiring any kind of mutation or selection to explain the emerging differences. Although in some cases, urbanization may facilitate dispersal and gene flow leading to increased connectivity and reduced differentiation (Miles, Dyer, and Verrelli 2018), a recent review revealed that a majority of urban landscape genetic studies support the typical urban fragmentation model in which the built environment increases isolation and reduces gene flow (Miles et al. 2019). Thus, urbanization is an important driver of rapid non-adaptive genetic differentiation.

A third way that urbanization influences the evolution of animals is by introducing novel selective pressures that result in adaptation to the urban environment. The most classic example of urban adaptive evolution is the case of industrial melanism in the peppered moth (Biston betularia; Kettlewell 1955). A more modern example utilized genomics to discover numerous gene-urbanization associations in loci with small-effect alleles, indicating a polygenic response to selection in great tits (Parus major) inhabiting cities (Perrier et al. 2018). There are many other examples of phenotypic divergence in morphological, physiological, or behavioral traits between urban and non-urban populations, but these trait changes are usually linked to phenotypic plasticity and not an evolutionary (genetic) response (Miranda et al. 2013; Donihue and Lambert 2015; Alberti, Marzluff, and Hunt 2017). Nonetheless, urbanization may be an important driver of rapid adaptive evolution for some species. A fourth way that urbanization may influence the evolution of animals is by inducing mutations or altering mutation rates, but this mechanism is less clear and has received little attention (Johnson and Munshi-South 2017).

Although some species thrive in both urban and natural environments, it is important to consider the impact that habitat specialization may have on populations in a fragmented landscape (Janecka et al. 2016). Many studies consider conspecific individuals as ecologically equivalent. However, one of the primary tenets of evolution is that individuals vary. Thus, the realized niche of each individual organism is slightly different, and some individuals are better suited to survive in certain environments within the theoretical niche of the species (Bolnick et al. 2003). Individual specialization is particularly important in species where individuals learn many behaviors from their parents. One consequence of this type of behavioral imprinting passed down from parents to offspring is the induction of habitat preferences. These preferences, in turn, can result in natal habitat-biased dispersal, the tendency of an animal to disperse to an area resembling its natal habitat as opposed to a more unfamiliar one due to the disperser being better behaviorally predisposed to that specific environment. This has long been recognized as a source of individual variation in habitat selection (Immelmann 1975; Davis and Stamps 2004). If natal habitat-biased dispersal occurs in a species that inhabits both urban and natural habitats, individual variation in habitat preferences may lead to genetic differences between urban and non-urban populations. Such differences imply a reduction of gene flow as a result of barriers to dispersal (e.g. physical, ecological or behavioral).

Coyotes (Canis latrans) thrive in many cities throughout North America (Gehrt 2007; Nagy et al. 2017). Urban coyotes generally occupy smaller home ranges, eat more anthropogenic food items, are more active at night and are bolder and more exploratory than their relatives in natural areas (Gehrt 2007; Breck et al. 2019). Greater Los Angeles is the second largest metroplex in the USA, encompassing several interconnected cities in southern California, with a combined human population exceeding 18.7 million (U.S. Census Bureau 2019). The landscape of Greater Los Angeles and its surroundings consists of highly urban, suburban, rural and wild mountainous habitats. Coyotes are found throughout all terrestrial habitats of Greater Los Angeles and are thus an exemplary model organism for investigating the effects of urbanization on wildlife.

The purpose of this study was to investigate whether urban coyotes differ genetically from coyotes in natural habitats. Specifically, we surveyed 10 microsatellite loci to quantify population structure and genetic diversity of coyotes inhabiting a mosaic of various habitat types differing in the degree of urbanization across Greater Los Angeles. We hypothesized that due to natal habitat-biased dispersal, coyotes living in natural habitats are more genetically related to other coyotes in similar natural habitats and less related to coyotes living in urban habitats, despite geographic distance. Further, we hypothesized that as urbanization increases, the genetic diversity of coyotes would decrease, likely as a consequence of habitat fragmentation, reduced gene flow and increased genetic drift. Our study investigates the fine-grained distribution of genetic variation in coyotes across a major urbanized area. More generally, this study contributes to our understanding of how urbanization affects the ecology and microevolution of wildlife.

Methods

Study location and sampling

The study site was located in the Greater Los Angeles area of southern California. We obtained genetic samples of coyotes (n = 134) from the National Park Service and the tissue depository at the University of California, Irvine. Coyotes were sampled from Los Angeles (n = 50), Orange (n = 56), Ventura (n = 16), San Bernardino (n = 4), Riverside (n = 6) and San Diego (n = 2) counties. The georeferenced localities of samples represent a spread along a continuum of land cover from naturally vegetated wildland to land developed for agriculture or urban/suburban use. All samples collected by the National Park Service were approved by the California Department of Fish and Wildlife under scientific collecting permit number SC-5636 and met all animal safety standards from the American Society of Mammalogists guidelines for the use of wild mammals in research (Sikes and Animal Care and Use Committee of the American Society of Mammalogists 2016). Other samples retrieved from UC Irvine came from animals killed for reasons other than research, such as roadkill or coyotes killed by private landowners.

Laboratory methods

We used the Qiagen DNeasy Blood and Tissue Kit to extract total genomic DNA from either coyote muscle, liver or skin tissue. We initially homogenized tissues with metal beads in a bullet-blender and then performed the actual DNA extraction according to the manufacturer’s instructions. We purified the extraction using the Quick-Precip protocol (EdgeBio) and quantified DNA concentrations using a NanoDrop ND-2000 Spectrophotometer (Thermo Fisher Scientific). To minimize variation of input concentrations, we diluted samples with >500 ng/µl 5-fold and samples with 200–500 ng/µl 3-fold. Average final concentration of DNA samples was 109 ng/µl.

We genotyped the following 10 microsatellite loci: CXX172, FH2001, FH2054, FH2062, FH2096, FH2137, FH2140, FH2226, FH2293 and FH2535 (Ostrander, Sprague, and Rine 1993; Francisco et al. 1996; Mellersh et al. 1997; Breen et al. 2001). Marker CXX172 contains a dinucleotide repeat, while all other markers have tetranucleotide repeats. All 10 microsatellite markers were designed for the domestic dog (Canis familiaris) and map to 10 distinct chromosomes of the canine genome (Supplementary Table S1). We used the polymerase chain reaction (PCR) to individually amplify each microsatellite locus across all samples. The total PCR volume of 25 µl included 12.5 µl Hot-Start Taq Blue Master Mix (Denville Scientific), 1 µl fluorescently tagged forward primer, 1 µl reverse primer, 2 µl Bovine Serum Albumin, 6.5 µl H₂O and 2 µl template DNA. We performed PCR amplifications in Bio-Rad T100 Thermal Cyclers with the following conditions: initial denaturation for 10 min at 95°C; 45 cycles of 1 min at 95°C, 1 min at 58°C (FH markers) or 1 min at 60°C (CXX marker) and 1 min at 72°C; and final extension for 10 min at 72°C. We visualized amplified fragments using gel electrophoresis on a 1.5% agarose gel stained with GreenGlo dye (Denville Scientific). We diluted the amplified DNA solution to a uniform concentration based on the intensities of the bands observed on the gel images. We then prepared a genotyping plate by adding 1 µl of diluted amplified DNA solution, 8.8 µl of Hi-Di Formamide (Thermo Fisher Scientific) and 0.2 µl of 400HD ROX size standard (Thermo Fisher Scientific) for a total volume of 10 µl per well. We performed DNA fragment analysis by capillary electrophoresis using an ABI Genetic Analyzer 3500 (Thermo Fisher Scientific) and made allele calls using the ABI GeneMapper software.

Analytical methods

Analyses of isolation by distance and population genetic structure may be sensitive to family structure. For example, unsupervised clustering algorithms may overestimate the number of genetic populations if there are highly related individuals in a sample (Anderson and Dunham 2008). Therefore, we first computed the Lynch and Ritland (1999) estimator of pairwise relatedness across all coyotes in GENALEX 6.5 (Peakall and Smouse 2006, 2012), and then removed one coyote from each pair of first-order relatives (r = 0.5 for parent–offspring or full-sib relationships). Second-order relatives are most likely grandparent–grandoffspring pairs; it is unlikely that they are half-siblings because coyotes form long-term monogamous pairs (Hennessy, Dubach, and Gehrt 2012; Henger et al. 2020). Either way, they would not be living together in a family group, as juvenile coyotes typically disperse from natal territories within their first year, before any grandoffspring or half-sibs would be born. Additionally, even a large proportion of half-sibs in a sample failed to induce spurious population structure (Anderson and Dunham 2008). Thus, we retained second-order relatives in analyses of population structure and all other downstream analyses.

We used a Mantel test of isolation by distance to examine whether coyotes that are farther apart are more distantly related. We used GENALEX to examine the correlation between linear codominant genotypic distance and linear geographic distance calculated from Lat/Long coordinates of individuals, applying 999 random permutations of the genetic distance matrix to test the significance of the correlation.

To examine the pattern of genetic structure more closely, we performed a maximum-likelihood model-based estimation of ancestry and population genetic structure to assign individuals to genetic populations using STRUCTURE 2.3.4 (Pritchard, Stephens, and Donnelly 2000). STRUCTURE implements a Bayesian algorithm to assign multilocus genotypes to genetic clusters by calculating the likelihood that a group of individuals constitutes a population. As geographic location is not incorporated in the STRUCTURE analysis, population assignments can be visualized on a map and accepted as an unbiased display of population genetic structure. For the STRUCTURE analyses, we used 10 replicate runs of 20 000 burn-in and 100 000 Markov chain Monte Carlo iterations and the admixture ancestry model with correlated allele frequencies (Falush, Stephens, and Pritchard 2003). We set the number of possible populations ranging from K = 1 to 10 and verified that likelihood statistics and alpha reached convergence during burn-in for each value of K. We used STRUCTURE HARVESTER (Earl and von Holdt 2012) to determine the most likely number of genetic populations, K, by plotting Ln P(D), the mean posterior probability of the data and Delta K, a value related to the second-order rate of change of the likelihood function with respect to K (Evanno, Regnaut, and Goudet 2005). We aligned and averaged the 10 replicate cluster membership coefficient matrices using CLUMPP (Jakobsson and Rosenberg 2007). We considered individuals with >70% ancestry assigned to any given population as members of that population, and those with <70% ancestry assigned to any given population as admixed. Higher thresholds of assignment (e.g. 75% and 80%) did not qualitatively alter results in downstream analyses but reduced the size of the genetic populations and increased the size of the admixed group.

We calculated the pairwise genetic differentiation (F_ST) among the STRUCTURE-defined populations and admixed group, and then tested the significance of the observed genetic differences by an analysis of molecular variance (AMOVA) using 999 permutations in GENALEX. To examine whether the STRUCTURE results reflect family units, we computed the Lynch and Ritland estimator of pairwise relatedness within populations, applying 999 random permutations and 1000 bootstraps in GENALEX.

We mapped the geographic distribution of highly assigned and admixed coyotes using Google Earth and QGIS. We used the raster (Hijmans 2017), maptools (Bivand and Lewin-Koh 2017) and rgdal (Bivand, Keitt, and Rowlingson 2018) R packages and the USGS 2011 Land Cover Database to extract land cover data in 30 × 30-m resolution from 5-km² buffers around each georeferenced sampling locality, corresponding to the average home range of a coyote in southern California (Riley et al. 2003). We reduced the proportions of the various USGS land cover classes to three composite categories: (i) urban/suburban developed, land with 20–100% impervious surface cover; (ii) agricultural/open developed, land used for cultivation or with <20% impervious surface cover; and (iii) natural vegetation, areas dominated by forest, shrub scrub, grassland and wetland. We extracted the total length of local roads and highways from 5-km² buffers around each georeferenced sampling locality as a second metric of urbanization. Road density is an appropriate and independent measure of urbanization as it encompasses the majority of open space available for a coyote to roam in an urban environment. We compared the extent of urban/suburban developed land cover, natural land cover, road density and highway density among the STRUCTURE-defined populations and admixed group using non-parametric Kruskal–Wallis tests.

We used GENALEX to calculate standard metrics of genetic diversity, such as observed number of alleles (A_N), effective number of alleles (A_E), number of private alleles (A_P), Shannon’s diversity information index (I), observed heterozygosity (H_O), unbiased expected heterozygosity (H_E) and inbreeding coefficient (F_IS). We used analysis of variance (ANOVA) to test for significant differences in A_N, A_E, I, H_O, H_E and F_IS among populations. Since genetic diversity metrics are highly sensitive to sample size, we also performed a rarefaction analysis for allelic richness to account for differences in population size using the program ADZE 1.0 (Szpiech, Jakobsson, and Rosenberg 2008).

Finally, we compared patterns of observed genetic diversity to several physiographic features including land cover and road density. We performed regression analyses to test for associations between indices of diversity and indices of urbanization. Response variables included metrics computed at the population level (e.g. observed number of alleles) and at the individual level (i.e. percentage of heterozygous loci). Explanatory variables included percent urban/suburban land cover, percent natural land cover, percent agricultural land cover, local road density and highway density within the 5-km² buffers around each sampling locality.

Results

We genotyped 134 samples at 10 microsatellite loci. Our quality control measures filtered out eight samples that performed poorly with either low amplification or genotyping rates. The analysis of pairwise relatedness revealed that one pair of coyotes were first-order relatives (r = 0.5). We removed one coyote from this pair and conducted all downstream analyses with the remaining 125 samples, which overall had a total genotyping rate of 97% (Supplementary Material).

The Mantel test showed a small and marginally significant correlation between genetic and geographic distances (r_M = 0.125, P = 0.050), suggesting a slight signal of isolation by distance. STRUCTURE revealed a strong signal that K = 4 was the most likely number of genetically distinct populations among the coyotes sampled (Fig. 1A). Other values of K had little support from the Delta K analysis and generated clusters with very few coyotes (Supplementary Fig. S1). Thus, we retained K = 4 as the most parsimonious and strongly supported population structure for downstream analyses. Of the 125 coyotes analyzed, 88 individuals assigned strongly to one of the four STRUCTURE-defined populations; the rest were admixed (Fig. 1B). Based on Q coefficients, the ancestry of the 37 admixed coyotes was assigned as follows: 23% to Population 1, 11% to Population 2, 40% to Population 3 and 26% to Population 4. This indicates that a large portion of these admixed individuals were more similar to coyotes in Population 3.

Figure 1:

(A) Evanno plot of the STRUCTURE analysis indicated that K = 4 is the most likely number of genetically distinct populations. (B) Genetic ancestry of 125 coyotes sampled in this study in relation to K = 4. Each vertical bar represents an individual coyote. In total, 88 individuals had >70% assignment to one of four genetically distinct populations.

AMOVA showed that a significant amount of genetic variance exists among the four populations and the admixed group (Global F_ST = 4%, P = 0.001). This value was relatively small; however, all four populations, as well as the admixed coyotes (grouped as an admixed ‘population’), were significantly differentiated from each other. Pairwise F_ST values among the four populations ranged from 0.042 to 0.076, and pairwise F_ST comparisons with the admixed group were lower, ranging from 0.011 to 0.025 (Table 1). The lowest F_ST value was between Population 3 and the admixed group. The average pairwise relatedness of individuals within the four populations was significantly >0 (P = 0.001) but still lower than expected for family units, ranging from r = 0.009 to 0.060; the average pairwise relatedness of admixed coyotes was not significantly different from 0 (P = 0.364; Supplementary Table S2 and Fig. S2).

Table 1:

Pairwise genetic differentiation among the four STRUCTURE-defined populations and the admixed coyotes (grouped together)

	Pop 1	Pop 2	Pop 3	Pop 4	Admixed
Pop 1	–
Pop 2	0.074	–
Pop 3	0.042	0.051	–
Pop 4	0.070	0.076	0.047	–
Admixed	0.022	0.023	0.011	0.025	–

Values of the fixation index F_ST are shown below the diagonal. All pairwise comparisons are statistically significant (P < 0.005) based on 999 AMOVA permutations.

Population 3 was associated with areas that have more naturally vegetated land near the Santa Monica Mountains, Santa Susana Mountains, Santa Ana Mountains, Santa Rosa Mountains and Simi Hills (Fig. 2). Coyotes in Population 3 were widespread, separated by up to 271 km. Population 3 had the most variability in land cover, with nine coyotes sampled from areas with <10% urban land cover and five coyotes sampled from areas with >90% urban land cover. Nonetheless, Population 3 had significantly less urban/suburban land cover (Kruskal–Wallis χ² = 13.426, df = 4, P = 0.009) and significantly more naturally vegetated land cover (χ² = 11.467, df = 4, P = 0.022) within buffers of 5 km² (Fig. 3A, B and Supplementary Table S3). Because Population 3 was geographically widespread and had a high variance of land cover, we conducted a STRUCTURE analysis of only this group to explore whether there is substructure or geographic clustering. In this analysis, K = 1 had the most support, thus there was no evidence of substructure or geographic clustering within Population 3 (Supplementary Fig. S3).

Figure 2:

Spatial distribution of the K = 4 genetic populations and admixed coyotes. Sampling localities mapped onto (A) Google Earth satellite imagery, and (B) USGS Land Cover satellite imagery (red = urban/suburban developed, yellow = agricultural/open developed, light brown = naturally vegetated). Colored markers represent individual coyotes assigned to one of four populations; black markers represent admixed individuals.

Figure 3:

Boxplots showing the distribution of urban/suburban land cover, natural land cover and road density within 5-km² buffers around georeferenced sampling localities for each genetic population and the admixed group. Differences among populations were analyzed by non-parametric Kruskal–Wallis tests. (A) Population 3 had significantly less urban/suburban land cover (χ² = 13.426, df = 4, P = 0.009). (B) Population 3 had significantly more naturally vegetated land cover (χ² = 11.467, df = 4, P = 0.022). (C) Population 2 had significantly greater local road density (χ² = 19.574, df = 4, P = 0.001). Landscape icons made by Freepik from flaticon.com.

In contrast to Population 3, the remaining Populations (1, 2 and 4) were associated primarily with highly urbanized habitats in Los Angeles and Orange Counties, and the admixed group was distributed in areas with an intermediate amount of urbanization. Populations 1 and 4 had very high levels of urban/suburban land cover and were primarily centered in Orange County. The roughly parallel California State Routes 55 and 261, and the intense commercial development between them, separate these nearby populations. Population 2 was located primarily in the highly urbanized areas of Downtown Los Angeles and San Fernando Valley that surround Griffith Park, the site of the iconic Hollywood Sign (Fig. 2). Three distant individuals assigned to Population 2 were sampled at Long Beach Airport, a suburban road in Agoura Hills and a golf course in Ontario. Coyotes in this population had significantly greater local road density compared with the other populations and admixed group (χ² = 19.574, df = 4, P = 0.001; Fig. 3C andSupplementary Table S3). Highway density did not differ among the four populations and admixed group (χ² = 4.986, df = 4, P = 0.289).

Population 3, the population associated with more naturally vegetated habitat, had the highest allelic richness (A_N, correlated samples ANOVA F = 12.07, P < 0.0001), the highest effective number of alleles (A_E, F = 4.73, P = 0.004), the highest number of private alleles (A_P) and the highest Shannon’s diversity information index (I, F = 6.17, P = 0.001; Table 2 and Supplementary Fig. S4). In comparison to the four populations defined by STRUCTURE, the admixed group showed the second highest values for allelic richness, effective number of alleles, number of private alleles and the Shannon diversity information index (Table 2). There were no significant differences among the four populations and admixed group in observed heterozygosity (H_O, F = 0.60, P = 0.665), expected heterozygosity (H_E, F = 0.91, P = 0.469) or coefficient of inbreeding (F_IS, F = 0.76, P = 0.558). Population 3 was the most genetically diverse population, but it was also the population with the largest sample size. However, the rarefaction analysis for allelic richness showed that the three urban populations (1, 2 and 4) had lower genetic diversity compared with Population 3, the more ‘natural’ population, regardless of sample size, and the admixed group’s allelic richness was intermediate (Fig. 4).

Figure 4:

Rarefaction analysis of the effect of sample size on estimates of genetic diversity. Populations 1, 2 and 4 were significantly more urban and had lower genetic diversity regardless of sample size.

Table 2:

Standard metrics of genetic diversity by population averaged across 10 loci, including number of samples (n), observed number of alleles (A_N), effective number of alleles (A_E), total number of private alleles (A_P), Shannon’s diversity information index (I), observed heterozygosity (H_O), unbiased expected heterozygosity (H_E) and inbreeding coefficient (F_IS)

	n	AN	AE	AP	I	HO	HE	FIS
Pop 1	26	5.9	3.8	2	1.409	0.703	0.704	−0.033
Pop 2	12	5.2	3.8	0	1.361	0.656	0.736	0.074
Pop 3	31	10.0	6.0	17	1.796	0.743	0.760	0.018
Pop 4	19	5.4	3.7	0	1.394	0.737	0.719	−0.067
Admixed	37	9.3	5.4	6	1.696	0.717	0.747	0.007

Population 3 had the highest A_N, A_E, A_P, I, H_O and H_E. Shaded cells represent significantly different values.

Genetic diversity (measured as the average number of alleles in a population) decreased as the degree of urbanization increased, with 1.6 alleles lost for every 10% increase in urbanization (R² = 0.790, P = 0.044; Fig. 5). Conversely, genetic diversity increased as the degree of natural land cover increased (R² = 0.772, P = 0.049) but was not associated with the degree of agricultural land cover (R² = 0.656, P = 0.097). Genetic diversity (measured as observed heterozygosity of individual coyotes) was not associated with local road density (R² = 0.012, P = 0.226) or with highway density (R² = 3.13 × 10⁻⁴, P = 0.847).

Figure 5:

Regression analysis of association between genetic diversity and urbanization. Genetic diversity (measured as the number of alleles in a population) decreased as the degree of urbanization increased (R² = 0.790, P = 0.044)

Discussion

This study demonstrates that fine-grained patterns of urbanization have a noticeable effect on the population genetic structure and diversity of coyotes. We found evidence of four distinct genetic populations within Greater Los Angeles. Three populations with low genetic diversity were associated primarily with urbanized habitats in Los Angeles and Orange Counties. The remaining geographically widespread and genetically diverse population was mostly associated with a mosaic of naturally vegetated mountainous areas throughout the entire study site. Nearly, 30% of the coyotes sampled did not assign strongly to any of the four genetic populations; these admixed individuals were mostly found in areas with intermediate degrees of urbanization and collectively had an intermediate amount of genetic diversity. These results support our hypothesis that the population genetic structure of coyotes in Greater Los Angeles is influenced by urbanization and that coyotes living in natural areas form a genetically distinct cluster despite long distances separating them. Additionally, these results support our hypothesis that genetic diversity in coyotes is negatively associated with urbanization.

Early studies of genetic variation in coyotes did not detect population structure, even at continental scales (Lehman and Wayne 1991; Roy et al. 1994). This apparent panmixia was attributed to the generalist and vagile nature of coyotes. However, more recent studies have detected population structure at smaller areas sampled in greater detail. In western North America, genetically distinct coyote populations were subdivided according to macrohabitat discontinuities. Population structure of coyotes in central California conformed to major contiguous ecoregions such as the Sierra Nevada Mountains, the Great Valley and the low-elevation coastal mountain ranges (Sacks, Brown, and Ernest 2004; Sacks et al. 2005). Similarly, coyotes in the arid Southwest and western California formed a genetic cluster despite vast distances separating them, and this cluster was distinct from coyotes in the more heterogeneous ecoregions of eastern California (Sacks et al. 2008). Genetic variation in northeastern coyotes was also highly structured, especially in the heterogeneous landscape of central New York and Pennsylvania, where two colonization fronts recently merged (Monzón 2014). The spatial distribution of genetic variation in northeastern coyotes did not correspond to contiguous ecoregions (Berkman et al. 2019) but was partly explained by heterogeneity in prey density and the relative proportions of forest cover, agriculture and urban/suburban land use (Monzón 2012; Monzón, Kays, and Dykhuizen 2014). This correspondence of genetic subdivisions with habitat heterogeneity may be a consequence of natal habitat-biased dispersal, suggesting that ecological factors drive genetic differentiation in coyotes and other highly vagile canids (Karlin and Chadwick 2012; Sanz-Pérez et al. 2018).

We observed a similar congruence of genetic structure with habitat heterogeneity, but our study examined population genetic structure at a much larger scale (i.e. smaller extent and finer grain) than previous investigations. Whereas the aforementioned studies surveyed genetic structure in coyote populations across the species range or a region or a state, we surveyed genetic structure in coyotes living in one urban area and its immediate surroundings. Moreover, our study explicitly considered urbanization as a quantitative variable in the description of coyote habitat. Our results showed that ecological variation can induce genetic differences when it occurs in a complex, fine-grained mosaic, and not just when it occurs in broad macrohabitat configurations. Part of this ecological variation involves urban land use, which has not been quantitatively considered in previous genetic studies of coyotes.

We detected higher levels of genetic differentiation compared with other studies of urban coyotes (Table 3). Coyotes that recently colonized New York City (NYC) were slightly differentiated from non-NYC coyotes (DeCandia et al. 2019). Urban coyotes in the small city of Auburn, Alabama were similarly slightly differentiated from coyotes in the surrounding rural areas (Damm et al. 2015). Coyotes separated by a freeway near Los Angeles had similar F_ST values (Riley et al. 2006). Coyotes in the suburbs surrounding Cleveland, OH did not freely intermix (Rashleigh, Krebs, and van Keulen 2008); however, F_ST in that study was calculated from mitochondrial DNA haplotypes and is, therefore, not comparable to values calculated from nuclear microsatellite genotypes. The F_ST values reported in this study of coyotes in Greater Los Angeles are even higher than those reported in regional comparisons. For example, coyotes in the northeastern USA were only slightly differentiated from coyotes in the southeastern USA: F_ST = 0.022 (Hinton et al. 2019). The high values of genetic divergence reported in our study are consistent with those reported for other vertebrate species, including a volant bird, inhabiting fragmented urban areas in southern California (Delaney, Riley, and Fisher 2010). We conclude that high levels of genetic differentiation among animal populations in Greater Los Angeles are evidence of a reduction in functional connectivity and permeability of the cityscape (Beninde et al. 2016).

Table 3:

Range of genetic differentiation values measured in various studies of urban coyotes

Urban context	Number of populations	Sample size	Genetic markers	FST range	References
Los Angeles, CA	4, defined analytically	125	10 microsatellites	0.042–0.076	This study
NYC, NY	2, pre-defined as NYC versus non-NYC	105	9 microsatellites	0.035	DeCandia et al. (2019)
Auburn, AL	2, pre-defined as urban versus rural based on land cover	58	10 microsatellites	0.030	Damm et al. (2015)
Los Angeles, CA	3, pre-defined based on roads	68	7 microsatellites	0.030–0.037	Riley et al. (2006)
Cleveland, OH	3, pre-defined based on rivers	57	mtDNA D-loop (531 bp)	0.07–0.17	Rashleigh, Krebs, and van Keulen (2008)

The presence of admixed coyotes suggests that despite the significant genetic structure observed, appreciable gene flow is still occurring between urban and rural populations. Based on the STRUCTURE and pairwise F_ST analyses, admixed coyotes appeared more similar to those in Population 3, which was widespread and mostly associated with natural areas. Interestingly, the admixed coyotes also were widespread, but they were distributed in areas with an intermediate amount of urbanization and possessed intermediate genetic diversity. There seems, therefore, to exist a gradient of gene flow from natural to urban areas in Greater Los Angeles. We propose that comparable gradients of gene flow exist with a distinct genetically admixed group linking natural and urban population of coyotes and other generalist species throughout the world’s cities.

Finding three genetically distinct populations with similar urbanization levels was surprising. Something other than differences in land cover must account for the genetic distinctions. One population near Downtown Los Angeles occupied a hyper-urban area with very high local road density. Two other populations were located primarily in Orange County and were separated by a putative geographic barrier, the high-intensity commercial and industrial development sandwiched between California State Routes 55 and 261. We speculate that this semi-permeable barrier made of two freeways reduces dispersal sufficiently to maintain the genetic distinctness between coyote populations on either side (genetic differentiation F_ST = 0.070 and the corresponding number of migrants per generation N_m = 3.306). In fact, freeways are physical dispersal barriers and restrict gene flow in urban carnivores of Greater Los Angeles (Riley et al. 2006, 2014b; Serieys et al. 2015; Benson et al. 2016). Future studies could investigate more exhaustively the factors causing the genetic divergence of adjacent urban populations.

The coyote is a generalist species with great movement potential; thus, it is unlikely that there are substantial topographic barriers to dispersal between natural and urban areas in southern California. Nonetheless, individual coyotes are not necessarily generalists but may instead specialize on resources specific to one habitat type or another (Sacks et al. 2008). A recent study of diet in coyotes of Greater Los Angeles revealed that individual resource specialization is more common in rural and suburban habitats (Larson et al. 2020). The populations in our study differed in the degree of natural land cover, urban/suburban land cover and local road density, suggesting that these environmental factors influence how individuals cue in on suitable habitat and move, reproduce, and utilize resources within it. However, the perception of suitable habitat is likely to be induced or imprinted by early life experiences (Davis and Stamps 2004). For example, a coyote born in an urban territory must learn crucial survival skills from its parents to become adept at navigating areas with many roads, cars, people, pollutants and other hazards to wild animals. These skills uniquely equip the young coyote to gain the benefits of urban living while avoiding its risks (Gehrt, Brown, and Anchor 2011). Adult coyotes can become increasingly habituated to anthropogenic disturbance and pass on their reduced fear of humans to their offspring (Schell et al. 2018). Thus, urban coyotes presumably transmit urban survival skills to their pups, and once the pups grow and disperse, they settle in an urban habitat similar to the one that they were trained in by their parents. Over time, natal habitat-biased dispersal would result in urban coyote populations becoming genetically distinct from populations in natural habitats.

This congruence of population genetic structure with habitat type would most likely emerge as a response to neutral evolutionary forces such as genetic drift and limited gene flow, but local adaptation to human-modified environments is also a possibility. High dispersal and gene flow are often thought to prevent local population divergence, but if dispersal is nonrandom in relation to environmental heterogeneity it may actually reinforce evolutionary differentiation and promote the emergence of locally adapted ecotypes (Garant et al. 2005; Kays and Monzón 2017). Detailed studies of coyote movements in Greater Los Angeles are necessary to confirm whether juvenile dispersal is influenced by levels of urbanization in the natal habitat. Future genetic studies may focus on surveying non-neutral markers to explore adaptive evolution in urban coyotes (e.g. DeCandia et al. 2019). Also, a landscape genomics approach can help to tease apart the signatures of genetic drift in small urban populations, reduced gene flow from habitat-biased dispersal and natural selection of urban-favored genotypes (Schell 2018; Miles et al. 2019; Rivkin et al. 2019).

Urbanization significantly reduced the genetic diversity of coyotes. Our results support the hypothesis that coyote populations experienced fragmentation and subsequent reduction in gene flow between isolated groups during the extremely rapid expansion of Greater Los Angeles beginning in the late 1800s. Urban coyotes in NYC and in Auburn, Alabama also had lower genetic diversity compared with their non-urban counterparts (Damm et al. 2015; DeCandia et al. 2019). However, it is important to note that coyotes inhabiting NYC and Auburn began to colonize these cities only in the past few decades (Damm et al. 2015; Nagy et al. 2017). Although these studies did not quantify levels of urbanization for each coyote as we did, the results are consistent. Urbanization erodes genetic diversity of both recent and long-established coyote populations at the extremes of the species range. This is not surprising; urbanization and other forms of anthropogenic landscape modifications represent loss and fragmentation of habitat and are usually negatively associated with genetic diversity in various taxa (Janecka et al. 2014; Munshi-South, Zolnik, and Harris 2016; Johnson and Munshi-South 2017).

It is very important to monitor the effects of human activities on the genetic diversity of animals because it is an important indicator of the health of a population. We found that urbanization negatively impacts genetic diversity even in a generalist carnivore that is thought to be minimally impacted by land development. The combined effects of urbanization (loss of genetic diversity and ensuing population genetic structure) are demonstrably impacting the microevolution of coyotes and presumably also many other species living in and around the world’s cities. Further studies should be conducted on other urban species to better understand the ecological and evolutionary effects of urbanization.

Data availability

Spatial and genetic data for 125 coyotes used in this study and genetic diversity metrics summarized for each locus are available as Supplementary Material.

Supplementary data

Supplementary data are available at JUECOL online.

Acknowledgements

We thank Lee Kats, Niamh Quinn, Ted Stankowich, Ben Sacks, Michael Gribble, Morgan Clark and Nicole Lopez for their assistance and advice with this project.

Funding

This research was funded by the National Science Foundation Research Experiences for Undergraduates award number 1560352, the Natural Science Division of Pepperdine University, the Pepperdine Academic Year Undergraduate Research Initiative (AYURI) grant to J. Monzón, and the Dean’s Research Grant to R. Honeycutt.

Conflict of interest statement. None declared.