https://orcid.org/0000-0002-5604-2942
Abstract
Island spotted skunks (Spilogale gracilis amphiala) are a rare subspecies endemic to the California Channel Islands, currently extant on Santa Cruz and Santa Rosa islands. How and when skunks arrived on the islands is unknown, hindering decision-making about their taxonomic status and conservation priority. We investigated these questions by sequencing the complete mitochondrial genomes of 55 skunks from the two islands and mainland (California and Arizona) and examining phylogenetic patterns and estimations of isolation times among populations. Island spotted skunks grouped in a single monophyletic clade distinct from mainland spotted skunks. A haplotype network analysis had the most recent common ancestral haplotype sampled from an individual on Santa Rosa, suggesting both islands were colonized by a single matriline. Additionally, no haplotypes were shared between skunk populations on the two islands. These patterns imply that both island populations were derived from a common ancestral population shortly after establishment and have remained isolated from each other ever since. Together with divergence estimates from three methods, this topology is consistent with colonization of the super-island, Santarosae, by a single ancestral population of spotted skunks in the early Holocene, followed by divergence as the sea level rose and split Santarosae into Santa Cruz and Santa Rosa islands 9,400–9,700 years ago. Such a scenario of colonization could be explained either by rafting or one-time transport by Native Americans. Given their distinct evolutionary history, high levels of endemism, and current population status, island spotted skunks may warrant management as distinct evolutionarily significant units.
Islands contribute substantially to global biodiversity because they support high levels of endemism (Whittaker et al. 2017; Veron et al. 2019). Consequently, understanding the processes that underlie endemism and speciation on islands is of particular interest (Rosauer and Jetz 2015; Palombo 2018; Veron et al. 2019). Studying the evolutionary and biogeographical mechanisms that give rise to island endemism can help explain global distributions of biodiversity (Veron et al. 2019). Further, investigating the evolutionary history of insular endemic species may help resolve questions of taxonomic status and guide decisions about conservation priorities (Cadotte and Davies 2010). However, island faunas are often depauperate and unbalanced, with fewer terrestrial, flightless mammals as compared with the mainland (MacArthur and Wilson 1967; Alcover et al. 1998; Palombo 2018), limiting opportunities to study island endemism for some taxa. For example, carnivores rarely occur naturally on islands (Williamson et al. 1989; Alcover and McMinn 1994; Palombo 2018), highlighting the importance of understanding the biogeography and evolutionary histories of insular endemic carnivores.
The California Channel Islands support two endemic carnivores, the island fox (Urocyon littoralis) and the island spotted skunk (Spilogale gracilis amphiala). The island fox is a species that occurs on six of the eight largest islands with each of the island populations considered a separate subspecies (Gilbert et al. 1990; Wayne et al. 1991; Collins 1993; Rick et al. 2009). Island spotted skunks are currently found on only two of the northern Channel Islands, Santa Cruz and Santa Rosa, where they co-occur with island foxes, although historically they may have also occurred on San Miguel Island (Walker 1980). Based on morphological measurements, island spotted skunks are considered a subspecies of western spotted skunks (S. gracilis) found on mainland California (Van Gelder 1959). For both carnivores, the evolutionary history, including timing and potential method of island colonization, is of great interest. However, while island foxes have been the focus of multiple studies of morphology, phylogeny, and genetic composition, island spotted skunks remain understudied.
Foxes and skunks are thought to have coexisted on Santa Cruz and Santa Rosa islands for at least several thousand years (Floyd et al. 2011; Goddard et al. 2015; Hofman et al. 2015), an idea supported by extensive investigation of island foxes but based on more speculative evidence for island skunks. Analysis of mitochondrial fragments suggests that foxes colonized the northern islands approximately 7,300 to 19,700 years ago (Goddard et al. 2015). This time frame agrees with radiocarbon dating of bones found at archeological sites that estimates foxes have been on the northern islands for at least 7,000 years (Hofman et al. 2015, 2016). In contrast to foxes, the timing of spotted skunk arrival on the Channel Islands is unknown. Based on two recent phylogenies, island and mainland spotted skunks are considered a part of the western clade of western spotted skunks, which diverged from other groups of western spotted skunks 480,000–820,000 years ago (McDonough et al. 2022; but see Ferguson et al. [2017] for a different time frame of divergence). However, both the Ferguson et al. (2017) and McDonough et al. (2022) phylogenies lacked the resolution or sample sizes necessary to investigate the divergence time of island spotted skunks from the nearest western spotted skunks on mainland California. Furthermore, there is sparse evidence of spotted skunks in the paleontological or archeological record for California or the Channel Islands (Orr 1968; Rick 2013) that might help resolve uncertainty surrounding the timing of island colonization.
The ancestors of island foxes and island spotted skunks likely arrived on the northern Channel Islands either by natural dispersal on floating debris from the mainland or via introduction by humans (Walker 1980; Wenner and Johnson 1980; Floyd et al. 2011; Goddard et al. 2015; Hofman et al. 2015). The Channel Islands currently group into two geographic clusters, with a northern group (Anacapa, Santa Cruz, Santa Rosa, and San Miguel) off the coast of Santa Barbara, California, and a southern group (San Clemente, San Nicolas, Santa Barbara, and Santa Catalina) off the coast of Los Angeles and San Diego, California. Throughout the Quaternary Period, the northern Channel Islands have not been connected to the mainland via a land bridge (Johnson 1983); hence, overland colonization by foxes or skunks is unlikely. However, during the last glacial maximum, lowered sea levels connected the northern Channel Islands into one large island, Santarosae (Orr 1968). Santarosae was at its maximum size approximately 20,000 years ago and was separated from the mainland by ca. 7 km (Orr 1968; Reeder-Myers et al. 2015). Subsequently, rising sea levels resulted in the gradual separation of Santarosae into the four northern islands 9,000–11,000 years ago, with Santa Cruz Island separating from Santa Rosa Island 9,400–9,700 years ago (Reeder-Myers et al. 2015). During this time, distance from the mainland was <12 km (Reeder-Myers et al. 2015). By approximately 7,000 years ago, the islands reached their modern configuration (Reeder-Myers et al. 2015), leaving Santa Cruz and Santa Rosa islands 30–45 km from the mainland and ca. 9 km apart. Many terrestrial carnivores are not proficient oceanic swimmers (Palombo 2018), but if foxes or skunks colonized via rafting on floating debris, this would have most likely occurred prior to the separation of Santarosae, when the distance between the mainland and super-island was at its lowest (Pergams and Ashley 1999; Floyd et al. 2011; Reeder-Myers et al. 2015). Alternatively, one or both carnivores may have been deliberately or inadvertently translocated to the northern Channel Islands by humans (Wenner and Johnson 1980). Humans reached the Channel Islands approximately 13,000 years ago, and Native Americans inhabited the northern islands consistently while traveling back and forth to the mainland via canoes until Europeans arrived on the islands and began ranching in the 1800s (Rick et al. 2014). Given the limited skunk remains in the archeological record and the similar morphology of island and mainland spotted skunks (Van Gelder 1959), it also remains a possibility that skunks were brought to the islands more recently, such as by people of European descent during or after the 17th century.
Investigating the evolutionary history of island spotted skunks is important from a conservation standpoint because they are considered a “Species of Special Concern” by the state of California due to a restricted distribution and small populations (Collins 1998). Additionally, skunk numbers appear to have recently declined on both islands (Bolas et al. 2020). Hence, uniting phylogenetic information on island spotted skunks with knowledge about their status is important for guiding ongoing management and conservation efforts (Cadotte and Davies 2010). Here, we estimate the timing of genetic divergence between spotted skunks on the mainland and islands by sequencing and analyzing the complete mitochondrial genome (mitogenome). We then apply our analysis of genetic divergence and patterns to explore the timing and methods of spotted skunk colonization of the Channel Islands.
Materials and Methods
Sample collection
We used DNA from 36 spotted skunks sampled from Santa Cruz (n = 22) and Santa Rosa (n = 14) islands obtained from tissues (blood, hair follicles, skin, or muscle) collected from island spotted skunks during two periods of trapping on both islands, 2001–2004 and 2017–2018. All trapping and sampling of live skunks followed guidelines for use of mammals in research (Sikes et al. 2016). Samples from 2001 to 2004 were collected under animal care and use protocols approved by the University of California, Davis and the University of California, Santa Barbara. Samples from 2017 to 2018 were collected by the National Park Service using approved study plans for annual monitoring of vertebrates. Additionally, we included DNA from museum samples (thin slice of toepad epidermis from specimens collected between 1919 and 1990) from the islands and mainland (California and Arizona, n = 34) that were either previously extracted (by Floyd et al. 2011, n = 33) or were part of a sample set of mainland spotted skunks for which we sought preliminary genome-wide sequencing data (n = 23) for a broader, ongoing study (Fig. 1; Supplementary Data SD1). Samples from mainland California included Santa Barbara County (n = 16), Los Angeles County (n = 8), Riverside County (n = 4), San Bernardino County (n = 2), and Kern County (n = 1). We used samples from Arizona (n = 3) as an outgroup representing the Arizona clade of western spotted skunks (Ferguson et al. 2017; Sonoran clade or S. leucoparia, McDonough et al. 2022). Hair samples were stored in envelopes in a dry area, and all blood, tissue, and extracted DNA samples were stored at −20°C or −80°C.
Localities for samples of 89 western spotted skunks (Spilogale gracilis) from Santa Cruz Island (n = 27, blue), Santa Rosa Island (n = 31, orange), and mainland California, United States (n = 31, black; inset). Samples were collected in three time intervals: 1919–1990, 2000–2004, and 2017–2018. For each locality, proportions indicate the number of successfully recovered mitogenomes (color fill in pie chart) over the number of samples sequenced. We also sequenced three samples from western spotted skunks in Arizona, not shown.
Laboratory methods
We extracted total genomic DNA from hair samples using a DNeasy Blood and Tissue kit (Qiagen Inc., Valencia, California) according to manufacturer’s instructions, with two 100-μl elutions of buffer AE to maximize DNA concentration and yield. Prior to DNA extraction, hair samples (5–10 hairs with follicles, 2 cm in length) were digested overnight at 56°C in proteinase-k and dithiothreitol. We extracted DNA from tissue using the Gentra Puregene Tissue kit (Qiagen Inc., Valencia, California). For museum samples, all extractions were conducted in a space dedicated to ancient DNA, using both Puregene and DNeasy Blood and Tissue kits. Some samples used DNA extracted for a previous study, described in Floyd et al. (2011). We sonicated DNA samples (except those from museum specimens) with a Covaris E220 machine (Covaris, Woburn, Massachusetts) to produce fragments averaging ~400 bp and prepared total genomic libraries using a NEBNext Ultra II DNA Library Prep kit for Illumina (E7645L), following manufacturer’s recommendations for amounts of DNA <50 ng but with reduced reactions to half the volume. We included negative controls during extractions and library preparation steps. We pooled libraries in equimolar ratios and sequenced in paired-end 150-bp mode on an Illumina HiSeq4000 platform (Illumina, San Diego, California) at the DNA Technologies and Expression Analysis Core Laboratory at the University of California, Davis Genome Center.
Alignment
We trimmed sequencing reads using default settings in ngsShoRT (Chen et al. 2013) and assembled them into circularized mitogenomes with Novoplasty v.2.7.2 (Dierckxsens et al. 2017), using a previously published Spilogale mitogenome as a starting seed and reference (GenBank accession number NC_010497.1). The reference genome was erroneously designated as S. putorius in GenBank, but comparison of the cytochrome b gene to other GenBank accessions unambiguously grouped it with S. gracilis, allowing incorporation as an additional sample in our alignment. For samples that did not produce circularized assemblies, we aligned reads to the reference mitogenome using BWA MEM v.0.7.17 (Li and Durbin 2009). After using SAMTools v.1.19 (Li et al. 2009) to remove PCR duplicates and only retain properly paired reads, we manually assembled aligned reads in Sequencher v.5.4.6 (Gene Codes Corporation, Ann Arbor, Michigan). In the final aligned data set, we removed insertions and deletions, and excluded the D-loop. We chose to exclude the D-loop because preliminary Sanger sequencing revealed little variation and the existence of large repetitive sections prevented its alignment in lower-coverage samples. We used this final alignment for all following analyses.
Network analysis
We constructed haplotype networks using a median-joining network (Bandelt et al. 1995; Bandelt et al. 1999) implemented in PopART v.1.7 (Leigh and Bryant 2015). We calculated rho statistics, the average number of mutations separating ancestral and descendent haplotypes in a clade (Forster et al. 1996; Macaulay et al. 2019) using Network v.5.0 (Bandelt et al. 1995). To scale rho to an estimate of the clade age in years, we applied the mutation rate estimated in our phylogenetic BEAST analysis (see below). Because Network treats uncalled or ambiguous bases, normally coded as “N,” as a distinct base, we imputed all such sites, most of which were invariant across our data set. We computed diversity statistics using DnaSP v.6 (Rozas et al. 2017).
BEAST analysis
We estimated divergence time between mainland and island mitogenomes in BEAST v.2.6 (Bouckaert et al. 2014). We used the Arizona clade of western spotted skunks as an outgroup and calibrated the root age using a splitting time estimated from a recently created Spilogale mitochondrial phylogeny (McDonough et al. 2022). We applied a lognormal prior to the root node (mean = 655,000 in real space, SD = 0.25), such that the 2.5 and 97.5 quantiles (389,000–1,040,000) captured the uncertainty indicated by the 95% highest posterior density (HPD) of the estimate. An alternative mitochondrial phylogeny for Spilogale estimated the divergence between the Arizona and western clades to be more ancient (1.36 Ma, 0.80–2.11 95% HPD; Ferguson et al. 2017), so we also estimated divergence times using this date as an alternate prior on the root node. For all BEAST analyses, we used a strict clock model appropriate for intraspecific trees, applying a lognormal prior with a broad SD (mean = −18.01521564, SD = 1.5) that spans previously published carnivore mitochondrial rates (quantiles for the clock rate prior in real space [substitutions/site/year]: median 1.50 × 10−8, 2.5% 7.93 × 10−10, 97.5% 2.48 × 10−7; Goddard et al. 2015). We partitioned sequence data by gene codon position and assigned all other positions (tRNA, rRNA, origin of the light strand) to a fourth partition. The ND6 gene, which is transcribed in the reverse direction, was grouped with the fourth partition. We linked trees across partitions and unlinked substitution models, using bModelTest v1.21 to average out uncertainty over possible models of evolution (Bouckaert and Drummond 2017). To evaluate the sensitivity of estimates to demographic assumptions, we considered two tree priors, the constant coalescent and the coalescent Bayesian skyline (Drummond et al. 2005). In both cases, Markov chains were run for 100 million generations with 10% discarded as burn-in, at which point we confirmed convergence and effective sampling sizes >1,000 for all parameters in Tracer v.1.7 (Rambaut et al. 2018). We then used Multithreaded Nested Sampling (Russel et al. 2019) to compare marginal likelihoods between tree priors, using five threads, 10 particles, and a subchain length of 200,000. The log Bayes factor (Kass and Rafferty 1995) was estimated as 7.3 in strong favor of the Bayesian skyline. Results therefore correspond to the model using the skyline tree prior, constructed from a maximum clade credibility tree based on median heights in TreeAnnotator v.2.6 (Bouckaert et al. 2014) and visualized using the ggtree R package (Yu et al. 2017).
IMa2 analysis
In addition to rho and BEAST estimates, we also estimated splitting times between the skunk populations on the two islands using coalescent MCMC simulations in IMa2 (Hey and Nielson 2007). We used an HKY substitution model and, upon finding no support for postdivergence gene flow between islands, only considered the no-migration model. We conducted five independent runs with different random seeds for 50 million steps following 10 million steps of burn-in, based on 40 independent chains that followed a geometric heating model (-hfg -hn40 -ha0.999 -hb0.3). We compared effective sample sizes and consistency of results across runs with different random number seeds to evaluate convergence. To convert scaled parameters to absolute years, we applied the median mutation rate estimated in the BEAST analysis.
Results
We obtained complete mitogenome sequences from 55 of 92 samples. Of samples collected before 2000, 23% successfully sequenced, and 87% of samples collected in or after 2000 were successful; all successful island samples were modern (Supplementary Data SD1). Combining our mitogenomes with the previously published reference sequence yielded a total data set of 46 island spotted skunks (15 from Santa Cruz, 31 from Santa Rosa), eight western spotted skunks from mainland California (seven from California, one from GenBank of unknown origin), and two western spotted skunks from Arizona (Table 1; Fig. 1; Supplementary Data SD1). After excluding the D-loop and removing insertions, the final data set consisted of 15,441 bp of aligned sequences based on a mean read depth of 55× (range 8–102×; Supplementary Data SD1). All sequences were deposited in GenBank (accession numbers OK561938).
Number of samples yielding mitogenomes (n), along with number (No.) of distinct haplotypes, haplotype diversity, and nucleotide diversity from island spotted skunks (Spilogale gracilis amphiala) collected on Santa Cruz and Santa Rosa islands, and western spotted skunks (S. gracilis) from mainland California and Arizona, United States.
| Locality | n | No. haplotypes | Haplotype diversity (SD) | Nucleotide diversity (SD) × 10−3 |
|---|---|---|---|---|
| Santa Cruz Island | 15 | 5 | 0.733 (0.089) | 0.10 (0.02) |
| Santa Rosa Island | 31 | 4 | 0.385 (0.099) | 0.06 (0.02) |
| Both Islands | 46 | 9 | 0.697 (0.064) | 0.17 (0.02) |
| California mainlanda | 8 | 8 | 1.000 (0.063) | 2.52 (0.29) |
| Arizona | 2 | 2 | 1.000 (0.5000) | 1.30 (0.65) |
| Total | 56 | 19 | 0.796 (0.050) | 2.83 (1.33) |
| Locality | n | No. haplotypes | Haplotype diversity (SD) | Nucleotide diversity (SD) × 10−3 |
|---|---|---|---|---|
| Santa Cruz Island | 15 | 5 | 0.733 (0.089) | 0.10 (0.02) |
| Santa Rosa Island | 31 | 4 | 0.385 (0.099) | 0.06 (0.02) |
| Both Islands | 46 | 9 | 0.697 (0.064) | 0.17 (0.02) |
| California mainlanda | 8 | 8 | 1.000 (0.063) | 2.52 (0.29) |
| Arizona | 2 | 2 | 1.000 (0.5000) | 1.30 (0.65) |
| Total | 56 | 19 | 0.796 (0.050) | 2.83 (1.33) |
aIncludes one sample from GenBank of unknown origin (accession no. NC_010497.1).
Number of samples yielding mitogenomes (n), along with number (No.) of distinct haplotypes, haplotype diversity, and nucleotide diversity from island spotted skunks (Spilogale gracilis amphiala) collected on Santa Cruz and Santa Rosa islands, and western spotted skunks (S. gracilis) from mainland California and Arizona, United States.
| Locality | n | No. haplotypes | Haplotype diversity (SD) | Nucleotide diversity (SD) × 10−3 |
|---|---|---|---|---|
| Santa Cruz Island | 15 | 5 | 0.733 (0.089) | 0.10 (0.02) |
| Santa Rosa Island | 31 | 4 | 0.385 (0.099) | 0.06 (0.02) |
| Both Islands | 46 | 9 | 0.697 (0.064) | 0.17 (0.02) |
| California mainlanda | 8 | 8 | 1.000 (0.063) | 2.52 (0.29) |
| Arizona | 2 | 2 | 1.000 (0.5000) | 1.30 (0.65) |
| Total | 56 | 19 | 0.796 (0.050) | 2.83 (1.33) |
| Locality | n | No. haplotypes | Haplotype diversity (SD) | Nucleotide diversity (SD) × 10−3 |
|---|---|---|---|---|
| Santa Cruz Island | 15 | 5 | 0.733 (0.089) | 0.10 (0.02) |
| Santa Rosa Island | 31 | 4 | 0.385 (0.099) | 0.06 (0.02) |
| Both Islands | 46 | 9 | 0.697 (0.064) | 0.17 (0.02) |
| California mainlanda | 8 | 8 | 1.000 (0.063) | 2.52 (0.29) |
| Arizona | 2 | 2 | 1.000 (0.5000) | 1.30 (0.65) |
| Total | 56 | 19 | 0.796 (0.050) | 2.83 (1.33) |
aIncludes one sample from GenBank of unknown origin (accession no. NC_010497.1).
Relationship of island to mainland populations
The Bayesian phylogenetic tree grouped island spotted skunks as a monophyletic clade nested within the clade of western spotted skunks from mainland California (Fig. 2). The network analysis also supported monophyly of island skunks, depicting a single basal haplotype sampled on Santa Rosa from which all other closely related island haplotypes were derived (Fig. 3). This basal island haplotype was 37 mutational steps away from the nearest sampled mainland haplotype, but our sampling of mainland diversity was cursory and unlikely to represent the high extant mitochondrial diversity on the mainland suggested by the observed haplotype and nucleotide diversity (Table 1). Consequently, we could not assume that the divergence between mainland and island samples accurately reflected the timing of island colonization. However, the monophyly and tight clustering of island haplotypes (i.e., no separation greater than two mutations) relative to average distance among haplotypes on the mainland strongly suggested that all skunks sampled from both islands derive from a single maternal ancestor (Fig. 3). Based on our calibrated Bayesian tree using the root prior of the McDonough et al. (2022) phylogeny, we estimated the time to most recent common ancestor (TMRCA) for the island matriline at 5,690 years (95% HPD: 1,910–11,320; Table 2; Fig. 2). When we applied the BEAST-estimated mutation rate to estimates of the TMRCA based on rho and coalescent IMa2 analyses, point estimates were broadly similar (within 1,200 years of the BEAST estimate) with overlapping uncertainty intervals (Table 2; Fig. 4A). When applying the root prior from Ferguson et al. (2017), the mutation rate halved and the timing of divergence approximately doubled in all three analyses (Table 2; Figs. 2 and 4A; Supplementary Data SD3).
Point estimates of the time to most recent common ancestor (TMRCA) for island spotted skunks (Spilogale gracilis amphiala) based on mitogenome analyses. Mutation rates and TMRCA were estimated in BEAST v2.6 using priors for the root node based on the divergence time of the Arizona and western clades of western spotted skunks (S. gracilis), estimated by two recent spotted skunk phylogenies. BEAST-estimated mutation rates were then applied to IMa2 and rho analyses. Our calibration based on the McDonough phylogeny (McDonough et al. 2022) used a root node of 655,00 years, which resulted in a BEAST-estimated mutation rate of 2.98 × 10−8 substitutions per site per year (Fig. 2). Our calibration based on Ferguson et al. (2017) used a root node of 1.36 Ma, which resulted in a BEAST-estimated mutation rate of 1.46 × 10−8 substitutions per site per year (Supplementary Data SD3).
| Analysis | Point estimate (error interval) | TMRCA—McDonough et al. (2022) calibration | TMRCA—Ferguson et al. (2017) calibration |
|---|---|---|---|
| Beast | Median (95% HPD) | 5,690 (1,910–11,320) | 11,560 (3,980–23,160) |
| IMa2 | High point (95% HPD) | 4,900 (2,500–10,340) | 9,980 (5,100–21,080) |
| Network rho | Mean (1.96 SD) | 4,450 (640–8,260) | 8,770 (1,260–16,280) |
| Analysis | Point estimate (error interval) | TMRCA—McDonough et al. (2022) calibration | TMRCA—Ferguson et al. (2017) calibration |
|---|---|---|---|
| Beast | Median (95% HPD) | 5,690 (1,910–11,320) | 11,560 (3,980–23,160) |
| IMa2 | High point (95% HPD) | 4,900 (2,500–10,340) | 9,980 (5,100–21,080) |
| Network rho | Mean (1.96 SD) | 4,450 (640–8,260) | 8,770 (1,260–16,280) |
HPD = highest posterior density.
Point estimates of the time to most recent common ancestor (TMRCA) for island spotted skunks (Spilogale gracilis amphiala) based on mitogenome analyses. Mutation rates and TMRCA were estimated in BEAST v2.6 using priors for the root node based on the divergence time of the Arizona and western clades of western spotted skunks (S. gracilis), estimated by two recent spotted skunk phylogenies. BEAST-estimated mutation rates were then applied to IMa2 and rho analyses. Our calibration based on the McDonough phylogeny (McDonough et al. 2022) used a root node of 655,00 years, which resulted in a BEAST-estimated mutation rate of 2.98 × 10−8 substitutions per site per year (Fig. 2). Our calibration based on Ferguson et al. (2017) used a root node of 1.36 Ma, which resulted in a BEAST-estimated mutation rate of 1.46 × 10−8 substitutions per site per year (Supplementary Data SD3).
| Analysis | Point estimate (error interval) | TMRCA—McDonough et al. (2022) calibration | TMRCA—Ferguson et al. (2017) calibration |
|---|---|---|---|
| Beast | Median (95% HPD) | 5,690 (1,910–11,320) | 11,560 (3,980–23,160) |
| IMa2 | High point (95% HPD) | 4,900 (2,500–10,340) | 9,980 (5,100–21,080) |
| Network rho | Mean (1.96 SD) | 4,450 (640–8,260) | 8,770 (1,260–16,280) |
| Analysis | Point estimate (error interval) | TMRCA—McDonough et al. (2022) calibration | TMRCA—Ferguson et al. (2017) calibration |
|---|---|---|---|
| Beast | Median (95% HPD) | 5,690 (1,910–11,320) | 11,560 (3,980–23,160) |
| IMa2 | High point (95% HPD) | 4,900 (2,500–10,340) | 9,980 (5,100–21,080) |
| Network rho | Mean (1.96 SD) | 4,450 (640–8,260) | 8,770 (1,260–16,280) |
HPD = highest posterior density.
Bayesian phylogeny of island (Spilogale gracilis amphiala) and mainland California western (S. gracilis) spotted skunks based on 56 mitogenomes (15,441 bp). The median and 95% highest posterior density (HPD) intervals for the time to the most recent common ancestor (TMRCA) for island spotted skunks was estimated in BEAST v.2.6, calibrated using the most recent estimate of divergence time between the western and Arizona clades of western spotted skunks (McDonough et al. 2022; Table 2). When the estimate from Ferguson et al. (2017) was used to calibrate the divergence time between the western and Arizona clades, the root node was 1.36 Ma and the TMRCA for island skunks was approximately twice as ancient at 11,560 years (95% HPD 3,980–23,160; Table 2 and Supplementary Data SD3). Filled circles indicate nodes strongly supported by Bayesian posterior probabilities (BPP = 1.0).
Median-joining network of island (n = 46, Spilogale gracilis amphiala) and mainland California western (n = 8, S. gracilis) spotted skunk mitogenomes (15,441 bp), rooted to the Arizona clade (n = 2). Circles represent haplotypes; sample size is indicated by the diameter and numbers inside the circle (for n > 2), and arrow indicates the direction of rooting for island samples, as inferred from the phylogenetic tree (see Fig. 2). Branch lengths are proportional to numbers of substitutions between haplotypes, indicated by perpendicular lines, except where condensed to save space, in which case numerals indicate numbers of substitutions. One haplotype from GenBank (*) was from an unknown location (accession no. NC_010497.1).
Coalescent analysis of mitogenomes (15,441 bp) of 46 island spotted skunks (Spilogale gracilis amphiala) in IMa2, showing probability (posterior density) distributions corresponding to (A) estimates of the time to most recent common ancestor (TMRCA) and population splitting time (divergence) for Santa Cruz and Santa Rosa islands, assuming an average mutation rate of 2.98 × 10−8 substitutions/site/year (McDonough et al. 2022 calibration) and 1.46 × 10−8 substitutions/site/year (Ferguson et al. 2017 calibration), and (B) unscaled estimates of historical effective population sizes (ϴ).
Relationship between island populations
None of the nine island haplotypes observed (Santa Cruz n = 5; Santa Rosa n = 4) were shared between islands (Fig. 3). This clustering of haplotypes of each island into distinct haplogroups indicated a lack of gene flow between islands following divergence. Coalescent IMa2 analysis estimated that divergence between the two island populations occurred approximately 1,500 years after the TMRCA, but broadly overlapping 95% HPD intervals suggest that the difference between estimates of splitting time and TMRCA was not significant (Fig. 4A; Supplementary Data SD2.1). Similarly, rho estimates of the TMRCA for Santa Rosa, Santa Cruz, and both islands together were essentially equivalent when the SD was accounted for (Supplementary Data SD2.2). Together, the sorting of haplotypes into nearly reciprocally monophyletic clades between islands and molecular dating suggested the skunk populations on Santa Rosa and Santa Cruz diverged soon after the time of colonization and have remained on independent trajectories for most of the Holocene.
Diversity in island populations
Haplotype and nucleotide diversity estimates were substantially lower on the islands compared to mainland California (Table 1). Between the islands, the Santa Cruz population possessed higher levels of haplotype and nucleotide diversity relative to the Santa Rosa population. According to the IMa2 analysis, matrilineal effective population size was slightly higher on Santa Cruz than on Santa Rosa, although error ranges overlapped (Fig. 4B; Supplementary Data SD2.1). There were no significant differences in diversity estimates between samples collected during 2001–2004 and 2017–2018 on either island.
Discussion
Timing of island colonization
Evidence of monophyly for island spotted skunks is consistent with a single colonization event from the mainland. Given the high matrilineal diversity in the mainland population, it is exceedingly unlikely that we would have observed monophyly and such tight clustering of haplotypes in the island populations if they were derived from more than one founding matriline from the mainland. Additionally, there were no haplotypes shared between skunk populations on the two islands and the haplotypes on each island formed distinct clades, albeit with Santa Cruz nested within Santa Rosa. This pattern suggests that the two island populations became isolated from one another soon after initial establishment of the ancestral population. The most parsimonious explanation is that skunks colonized in one event, shortly before the breakup of Santarosae Island. For example, a single pregnant female could have been the sole founder for all island skunks. Alternatively, the founder population could have included multiple individuals and haplotypes from the mainland, followed by stochastic loss of all but one matriline before populations diverged. Although less likely, spotted skunks could have colonized the islands after the breakup of Santarosae. The most plausible such scenario would be a stepping-stone model, whereby skunks arrived first on Santa Rosa Island and then, soon after, founded a population on Santa Cruz Island. Regardless, the near reciprocal monophyly observed between matrilines of the two islands indicates that these populations have existed in isolation for enough time for their private mutations to accumulate. The estimated range of divergence time between mainland and island spotted skunks from our molecular dating analysis in BEAST using the calibration from the McDonough et al. (2022) phylogeny was 1,910–11,320 years ago, while using the calibration from Ferguson et al. (2017) doubled the time range. Both time frames suggest a Holocene colonization of the Northern Channel Islands by spotted skunks that encompasses the time of the separation of Santa Cruz and Santa Rosa islands (9,400–9,700 years ago; Reeder-Myers et al. 2015), and postdates human arrival but predates the arrival of European settlers (13,000 years ago and 250 years ago, respectively; Rick et al. 2014). Given the support for a colonization prior to the breakup of Santarosae indicated by the network topology, we suggest that the most plausible timing of arrival for spotted skunks to the islands was during the early Holocene, which is consistent with the midpoint of the Ferguson calibration and within the range of the McDonough calibration.
Our findings and proposed time frame are consistent with both rafting on floating debris and Native American transport as explanations for how spotted skunks colonized the Channel Islands but refute translocation of skunks by European settlers. If skunks arrived on the Channel Islands by rafting, the most plausible timing would be 9,000–20,000 years ago, when the northern Channel Islands were grouped together as Santarosae and were only 7–12 km from mainland California. Rafting seems feasible for small-bodied animals (King 1962; Ali and Vences 2019), including small carnivores (i.e., Euplerida [1.5–1.8 kg], Yoder et al. 2003; Paradoxurus hermaphroditus [2–5 kg], Ali and Vences 2019), and may be facilitated by storms (Censky et al. 1998). Coastal California can experience winter rain and flooding capable of propelling debris from the mainland to the islands. For example, in 1955, a live black-tailed jackrabbit (Lepus californicus; 1.3–3 kg) was recovered from floating debris 63 km from mainland Southern California, beyond Santa Catalina Island; no jackrabbits inhabited the Channel Islands at that time (Prescott 1959). Western spotted skunks (0.5–0.9 kg) mate in the fall and experience delayed implantation before giving birth in the spring (Mead 1968), so it is plausible that a single pregnant female could have been swept to sea by a winter storm and survived to reach the northern Channel Islands without losing the pregnancy.
Similar to rafting, a pregnant skunk presumably would be able to survive an unintentional translocation as a stowaway in a Native American boat, although transport may have also been intentional. Some indigenous groups in North America appear to have valued skunk species for religious purposes or as food (Wenner and Johnson 1980), and striped skunks (Mephitis mephitis) have been recovered from archeological sites along the southern California coast (Landberg 1965). However, archeological evidence of spotted skunks on the Channel Islands is limited and radiocarbon dating is needed for the few skunk remains that have been found: spotted skunk bones were recovered from one midden site each on San Miguel Island (strata dated to be roughly 700–3,000 years old; Walker 1980; Erlandson et al. 1996), Santa Rosa Island (strata dated at 300–400 years old; Rick 2011), and Santa Cruz Island (strata dated 170–650 years old; Rick 2013). Additionally, oral histories describing associations between island spotted skunks and Native Americans are sparse enough that the relationship between island skunks and Native Americans on the Channel Islands is unknown. However, an account of a skunk dance performed by Native Americans living on the mainland of California to the east of the Channel Islands suggests familiarity with skunks (Rick 2013). In contrast, island foxes have been found in ceremonial burial sites across the northern and southern islands, and the ethnographic record suggests foxes were used as pets, food, and in religious ceremonies by Native Americans on the Channel Islands (Collins 1991a; Rick et al. 2019). Additionally, Native Americans likely brought island foxes to the southern Channel Islands, presumably because of their cultural or dietary value (Collins 1991b; Wayne et al. 1991; Rick et al. 2009; Goddard et al. 2015; Hofman et al. 2015, 2016). If Native Americans valued skunks and deliberately brought them to the Channel Islands, it seems likely that such an introduction would have involved multiple individuals, rather than the single matrilineal ancestor, and occurred on multiple islands.
Among the four land mammals native to Santa Cruz and Santa Rosa islands, there is evidence for colonization both by rafting and human transport over a wide time frame. Island foxes derive from a single matrilineal founder associated with the northern islands approximately 7,300–19,700 years ago (Goddard et al. 2015), although it is unknown whether foxes arrived on the northern Channel Islands via rafting or human introduction. Deer mice (Peromyscus maniculatus), of which there is an endemic subspecies found on each of the eight islands, appear to have colonized the Channel Islands multiple times over the last 500,000 years (Ashley and Wills 1987; Shirazi et al. 2018). Several of these colonization events were probably inadvertent translocation by humans (Rick 2013; Shirazi et al. 2018), but at least one colonization would have been independent of humans and was likely a rafting event prior to the breakup of Santarosae, as evidenced by a single shared haplotype among deer mice on the northern islands (Ashley and Wills 1987). Western harvest mice (Reithrodontomys megalotis) on the northern islands are weakly differentiated from western harvest mice on the mainland based on both morphological and genetic analysis, suggesting a late-Holocene introduction, perhaps as inadvertent stowaways in a Native American canoe (Ashley 1989; Collins and George 1990). Taken together, the pattern and timing of colonization by spotted skunks to the northern Channel Islands is very similar to that of gray foxes but differs from that of the two rodents.
Diversity of island spotted skunk populations
As expected, island spotted skunk populations on both islands evidenced substantially lower mitogenomic haplotype and nucleotide diversity than found in the mainland population. Microsatellite data exhibited the same pattern, with lower heterozygosity in the two island populations (He = 0.58) relative to the mainland (He = 0.78; Floyd et al. 2011). Low genetic diversity is a common attribute of island populations due to a combination of founder effects and inherent constraints on population size (Frankham 1997; Whittaker and Fenandez-Palacios 2007; Furlan et al. 2012). The isolation of island populations also makes them susceptible to severe population fluctuations, further driving down genetic diversity and potentially leading to genetic drift, as seen with island foxes (Habel and Zachos 2013; Funk et al. 2016). We found that skunks on Santa Cruz Island had slightly higher diversity and larger effective population sizes than did skunks on Santa Rosa Island, but sample sizes were insufficient for robust comparison. Analysis of the nuclear genome could clarify the demographic histories of skunk populations on each island and help evaluate the degree to which inbreeding and genetic load may be impacting population fitness (e.g., Robinson et al. 2016, 2018).
Conservation implications
Evolutionary histories and genetic lineages should be carefully considered in defining conservation priorities (Cadotte and Davies 2010; Rosauer and Jetz 2015) and assigning taxonomic designations (Crandall et al. 2000; Fraser and Bernatchez 2001; Sacks et al. 2010). For island spotted skunks, their biogeographic history, as well as the potential similarity with the biogeographic histories of other endemic mammals on the Channel Islands, may be informative for guiding management and conservation planning. Our results demonstrate that spotted skunks arrived on the northern Channel Islands in a similar time frame as gray foxes. Additionally, rapid evolution is a common pattern in island mammals (Foster 1964; Millien 2006) that appears to have occurred in both island foxes and island deer mice, and thus may hold true for island spotted skunks as well. Island foxes are characterized by island dwarfism, which may have evolved over only a few thousand years (Hofman et al. 2015; Funk et al. 2016), and the subspecies of deer mice that are endemic to three islands (Anacapa, Santa Cruz, and Santa Barbara islands) evidenced rapid morphological change over a 90-year period, likely genetically based (Pergams and Ashley 1999). Further research is needed to assess morphological or cryptic functional differences that may have arisen between skunks on the islands and mainland, and additional evaluation of mitochondrial and nuclear DNA will continue to add information about the evolutionary history of island spotted skunks. Nonetheless, we found evidence of several thousand years of separation between island and mainland populations of spotted skunks and between the populations on the two islands. Moreover, each island supports an endemic lineage of spotted skunks, with genetic diversity that is spatially restricted and has a distinct evolutionary history. Taken together, evidence of biogeographic and phylogenetic endemism and isolation over a long time period for spotted skunks on both Santa Cruz and Santa Rosa islands suggests that reevaluation of the taxonomic status of island spotted skunks may be warranted. Given the evidence of long-term endemism and current rarity of island spotted skunks, island spotted skunks on each island should be managed as distinct evolutionarily significant units and prioritized for conservation.
Supplementary Data
Supplementary data are available at Journal of Mammalogy online.
Supplementary Data SD1.—Sample and sequencing information for DNA of 92 western spotted skunks (Spilogale gracilis), including museum catalogue number (NA if not a museum specimen), collection site and region, latitude and longitude where available (or the closest known location), year collected, and sample type. Abbreviations in museum catalogue numbers refer to the Natural History Museum of Los Angeles County (LACM); Santa Barbara Museum of Natural History (SBCM); and Dickey Collection of Birds and Mammals, University of California at Los Angeles (UCLA). Table also includes the mean sequence quality score, total number of sequenced reads, median insert size of total sequenced reads, the number and percentage of reads that aligned to the nuclear Spilogale gracilis reference genome (GenBank accession number GCA_004023965.1), the number of reads that aligned to the mitochondrial Spilogale reference genome (GenBank accession number NC_010497.1) according to both BWA and Novoplasty methods, the percentage of BWA aligned mitochondrial reads that were PCR duplicates, whether a full mitogenome was successfully recovered, method of final assembly, mean depth of assembled mitogenome, and the GenBank accession number.
Supplementary Data SD2.—Results from analysis of mitogenomes for island (n = 46; Spilogale gracilis amphiala) and mainland western spotted skunks (n = 10; S. gracilis), from Santa Cruz and Santa Rosa islands, and mainland California and Arizona, United States, that examine relationships between island populations. Table S2.1.—Summary of maternal effective population sizes and divergence estimates from IMa2 for island spotted skunk mitogenomes (15,441 bp), showing estimates scaled to coalescent units and, in the case of time estimates, scaled to absolute years using an average mutation rate of 2.98 × 10−8 substitutions/site/year (McDonough et al. 2022 calibration) and 1.46 × 10−8 substitutions/site/year (Ferguson et al. 2017 calibration). Table S2.2.—Rho statistics of average (and standard deviation, SD) numbers of mutations separating descendant haplotypes from ancestral haplotypes and corresponding estimates of time to the most recent common ancestor (TMRCA), assuming an average mutation rate of 2.98 × 10−8 substitutions/site/year (McDonough et al. 2022 calibration) and 1.46 × 10−8 substitutions/site/year (Ferguson et al. 2017 calibration).
Supplementary Data SD3.—Bayesian phylogeny of island (Spilogale gracilis amphiala) and mainland California western (S. gracilis) spotted skunks based on 56 mitogenomes (15,441 bp). The median and 95% HPD intervals for the time to the most recent common ancestor (TMRCA) for island spotted skunks was estimated in BEAST v.2.6, calibrated using an estimate of divergence time between the western and Arizona clades of western spotted skunks from Ferguson et al. (2017; Table 2). Filled circles indicate nodes strongly supported by Bayesian posterior probabilities (BPP = 1.0). See Fig. 2 for a phylogeny of island spotted skunks calibrated using the most recent estimate of divergence time between the western and Arizona clades of western spotted skunks (McDonough et al. 2022 calibration).
Acknowledgments
We thank curators from the Natural History Museum of Los Angeles County; the Santa Barbara Museum of Natural History; and the Dickey Collection of Birds and Mammals, University of California at Los Angeles for generously allowing us to sample toe pads from their spotted skunk specimens. We also thank the National Park Service for assistance from staff and volunteers with logistics and field work, and specifically P. Powers, J. Schamel, S. Baker, M. Marshall, L. McMahon, L. Lee, and J. Fitzgerald. We thank members of the Sacks Lab for their assistance with lab work.
CONFLICT OF INTEREST
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
FUNDING
This work was supported by grants from the National Park Service Southern California Research Learning Center; the University of California, Davis, Graduate Group in Ecology Student Endowment and Jastro Award; and the University of California, Davis, Department of Wildlife, Fish, and Conservation Biology Howard Wildlife Management Fund. ECB was supported by a University of California, Davis Graduate Group in Ecology Fellowship and a National Science Foundation Graduate Research Fellowship.
Literature Cited
Author notes
These authors contributed equally.
