Supplementary Information

Appendix S1 – Moment-based estimation of genetic effective population size and approximate 95% confidence intervals.

In principle, the variance effective population size can be measured via comparison of allele frequencies at 2 points in time (Waples 1989). However, estimators available for such calculations require large sample sizes of individuals in each time period because they sum covariance terms across individual alleles, necessitating the subtraction of a sampling variance component dependent on the numbers of individuals in the 2 temporal samples but unaffected by the number of loci used (e.g., Waples 1989). If sample size is small, this variance component will be large relative to true changes in allele frequencies. Not only does this result in imprecise estimates, but Ne estimates will be biased upwards in small populations. Therefore, we took advantage of the fact that unbiased estimates of He can be calculated in each time period with sampling variance inversely proportional to the combined number of individuals and loci. This method also can be used for mitochondrial-based estimates of gene diversity. Point estimates were calculated from the decline in He as follows (Hedrick 2000):

Ne = 1 / {2*[1-( He(Mod)/He(Hist))1/t]},

Where Ne is the variance effective population size, He(Mod) is the modern expected (under Hardy-Weinberg equilibrium) heterozygosity (or gene diversity), He(Hist) is the historical heterozygosity (or gene diversity), and t is the number of generations between modern and historical samples, which would be 75, conservatively assuming a 1-year generation time, or 37.5 assuming a 2-year generation time. Although red foxes have overlapping generations and the formal proof of this equation assumes nonoverlapping generations, the large number of generations separating historical and modern samples likely rendered the violation of this strict assumption negligible relative to uncertainty in the generation time.

Upper 95% confidence limits were approximated using the heterozygosity (gene diversity) estimates and their estimated standard deviations (SD) in each time period as follows:

Upper CL Ne = 1 / {2*[1-( [He(Mod)+1SD]/[He(Hist) -1SD])1/t]}.

This confidence limit estimate for Ne assumes a 65% probability that the true modern (and historical) heterozygosity is within +1 SD of the point estimate, which implies a 17.5% probability that the true value is >1 SD greater than the point estimate and an equal (17.5%) probability that the true value is <1SD less than the point estimate. The upper confidence limit for Ne then is based on the co-occurrence of an underestimated modern heterozygosity and an overestimated historical heterozygosity. If the magnitude of the over- and underestimates is 1 SD, this implies a co-occurrence probability of 0.1752 = 0.03 (i.e., α/2), which corresponds to a 94% (1 - α) confidence limit.


Appendix S2 – Microsatellite Hardy-Weinberg and linkage equilibria

Heterozygote deficiency or Hardy-Weinberg disequilibrium (HWD) can indicate allelic dropout (random or null alleles), which represent unwanted artifacts, or from inbreeding, which is a characteristic of the population. We performed analyses of variance (ANOVAs) to differentiate between these causes, using locus- and population-specific estimates of FIS (a measure of heterozygote deficit) as a response variable. Initially, we included all loci; a significant main effect of locus was followed-up with post-hoc tests (Fisher least significant difference [LSD]) to identify a high-HWD grouping of loci, likely to reflect higher allelic dropout or null alleles. Next, we repeated the analysis excluding the high-HWD loci to assess differences among populations, likely to correspond to differences in inbreeding or population substructure.

Based on a 1-way ANOVA (F13, 84 = 4.27, P <0.001) and post-hoc Fisher’s LSD test on locus-specific FIS across populations, 3 significantly distinct clusters of loci were identified: (1) little or no heterozygote deficiency, average FIS = 0.12 (+SD = 0.08; AHT 133, AHT140, FH2054, C01.424, CPH2b, FH2010, FH2289, FH2328, FH2004), (2) moderate heterozygote deficiency, average FIS = 0.28 (+0.03;FH2380, FH2457), and (3) high heterozygote deficiency, average FIS = 0.39 (+0.06; C08.618, FH2088, FH2001). Failure to amplify any alleles in a small number of individuals suggested that 2 of these loci (FH2088, FH2001) had null alleles. To assess population-specific differences in FIS with minimal influence of allelic dropout, we performed a 1-way ANOVA excluding loci in the high-heterozygote-deficiency cluster. This comparison was significant (F6, 70 = 2.26, P = 0.047), indicating differences among populations in inbreeding (Table 3 in text). Using these same loci, the heterozygote deficit was slightly higher in the historical samples ( = 0.28) than the modern samples ( = 0.15), although the difference was marginally significant (t4 = 2.52, P = 0.07).

We used Genepop on the web (Raymond & Rousset 1995) to assess linkage disequilibrium (LD) for each locus pair in each population and tested the frequency distribution of locus-pairs in terms of numbers of populations exhibiting significant LD against Poison expectations using a Chi-square Goodness-of-fit test (Zar 1999; Sacks et al. 2004). This approach allowed us to identify particular locus-pairs with an unexpected degree of LD, possibly due to physical linkage. To test for population differences in LD, we conducted a Chi-square test of independence in frequency of significant LD (uncorrected to maximize statistical power to detect differences) among populations (Zar 1999).

The occurrence of significant (uncorrected) linkage disequilibrium (LD) associated with particular pairs of loci varied from 0-4 of the 7 populations and did not differ significantly from random expectations (= 1.25, P = 0.87), suggesting that no pair of loci was especially prone to LD (e.g., as would be expected if physically linked). However, the proportion of locus-pairs with significant LD varied among populations (= 33.0, P < 0.001), with the Sacramento Valley exhibiting the greatest LD (Table 3 in text).

Figure S1. Temporal distribution of red fox binned by decade (dates indicate beginning) (A) cytochrome-b sequences (n = 210), (B) D-loop sequences (n = 206), and (C) microsatellite genotypes (n = 211) used in this study, shaded according to where sequences were first published (A,B). Shadings correspond to the following sources: present study (black); Aubry et al. 2009 (grey); and Perrine et al. 2007 (white).

Figure S2. Relationship between elevation and latitude of red foxes in the Montane group, illustrating distinctive distribution of Sacramento Valley red foxes relative to montane populations (Northern Cascades, Southern Cascades, Sierra Nevada, Rocky Mountains).

Figure S3. FST and 95% CI estimated from simulations for 2 equal-sized populations, Ne = 100 (red) and Ne = 30 (blue), diverging without gene flow from 10 coalescent simulations conducted in Easypop (Balloux 2001). Parameters of the simulation were as follows: u = 0.0006, 14 loci, 10 alleles each, n = 20 individuals sampled. Shaded area corresponds to uncertainty in generation time and indicates that observed change in FST corresponds to 30 < Ne
Table S1. List of Vulpes vulpes specimens used in mitochrondrial sequencing analyses, along with additional samples used to augment microsatellite samples. Accessioned samples were from the following institutions: National Museum of Natural History, Washington, D.C. (USNM); Museum of Vertebrate Zoology, University of California at Berkeley (MVZ); University of Washington Burke Museum, Seattle (UWBM); Museum of Wildlife and Fisheries Biology, University of California at Davis (WFB); California State University, Chico (CSU); Los Angeles County Museum of Natural History (LA). All other abbreviations on specimen identifiers (ID) indicate modern samples not currently accessioned in a museum collection. Period indicates that samples were collected in either historical (H; 1880-1950) or modern (M; 1951-2008) time periods. Priming regions are described in text; names of major clades or subclades follow Aubry et al. (2009) except that the “O-24 subclade” is used here to indicate a subclade stemming from composite haplotype O-24 nested within the Mountain subclade.

Specimen ID / Population6 / Period / Year collected / Mountain subspecies6 / Major clade/subclade / Cyt-b / D-loop / Composite haplotype4
USNM-154153 / Eurasia / H / Holarctic / U6 / 57 / U6-57
USNM-175597 / Eurasia / H / --5
USNM-175598 / Eurasia / H / 1912 / Holarctic / U / 48 / U-48
USNM-200665 / Eurasia / H / --5
USNM-200667 / Eurasia / H / Holarctic / U4 / 55 / U4-55
USNM-200668 / Eurasia / H / 1915 / Holarctic / U4 / 55 / U4-55
USNM-200669 / Eurasia / H / 1915 / Holarctic / U / 70 / U-70
USNM-200671 / Eurasia / H / 1914 / Holarctic / U / 70 / U-70
USNM-200676 / Eurasia / H / 1915 / Holarctic / U / 70 / U-70
USNM-200677 / Eurasia / H / 1915 / Holarctic / U4 / 54 / U4-54
USNM-200681 / Eurasia / H / 1915 / Holarctic / U / 70 / U-70
USNM-200682 / Eurasia / H / 1915 / Holarctic / U4 / 55 / U4-55
USNM-200684 / Eurasia / H / 1915 / Holarctic / U / 70 / U-70
USNM-200684 / Eurasia / H / --5
USNM-200686 / Eurasia / H / 1914 / Holarctic / U / 69 / U-69
USNM-200677 / Eurasia / H / --5
USNM-319220 / Eurasia / M / 1960 / Holarctic / U4 / 56 / U4-56
M5857 / Eurasia / M / --5
M5858 / Eurasia / M / --5
M5859 / Eurasia / M / --5
M5860 / Eurasia / M / --5
M5861 / Eurasia / M / --5
S08-0428 / Eurasia / M / --5
S08-0433 / Eurasia / M / --5
S08-0435 / Eurasia / M / --5
S08-0436 / Eurasia / M / --5
S08-0440 / Eurasia / M / --5
S08-0442 / Eurasia / M / --5
S08-0444 / Eurasia / M / --5
S08-0452 / Eurasia / M / --5
ACVC-98-02 / S.J. Valley / M / 1998 / Holarctic / 7 / N-7
MVZ-207040 / S.J. Valley / M / 2002 / Eastern / 9 / E-9
ESRP-R002 / S.J. Valley / M / 2003 / Eastern / E / 9 / E-9
ESRP-R003 / S.J. Valley / M / 2003 / Eastern / E / 9 / E-9
ESRP-R004 / S.J. Valley / M / 2003 / Eastern / E / 9 / E-9
ESRP-R005 / S.J. Valley / M / 2004 / Eastern / E / 9 / E-9
ESRP-R006 / S.J. Valley / M / 2004 / Eastern / E / 9 / E-9
ESRP-R007 / S.J. Valley / M / 2004 / Eastern / E / 9 / E-9
ESRP-R008 / S.J. Valley / M / 2004 / Eastern / E / 9 / E-9
ESRP-R009 / S.J. Valley / M / 2004 / Eastern / E / E-9
ESRP-R010 / S.J. Valley / M / 2002 / Eastern / E / 9 / E-9
ESRP-R011 / S.J. Valley / M / 2002 / Eastern / E / 9 / E-9
ESRP-R012 / S.J. Valley / M / 2003 / Eastern / F / 12 / F-12
EV37 / S.J. Valley / M / 2002 / Holarctic / 7 / N-7
EV38 / S.J. Valley / M / 2002 / Holarctic / 7 / N-7
EV40 / S.J. Valley / M / 2002 / --5
FRC-027 / S.J. Valley / M / 1997 / Holarctic / N / 7 / N-7
FRC-061 / S.J. Valley / M / 1999 / Holarctic / G / 38 / G-38
FRC-087 / S.J. Valley / M / 2000 / Holarctic / G / 38 / G-38
MAM-1937 / S.J. Valley / M / 2004 / Holarctic / N / 7 / N-7
MAM-1938 / S.J. Valley / M / 2004 / Holarctic / N / 7 / N-7
MAM-1974 / S.J. Valley / M / 2004 / Holarctic / N / 7 / N-7
MAM-2076 / S.J. Valley / M / 2004 / Holarctic / N / 7 / N-7
MAM-2703 / S.J. Valley / M / 2005 / Eastern / F / 14 / F-14
MAM-2704 / S.J. Valley / M / 2005 / Eastern / F / 14 / F-14
MAM-2705 / S.J. Valley / M / 2005 / Eastern / F / 14 / F-14
MAM-2706 / S.J. Valley / M / 2005 / Eastern / F / 14 / F-14
MAM-2748 / S.J. Valley / M / 2005 / Eastern / G / 38 / G-38
MAM-2793 / S.J. Valley / M / 2005 / Eastern / E / 9 / E-9
MAM-2794 / S.J. Valley / M / 2005 / Holarctic / G / 38 / G-38
MVZ-208595 / S.J. Valley / M / 1997 / Eastern / F / 12 / F-12
V01-0197 / S.J. Valley / M / 2001 / --5
V02-1003 / S.J. Valley / M / 2002 / --5
MVZ-90938 / Rockies / H / 1938 / macroura / Widespread / A6 / 37 / A6-37
USNM-120523 / Rockies / H / 1902 / macroura / -- / A / --
USNM-128536 / Rockies / H / 1903 / macroura / Mountain / A / 19 / A-19
USNM-130371 / Rockies / H / 1903 / macroura / Mountain / A / 19 / A-19
USNM-130372 / Rockies / H / 1903 / macroura / Widespread / A / 68 / A-68
USNM-159348 / Rockies / H / 1909 / macroura / O-24 subclade / A3 / 59 / A3-59
USNM-168812 / Rockies / H / 1910 / macroura / O-24 subclade / A3 / 59 / A3-59
USNM-210204 / Rockies / H / 1915 / macroura / -- / A / --
USNM-211097 / Rockies / H / 1916 / macroura / O-24 subclade / A3 / 59 / A3-59
USNM-213097 / Rockies / H / 1915 / macroura / Mountain / A / 19 / A-19
USNM-213098 / Rockies / H / 1915 / macroura / Mountain / A / 19 / A-19
USNM-213110 / Rockies / H / 1916 / macroura / Mountain / A / 19 / A-19
USNM-213111 / Rockies / H / 1915 / macroura / Mountain / A / 19 / A-19
USNM-214792 / Rockies / H / 1915 / macroura / Mountain / A / 19 / A-19