Estimating Trend in Occupancy for the Southern Sierra Fisher Martes pennanti Population

The fisher Martes pennanti (Figure 1) is a forest-dwelling mustelid carnivore whose primary habitat in the western United States is dense coniferous forest, usually with a deciduous component and abundant physical structure in the form of large trees, snags, and logs (Buskirk and Powell 1994; Powell and Zielinski 1994; Lofroth et al. 2010; Raley et al. 2012). In the West, the fisher occurred historically throughout the northern Rocky Mountains, Cascade and Coast ranges, and the Sierra Nevada (Gibilisco 1994). The range and abundance of this forest-dwelling carnivore have decreased in this region due to commercial trapping, changes in forest structure associated with logging and altered fire regimes, increased human access, and habitat loss to urban and recreational development (Powell and Zielinski 1994; Zielinski et al. 2005; Lofroth et al. 2010). In California, the fisher occupies less than half of its historical range, as described in the early 1900s (Grinnell et al. 1937), and a remnant population in the southern Sierra Nevada is separated by >400 km from the nearest population in northern California (Aubry and Lewis 2003; Zielinski et al. 2005). Recent genetic research has found, contrary to the perspective of Grinnell et al. (1937), that the northern California and southern Sierra Nevada populations may have been genetically isolated from each other prior to European settlement (Knaus et al. 2011; Tucker et al. in press). In 2004, the fisher was deemed “warranted but precluded by higher priority actions” for listing under the US Endangered Species Act (ESA 1973, as amended; USFWS 2004).

Although fisher conservation is one of competing concerns for land managers in the Sierra Nevada, the species has gained attention due to the perceived conflict between its habitat requirements and the growing need for fire and fuel management (USDA 2001, 2004; Scheller et al. 2011). Fisher habitat in the Sierra Nevada also occurs at elevations (1,500–2,900 m) where fire risk is also greatest (Miller et al. 2009; Scheller et al. 2011; Spencer et al. 2011) and where numerous mountain communities occur within the urban–wildland interface. Managers of public lands and stakeholders are uncertain about the trade-off between the loss of fisher habitat that may be affected by fuels treatments compared with the loss of habitat that will occur via uncharacteristically severe wildfire if forests are untreated. This uncertainty occurs as we also learn about the relatively high mortality rates of fishers in some portions of their range in California (Higley and Mathews 2009; SNAMP 2010) and new sources of mortality (e.g., rodenticide; Gabriel et al. 2012). The loss of occupied range, the vulnerability of fishers to wildfire and fuels treatments, and the potentially high rate of mortality in a species with a low reproductive rate, make it especially important that populations be monitored.

The size of the population in the southern Sierra Nevada is unknown, but various estimates have ranged from 100 to 400 individuals (Lamberson et al. 2000; Spencer et al. 2011). This small number is of concern, especially because the population is effectively isolated from the nearest population in northwestern California. It is likely that the latent effects of fur-trapping, the loss of structurally complex forests, the reduction in large-diameter trees, and the fragmentation of habitat by roads and residential development (McKelvey and Johnson 1992; Franklin and Fites-Kaufmann 1996; Campbell 2004) have contributed to a reduction of the fisher's range in the Sierra Nevada and the difficulty that dispersing animals may have recolonizing historically occupied areas.

In the 1990s, the U.S. Department of Agriculture Forest Service began a land-management planning process that resulted in Sierra Nevada Forest Plan Amendments (“Framework”) in 2001 and 2004. Both Sierra Nevada Forest Plan Amendments included plans to address the risk of uncharacteristically severe fire and to evaluate the effects of proposed treatments on potentially vulnerable species of wildlife such as the fisher. The Final Environmental Impact Statement for the Sierra Nevada Framework included monitoring for key species, including the fisher (USDA 2001, Appendix E). The adaptive management strategy specified monitoring to address the following question: “What is the status and change of the geographic distribution, abundance, reproductive success, and survivorship of the fisher population?” (USDA 2001:67, Appendix E-69). Initial attention focused entirely on the geographic distribution portion of this question, and monitoring the status and the change in the geographic range and relative abundance of fishers was the strategy to determine whether the actions prescribed in the Sierra Nevada Forest Plan Amendment would benefit the fisher population.

Fishers are easily detected using noninvasive survey methods (Long et al. 2008), and these methods have been used in standard protocols (Zielinski and Kucera 1995; Zielinski et al. 2005) to generate systematically collected, independently verifiable (McKelvey et al. 2008), and spatially precise detection data. Data from systematic surveys using track-plate and camera stations in the 1990s and early 2000s resulted in detection and nondetection data that were used to build statistical landscape-scale habitat models for the fisher in various portions of California (e.g., Carroll et al. 1999; Davis et al. 2007; Zielinski et al. 2010). Thus, when the circumstances arose to propose a Sierra-wide inventory and monitoring program for fishers, the survey methods available to do so had been well-tested.

The methods to detect presence or absence of fishers can be used to determine geographic distribution and occupancy, but they do not directly measure abundance or provide estimates of reproduction or survival. By monitoring occupancy (presence–absence), we make the key assumption that changes in occupancy reflect changes in population size (e.g., Noon et al. 2012). The survey design for fisher monitoring in the southern Sierra included a systematically arranged grid of sample units sampled for a fixed duration of time to result in either detection or no detection of at least one individual of the species (Zielinski and Mori 2001). The proportion of sample units with a detection (P) was proposed as an index of population status and trend, and therefore as a means to monitor the population. Zielinski and Mori (2001) proposed to evaluate the null hypothesis of no change in the trend in the proportion of occupied units using logistic regression (Hosmer and Lemeshow 2000). The sampling was designed to detect at least a 20% decline (1-sided alternative hypothesis) in the first 10 y of sampling (equivalent to a 2.45% annual decrease). The Type I error rate was set at 0.20, and a Type II error rate of ≥0.20 (i.e., statistical power [1 − β] of ≥80%) was selected. Sampling demands for a 2-sided alternative hypothesis led to projected costs that exceeded the expected annual budget, and detecting an increase in fishers was of less interest than detecting a decrease. A priori power analysis revealed that these standards could be achieved by sampling ≥288 sample units/y.

Monitoring following the basic approach described above occurred from 2002 to 2009. During this period, however, an important change occurred in the methods recommended for the analysis of binomial detection data (“occupancy” data). The logistic regression approach to analyzing binomial data of this nature was being eclipsed by modern occupancy modeling and estimation (MacKenzie et al. 2006), combined with multimodel methods of statistical inference (Burnham and Anderson 2002). A number of researchers had already been exploring the use of these approaches to monitor populations of species within the genus Martes (see Slauson et al. 2012). Thus, we adapted these newer analytical methods to the original design.

We report here the results of 8 y of sampling and analysis to determine occupancy and trends in occupancy for fishers in the southern Sierra Nevada from 2002 to 2009. We consider the effect of probability of detection, which was a component of the work from the beginning, and model the combined effects of detection and occupancy to estimate occupancy persistence rates (related to colonization and extinction rates; MacKenzie et al. 2003, 2005). The broad objective is to fit a model to the observed fisher data to assess the status and trend of fisher occupancy dynamics. Occupancy can be predicted by including environmental covariates (MacKenzie et al. 2006), but we do not include covariates in this analysis. Our primary goal is to report our approach to estimating occupancy for fishers and whether or not the results of 8 y of monitoring suggest a trend in occupancy rates.

The study area is, broadly, that portion of the historical range of the fisher from Yosemite National Park and Sierra National Forest south to the end of the Sierra Nevada (Figure 2). No surveys conducted prior to the beginning of the monitoring period (2002) had detected a fisher in the area bounded by 30 mi (48 km) east of Interstate 5 in Shasta County south to the Merced River in Yosemite National Park (Zielinski et al. 2005). Thus, the sampled population was defined to include largely public lands occurring within the approximate historical distribution of the fisher (as described by Grinnell et al. 1937) south of CA Highway 120 and west of U.S. 395. This area includes Sequoia National Forest (NF), Sierra NF, the southern portion of Yosemite National Park, small portions of Stanislaus NF, and Inyo NF. The annual, systematic survey data analyzed here were collected only within the general area presumed to be occupied by fishers (Zielinski et al. 2005); however, some sampling occurred in the presumably unoccupied region to the north and east.

We collocated our monitoring sample units with grid points included in the Forest Inventory and Analysis (FIA) system (Bechtold and Patterson 2005). The grid of “phase 2” FIA points was provided by the U.S. Department of Agriculture Forest Service Pacific Northwest Research Station and served as the sampling frame for monitoring. The FIA system is based on a system of hexagons, each 2,428 ha (6,000 ac) in size, distributed across the landscape, with the centers of each hexagon 5.3 km (3.27 mi) apart (Reams et al. 2005). Each FIA point, however, is randomly located within a hexagon such that adjacent points can, in practice, be closer or farther apart. To assure consistency with the 1985 Food Security Act and the FIA National Interim Privacy Policy, a Memorandum of Understanding was developed between the Pacific Northwest Research Station and the Pacific Southwest Region of the Forest Service, who administered the monitoring program that allowed sampling in close proximity to the FIA plots but not actually overlapping the FIA plots (offset from their true locations 100 m in a random direction).

Our first step in sample selection was to identify the FIA points occurring within the historical fisher distribution described by Grinnell et al. (1937). We eliminated FIA points from the sampling frame that occurred in cover types known to be unsuitable (e.g., sagebrush, grassland, urban) as well as those in the least suitable elevation range (below 800 and above 3,200 m). This resulted in 596 FIA points on all ownerships (388 on Forest Service lands). We originally planned to sample half the units in alternating years and did so the first 2 y (2002, 2003). Due to the random nature of this sample-unit selection, the reduced sample in these years was still considered a statistically valid annual sample for the purpose of estimating occupancy. Starting in 2004, however, we created a final set of 223 core sample units that were within the most likely elevations (1,140–2,560 m), were primarily on national forest land, and excluded units that were dangerous to access or extremely remote (>14 km from road access). None of the sample units was selected because of their known or suspected fisher occupancy status. The great majority of the core sample units were on national forests (92.3%), but national park, state, and private lands comprised 5.0%, 1.0%, and 1.7% of the sample units, respectively. The core 223 monitoring sites included an area of approximately 12,240 km². We planned to sample every sample unit annually, but this became financially and logistically difficult, so the total number of sites sampled varied somewhat from year to year. Whenever a sample unit was resampled, it was always located at the same geographic point. The mean minimum distance (SD) between pairs of points within the core set of sample units was 4.1(1.2) km.

In addition to conducting an analysis of the entire dataset, which we refer to as the combined analysis, we also identified three geographic zones within the study area where preliminary data (Zielinski et al. 2005) and local knowledge suggested that fisher abundance and occupancy may vary: the “southwestern” (primarily the west side of Sequoia NF; n  =  66 sample units), the “northwestern” (primarily the Sierra NF; n  =  126 sample units), and the “southeastern” (primarily the Kern Plateau; n  =  31 sample units) zones (Figure 2). The potential differences led us to evaluate whether probability of detection or occupancy varied among these zones (see below).

We used enclosed track stations (Ray and Zielinski 2008) as the primary detection method at sample units. In this baited device, a fisher enters the enclosure and leaves a sooted track impression on a track-receptive surface. Track stations have been used to detect fishers and martens Martes americana in California for years (Zielinski et al. 1995, 2005; Aubry and Lewis 2003) and the tracks of fishers can be quantitatively distinguished from those of marten (Zielinski and Truex 1995). Each sample unit had six track stations arrayed such that a central station (aligned with the offset FIA grid point as described above) was surrounded by five other stations positioned at 72° intervals, approximately 500 m from the center (Figure 3). From 2006 to 2008, all six track-plate stations were modified to include a barbed-wire hair snare (Zielinski et al. 2006a). In 2009, to reduce installation and servicing time, barbed-wire hair snares were only used at sample units that had previously detected a fisher, or were added midsurvey to units with a new fisher detection. Previous work has demonstrated that a hair snare does not affect visits by fishers to a station (Zielinski et al. 2006a).

We baited the track stations with a piece of raw chicken and applied a commercial trapping lure (Gusto; Minnesota Trapline Products, Pennock, MN). Each sample unit was deployed for 10 d, and we replaced bait and track plate every 2 d for five “visits” during each sample period. This 10-d sample period was assumed to be sufficiently short in duration to assure a reasonable degree of population closure. The lure was only applied once, at the beginning of the sample period. Multiple visits to sample units allowed us to estimate probability of detection and account for this parameter when estimating occupancy over time (MacKenzie et al. 2006; Slauson et al. 2012). A sample unit was registered as having detected a fisher if, after we checked the unit, at least one fisher detection was verified at any of the six stations within the unit.

The observed status (presence or absence) at a sample unit during a single visit is affected by the probability that a fisher is present at the site and the probability that the fisher is also detected. Thus, the probability of detection must be considered during the process of modeling occupancy (presence; MacKenzie et al. 2006). Models for detection probability and probability of occupancy were, ultimately, combined but we describe their components separately first.

We consider four models for “detection probability”:

1. A model where each visit to an occupied site has a common probability of detection and visits are independent of each other. We label the detection probability p and identify this model as p(.).

2. A model where each visit to an occupied site has, within a specific year, a common probability of detection and visits are independent of each other. This is the same as Model (1) but with the detection probabilities varying among years. We label this model as p(Year).

3. A model that postulates that there is a mixture of two types of sample units common among years: type A sample units have detection probability p_A and type B sample units have detection probability p_B. The proportion of units of type A is π and the proportion of type B is 1 − π. This model accounts for a particular type of overdispersion (greater variability in the number of detections than expected if the sample conformed to the distribution of the simpler model). In our case, overdispersion can result when there are two different kinds of sample units and p_A is much different from p_B. In turn, this results in a larger number of occupied sample units having no detections than expected under a model with a single constant detection probability. We label this model as p(2 group) and call it the “2-group model” (see Pledger 2000).

4. A model that postulates that there is a mixture of two types of sample units but the parameters potentially vary among years. This is the same probability structure as Model (3) but with separate sets of parameters for each year. We label this model as p(Year · 2 group).

For the probability of occupancy we considered a slightly different, but mathematically equivalent, formulation of the usual extinction–colonization parameterization. The usual parameterization over T seasons is for a Markov chain with the first year probability of presence (ψ₁) and the extinction probabilities (ε₁, ε₂,…,ε_{T − 1}) and colonization probabilities (γ₁, γ₂,…,γ_{T − 1}). However, because our main objective concerns examining linear trends in the probabilities of presence (ψ₁, ψ₂,…ψ_T; and potentially other smooth changes over time), we started with the probabilities of presence ψ₁, ψ₂,…,ψ_T. Once these are set, the potential range for the extinction and colonization probabilities are more restricted than just between zero and one. We used the directional proportion (labeled φ_t at year t) of the distance that the extinction probability is away from the value associated with independence. At the value of independence (ε_t  =  1 − ψ_{t + 1}) we have φ_t  =  0. At the largest possible value of

we have φ_t  =  −1. At the smallest possible value of

we have φ_t  =  +1.

We call φ the persistence factor (after the parameter used in Barton et al. 1962), and it ranges (like a correlation coefficient) from −1 to +1 with zero indicating independence. Values larger than zero indicate positive persistence with sites more likely not changing status from one season to the next than expected under independence. Values less than zero indicate negative persistence with sites more likely changing status than under independence.

Using the original (but again, equivalent) parameterization is useful when one wants to examine potential mechanisms for changes in probabilities of presence. When there is no change in the probabilities of presence over time, an examination of any internal changes in persistence can be performed equivalently with extinction or colonization probabilities or with the persistence factor.

We consider three models for the “probability of occupancy” (presence) and three models for “persistence.” First we describe the models for occupancy:

1. A model where the probability of occupancy is ψ and is constant for all years. We label this model ψ(.).

2. A model where the probability of occupancy varies among years. We label this model as ψ(Year).

3. A model where the logit of the probability of occupancy is a linear function of year. Here we have

Finally, we have three models for persistence:

1. A model of independence. Here we have φ_t  =  0, which corresponds to a model with ε_t  =  1 − ψ_t+1.

2. A model of constant persistence labeled φ(.). Here we have ε_t being a constant directional proportional distance from 1 − ψ_t+1, which is the value of ε_t under independence.

3. A model where the persistence varies by year. We label this model φ(Year).

Combining the models for detectability, occupancy, and persistence (Table 1), we considered 36 possible models for each zone and for the combined dataset (Table S1, Supplemental Material). The parameterization and model structure is unique and, thus, is described below for the following hypothetical observed detection history at a single sample unit (Table 2).

In this example, there are visits in 6 of the 8 y with detections in the last 3 y. To construct the likelihood contribution for this sample unit, we first need to consider all possible true annual status histories that could have resulted in such an observed history. The observed annual status is given as 0.0.0111 where “0” represents no detections, “.” represents no visits, and “1” represents at least one detection. There are 32 different possible true annual status histories consistent with the observed history, starting with 00000111 and ending with 11111111 (in lexicographic order), with the last 3 y always being a “1” (because we assume no false positives). That set of potential annual status histories is labeled H_i for sample unit i.

The probability of each potential true annual status history is calculated in the usual manner where we convert the ψ and φ values (“persistence model” parameterization) to ψ₁ and the ε and γ values (“colonization–extinction model” parameterization) because it is more straightforward to calculate the probability of the true annual status history in a hierarchical manner with the colonization–extinction parameterization. For example, if we consider the model labeled ψ(logit), φ(.), p(2 groups), the probability of the true annual status history only depends on the ψ(logit) and φ(.) parts of the model. We have a constant value of φ, and for the annual probabilities of presence we have two parameters, a and b, where

Given a and b, we calculate ψ₁, ψ₂, ψ₃, ψ₄, ψ₅, ψ₆, and ψ₇ (where ψ₁ is the probability of occupancy for the year 2002, in this example). Then we determine the corresponding values of ε₁, ε₂, ε₃, ε₄, ε₅, ε₆, and ε₇ and γ₁, γ₂, γ₃, γ₄, γ₅, γ₆, and γ₇:

For the hypothetical observed history (h  =  00110111; 1 of the 32 possible true histories), the following is the occurrence probability:

We then multiply the above by the probability of the particular detection history for each year. With the 2-group model, we allow for a site to be one of two types with respective probabilities of detection p₁ and p₂. These two types occur with probability π and 1 − π, respectively. If, for true detection history h, there are n_h total visits in years with presence and x total detections across all years, then the probability of that complete detection history is

Finally, the log-likelihood contribution for site i with specific values of a, b, φ,π, p₁, and p₂ is given by

All programming for the maximum likelihood estimators was performed using SAS software, Version 9.3 of the SAS System for Windows (SAS Institute, Inc., Cary, NC; Table S2 and Table S3, Supplemental Material). PROC NLMIXED in SAS was used as the optimizing routine. We compared the relative fit of the models using Akaike's Information Criteria (AIC; Akaike 1973; Burnham and Anderson 2002).

Models characterized with ψ(.) assume no year-to-year change in occupancy. Thus, if the data best fit a model with this single occupancy parameter (i.e., models containing ψ(.)), then we conclude that there is no evidence of a change in occupancy over time. However, we also evaluated trend an alternative way, by fitting a model that is linear in the logit of the probability of presence. If the slope parameter of this model is different from zero, it would suggest a significant trend in the probability of presence. To characterize an average annual rate of change, we estimated the following:

where ψ₂₀₀₂ and ψ₂₀₀₉ are the respective probabilities of presence in 2002 and 2009. Standard errors were then estimated using the delta method (Bishop et al. 1975).

A total of 223 sample units were sampled during 2002–2009, but not all were surveyed each year (Table S1, Supplemental Material). An average of 139.5 units was sampled/y, and the average sample unit was surveyed for 5 of the 8 y (Table 3; Figure 4). The greatest number of sample units was in the northwestern zone and the fewest in the southeastern zone (Table 3; Figure 2). We recorded an average of 0.255 sample units with fisher presence per year, ranging from 0.203 in 2004 to 0.279 in 2005 (Table 3). Across all years, the average (SD) number of stations per sample unit with a fisher detection was 2.05 (1.35), the average (SD) number of detections per sample unit was 3.03 (3.29) and the average (SD) number of days to the first detection at a sample unit was 4.90 (2.85).

Model ψ(.), φ(.), p(2 groups), which assumes a constant probability of occupancy, constant persistence factor, and the 2-group detection model clearly outperformed the other models (AIC weight  =  0.707; Table 4; Table S1, Supplemental Material) for all zones except the northwestern zone (where this model was ranked second). This indicates no trend or statistically significant variation in occupancy over the 8-y period (Figure 5). Moreover, the same model accounted for the majority of the AIC weight in two of the three zones and in the combined dataset, suggesting that there was no change in fisher occupancy rates in any of the three zones (Table 4; Figure 5). The top model for the northwestern zone was the ψ(logit), φ(.), p(2 groups) model and the slope parameter was not significantly different from zero, thus this line of evidence suggests no trend in the northwestern zone.

The probability of occupancy, adjusted to account for uncertain probability of detection, was 0.367 for the combined dataset and ranged from a low of 0.261 in the southeastern zone to 0.583 in the southwestern zone (Table 5). Positive persistence of sample unit status best fit the data, meaning that sample units tend not to change status (from absence to presence or from presence to absence) as often as expected if status was independent from year to year. As specified by the top model, for the combined dataset, there was a 74.4% chance of a fisher at sample units where there was a fisher the year before (one minus the extinction probability), and an 89.8% chance that sample units without a fisher also had no fishers the previous year (one minus the colonization probability). For comparison, a randomly chosen site had only a 63.3% chance of no fisher (one minus the probability of presence). Persistence factors were the lowest for the southeastern zone (52.7%; Table 5), indicating a more dynamic occupancy pattern over time for this region compared with other regions but still more persistent than if occupancy status was independent from year to year.

The values of the 2-group detection parameters suggest that there were more sample units where fishers were present, but undetected, than expected from the fit of the simpler detection model. This, in turn, increases the adjustment applied to the naïve (i.e., unadjusted for probability of detection) estimates of the proportion of sites with presence. The 2-group detection probability model resulted in a probability of detection of 0.778 and 0.179 for each group, and a weighted mean probability of detection for a single visit of 0.31 for all zones combined. Weighted mean probability of detection for a single visit was 0.17, 0.32, and 0.42 for the northwestern, southeastern, and southwestern zones, respectively. For the 5-visit protocol used here, this translates to an overall detection probability, if a fisher is present, of 0.71 with five visits for all zones combined, and overall detection probabilities of 0.51, 0.82, and 0.82 for the northwestern, southeastern, and southwestern zones, respectively. Adjusting for probability of detection in the models led to annual estimates of occupancy (Table 5) that were substantially higher than the naïve estimates (Table 3).

The lack of evidence for a trend was also supported by determining that the estimates of slope for the models where occupancy was allowed to vary (i.e.,ψ(logit)) did not differ from zero (Table 6; Figure 5). Estimates of annual rate of change of the probability of presence for the combined dataset varied from −0.021 (SE  =  0.013) to 0.029 (SE  =  0.022; Table 6) and none of the estimates were significantly different from zero.

Model evaluation and trend assessment revealed no evidence for a change or trend in fisher occupancy estimates over the 8-y period. There were no statistically significant slopes in any of the models with a trend parameter (specifically all models containing ψ(logit)). More important, perhaps, is the fact that the model that best fit the data for each zone did not include the effect of year for all but one zone. And while that one zone (northwestern) had a trend parameter estimated, it was not found significantly different from zero. Thus, we have no evidence that the fisher population, as indexed by our measure of occupancy, has changed in the southern Sierra Nevada, or any zone therein, from 2002 to 2009.

Occupancy in the southwestern zone of the study area was about twice that of the northwestern or southeastern zones (Table 5), which may mean that the southwestern portion of the study area is either closer to carrying capacity than is the remainder of the study area, or has a higher capacity. Persistence was only slightly higher in the southwestern zone than the northwestern zone, but both were higher than persistence in the southeastern zone. Additional analyses are needed to determine whether habitat or other factors, such as a large fire on the Kern Plateau in 2002, contribute to our findings. These differences in occupancy and persistence among zones may have important implications for the transferability of results from ongoing fisher population studies being conducted in the northwestern (SNAMP 2010) and southwestern (Thompson et al. 2011) zones to other parts of the southern Sierra fisher range.

The positive persistence factors that characterize the patterns of occupancy (combined estimate  =  0.81) suggest that sample units tend not to change status (from absence to presence or from presence to absence) as often as expected if status was independent from year to year. Fishers are territorial and relatively long-lived animals (5–10 y; Powell 1981); thus, our finding that sample units remain occupied (or unoccupied) from year to year is not surprising. Moreover, relative habitat suitability is heterogeneous in the southern Sierra Nevada (Davis et al. 2007; Spencer et al. 2011) and places where occupancy is consistent may be places of relatively high suitability. A confounding factor, however, is the possibility that positive persistence values are affected by our ignorance of how many individual fishers are available to be detected at each sample unit. Preliminary analysis of hair snared at sample units revealed that 13.8–22.2% of the units had evidence of visits by at least two different fishers, depending on year (J. Tucker, unpubl. data). Thus, circumstances do arise when a site can maintain its status as “occupied” when all but one of the individuals that had occurred in the general vicinity die or relocate. We believe, however, that whether persistence is due to the annual detection at a sample unit of a single individual or several fishers, it does not affect the interpretation of persistence: sample units with a detection are regularly occupied by one or more fishers from year to year. Future research should investigate the environmental correlates of these locations, including predicted relative habitat suitability (Davis et al. 2007; Spencer et al. 2011). We also expect that a full analysis of the covariates associated with the probability of presence (MacKenzie et al. 2006) would shed light on the factors affecting high persistence rates, but was beyond the scope of this analysis.

The possibility that more than one fisher could be detected at each sample unit also affects the interpretation of occupancy as an index of abundance because, insofar as this occurs, there would not be a one-to-one relationship between the number of fishers detected at a sample unit and the number of sample units. The original design (Zielinski and Mori 2001) specified selecting sample units that were ≥9 km apart to assure spatial independence, but maintaining this spacing would have greatly reduced the sampling frame and implementation efficiency. During the first year of sampling, we relaxed this requirement to achieve a larger sample. The mean minimum distance was 4.1 km and about 30 sample units were approximately 3 km apart. Average summer home ranges in our northwestern zone are 15 km² for females and 55 km² for males (Sweitzer 2011). Thus, it may not be surprising that some sample units were visited by more than one individual. Indeed, spatial independence was probably an unrealistic expectation because fisher home ranges overlap, particularly those of males and females (Powell 1981). Because of sexual overlap, there will not be a one-to-one relationship between the number of sample units with a detection and the number of fishers. This does not mean, however, that sample units cannot still be independent, which requires only that the same fisher is not detected at more than one sample unit. We have less data on this possibility; however, genetic identification of individuals captured via hair samples from 2006 to 2009 indicates that the percent of individuals detected at multiple sample units each year was low, averaging 5.4%/y (J. Tucker, unpublished data). Given the general spatial structure of fisher populations, the relatively large distances between sample units, and the relatively short sample duration (10 d), we believe that our estimates of occupancy are a useful, practical, and economical index of fisher abundance. Future methods for fisher monitoring will probably focus on identifying individual fishers and will likely reduce or eliminate this problem.

A peculiarity of the favored models in this application is the consistent inclusion of the 2-group detection model instead of a simpler model with a single detection probability (Table 4). Although there are two very different detection probabilities, it is not possible to identify which sampling unit is assigned to which detection probability. The two groups are parameters that capture overdispersion in the data but the cause for the overdispersion is unknown. The method of maximum likelihood helps estimate the values of π, p_A, p_B but cannot assign a group to each sample unit. It is possible that the two groups relate to habitat covariates that are distributed differently among sample units, sex-specific detection probabilities, observer effects, nonrandom elements of site selection, or lack of independence among visits. However, we did not explore the reasons for the grouping in this exercise.

Finally, we realize that we were unable to sample as many units per year (i.e., 288) as our original simulation modeling (Zielinski and Mori 2001) suggested were necessary to detect a meaningful decline. However, this early modeling assumed independence in occupancy estimates from year to year, contrary to the strong persistence estimates that our data revealed. The predicted sample size necessary to detect a decline would have been substantially lower if dependency in persistence were assumed in the original analysis.

The fisher population in the Sierra Nevada is assumed to be at risk due to its small size, geographic and genetic isolation, and the fact that much of its historical range is unoccupied (Grinnell et al. 1937; Zielinski et al. 2005). It is encouraging that the small population in the southern Sierra does not appear to be decreasing. Given the degradation of forest condition that has occurred in the Sierra Nevada in the past century due to timber harvest and fire suppression, and the potential latent effects of fur trapping, we assume that the fisher population within its current extent has not reached equilibrium (i.e., is not at carrying capacity). Thus, conservationists would desire an increasing trend in the index of occupancy. This increase would suggest that changes in forest condition, due to either happenstance or to conservation measures taken by regulatory and management agencies, have been beneficial. We cannot conclude from these data whether management actions taken to protect fishers and their habitat are responsible for the stability of the occupancy measure over the 8-y period. Guidelines to protect den sites, reduce old forest fragmentation, and enhance canopy cover were instituted in 2001 (USDA 2001) and strengthened further in 2004 with restrictions on the time of year that forest management activities can occur near fisher dens and restrictions on the extent and magnitude of thinning forests to reduce fire threat, and with mitigation of road and recreational effects around dens (USDA 2004). However, these protective measures are unlikely to have affected the fisher population significantly during the period we collected data; a century of potential negative effects cannot be reversed in a few years of protective measures.

Ultimately, as monitoring continues, it will be important to link population monitoring data, such as those reported here, with habitat monitoring data. This will help us understand whether changes in the fisher population can be linked to changes in habitat suitability. For example, methods for assessing and monitoring fisher habitat at the resting habitat scale (Zielinski et al. 2006b), home range scale (Thompson et al. 2011), and landscape scale (Davis et al. 2007; Spencer et al. 2011) have been developed. The resting habitat model is the most refined for this purpose and it has revealed relatively stable predicted resting habitat values for the southern Sierra Nevada (Zielinski et al. 2010). We hasten to add, also, that fisher populations change for reasons other than changes in abundance or quality of habitat (e.g., predation, disease, competition), and that these elements should also be considered when evaluating future changes in the status of the fisher population.

Please note: the Journal of Fish and Wildlife Management is not responsible for the content or functionality of any supplemental material. Queries should be directed to the corresponding author for the article.

Table S1. Work sheets that include models and AIC values, parameter estimates, detection information, and additional estimates for modeling fisher Martes pennanti occupancy from detections from 2002 to 2009 in the southern Sierra Nevada mountains of California. Worksheet 1 (“AIC”) includes AIC, AIC weights, and ΔAIC values. Worksheet 2 (“Parameter Estimates”) includes the parameter estimates for all the parameters included in the models, standard errors (se), degrees freedom (df), t value (tvalue), probability (pi), lower and upper confidence interval values (lower; upper) and gradient value (gradient). Worksheet 3 (“Histories”) includes information about sample units, including the zone (zone), the Universal Transverse Mercator system coordinates (tx_u10 ty_u10), the quadrat identifier (tquad), the Forest Inventory and Analysis plot identifier (fiaid), the national forest (forest), and the number of years that the sample was included (nv1…nv8) and the number of fisher detections for each year (nd1….. nd8). Worksheet 4 includes parameter estimates for additional parameters, primarily extinction and colonization, necessary to estimate as intermediate steps in calculating likelihood.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S1 (613 KB XLSX).

Table S2. SAS (SAS Institute, Inc., Cary, North Carolina) code for analysis of fisher Martes pennanti occupancy (“fisher analysis code.sas”) from 2002 to 2009 in the southern Sierra Nevada mountains of California.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S2 (6 KB TXT).

Table S3. SAS (SAS Institute, Inc., Cary, North Carolina, USA) macros for analysis of fisher Martes pennanti occupancy (“fisher macros.sas”) from 2002 to 2009 in the southern Sierra Nevada mountains of California.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S3 (12 KB TXT).

Reference S1. Powell RA, Zielinski WJ. 1994. Fisher. Pages 38–66 in Ruggiero LF, Aubry KB, Buskirk SW, Lyon LJ, Zielinski WJ. editors. The scientific basis for conserving forest carnivores: American marten, fisher, lynx and wolverine in the Western United States. Fort Collins, Colorado: U.S. Department of Agriculture Forest Service, Rocky Mountain Research Station. General Technical Report RMRS-GTR-254.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S4; also available at http://www.fs.fed.us/rm/pubs_rm/rm_gtr254/rm_gtr254_038_073.pdf (837 KB PDF).

Reference S2. Lofroth EC, Raley CR, Higley JM, Truex RL, Yaeger JS, Lewis JC, Happe PJ, Finley LL, Naney RH, Hale LJ, Krause AL, Livingston SA, Myers AM, Brown RN. 2010. Conservation assessment for fishers (Martes pennanti) in south-central British Columbia, western Washington, western Oregon, and California. Volume I. Denver: U.S. Department of the Interior Bureau of Land Management.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S5; also available at http://www.fws.gov/yreka/PDF/Lofroth_etal_2010.pdf (28 MB PDF).

Reference S3. Higley JM, Mathews S. 2009. Fisher habitat use and population monitoring on the Hoopa Valley Reservation, California. Final Report USFWS TWG U-12-NA-1. Unpublished report.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S6 (1.1 MB PDF).

Reference S4. Lamberson RH, Truex RL, Zielinski WJ, Macfarlane DC. 2000. Preliminary analysis of fisher population viability in the southern Sierra Nevada. Arcata, California: U.S. Department of Agriculture.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S7 (168 KB PDF).

Reference S5. McKelvey KS, Johnson JD. 1992. Historical perspectives on the forests of the Sierra Nevada and the Transverse Ranges in southern California: forest conditions at the turn of the century. Pages 225–246 in Verner J, McKelvey KS, Noon BR, Gutierrez RJ, Gould GI, Jr, Beck TW, technical coordinators. The California spotted owl: a technical assessment of its current status. Albany California: U.S. Department of Agriculture Forest Service, Pacific Southwest Research Station. General Technical Report PSW-GTR-133.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S8; also available at http://www.treesearch.fs.fed.us/pubs/3536 (2.1 MB PDF).

Reference S6. Franklin JF, Fites-Kaufmann J. 1996. Assessment of late-successional forests of the Sierra Nevada. Pages 627–661 in Sierra Nevada Ecosystem Project, Final Report to Congress. Volume II, Assessments and Scientific Basis for Management Options. Davis: University California, Centers for Water and Wildlands Resources.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S9; also available at http://www.sierraforestlegacy.org/Resources/Conservation/FireForestEcology/ThreatenedHabitats/OldGrowthForests/OGF-SNEP1996Part1.pdf (3.1 MB PDF).

Reference S7. Zielinski WJ, Kucera TE. 1995. American marten, fisher, lynx, and wolverine: survey methods for their detection. Arcata, California: U.S. Department of Agriculture Forest Service, Pacific Southwest Research Station. General Technical Report PSW-GTR-157.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S10; also available at http://www.fs.fed.us/psw/publications/documents/gtr-157/ (5 MB PDF).

Reference S8. Zielinski WJ, Mori S. 2001. What is the status and change in the geographic distribution and relative abundance of fishers? Arcata, California: U.S. Department of Agriculture Forest Service, Pacific Southwest Research Station. Adaptive Management Strategy, Sierra Nevada Framework, Study Plan.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S11 (365 KB PDF).

Reference S9. Bechtold WA, Patterson PL. 2005. The enhanced forest inventory and analysis program: national sampling design and estimation procedures. Asheville, North Carolina: U.S. Department of Agriculture Forest Service, Southern Research Station. General Technical Report SRS-GTR-80.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S12; also available at http://www.srs.fs.usda.gov/pubs/gtr/gtr_srs080/gtr_srs080.pdf (2.7 MB PDF).

Reference S10. Reams GA, Smith WD, Hansen MH, Bechtold WA, Roesch F, Moisen GG. 2005. The forest inventory and analysis sampling frame. Pages 11–26 in Bechtold WA, Patterson PL, editors. The enhanced forest inventory and analysis program: national sampling design and estimation procedures. Asheville, North Carolina: U.S. Department of Agriculture Forest Service, Southern Research Station. General Technical Report SRS-GTR-80.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S13; also available at http://www.srs.fs.usda.gov/pubs/gtr/gtr_srs080/gtr_srs080-reams001.pdf (861 KB PDF).

Reference S11. [SNAMP] Sierra Nevada Adaptive Management Project. 2010. Fall 2010 newsletter: vol. 4 no. 2 — fisher team.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S14; also available at http://snamp.cnr.berkeley.edu/news/2010/sep/20/fall-2010-newsletter-vol-4-no-2-fisher-team/ (1 MB PDF).

Reference S12. [USDA] U.S. Department of Agriculture. 2001. Sierra Nevada forest plan amendment: final environmental impact statement. Vallejo, California: U.S. Forest Service, Pacific Southwest Region.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S15; also available at http://www.sierraforestlegacy.org/Resources/Conservation/LawsPoliciesRegulation/KeyForestServicePolicy/SierraNevadaFramework/Framework-FSROD01.pdf (805 KB PDF).

Reference S13. [USDA] U.S. Department of Agriculture. 2004. Sierra Nevada forest plan amendment: final supplemental environmental impact statement. Vallejo, California: U.S. Forest Service, Pacific Southwest Region.

Found at DOI: http://dx.doi.org/10.3996/012012-JFWM-002.S16; also available at http://www.fs.usda.gov/Internet/FSE_DOCUMENTS/fsbdev3_046095.pdf (775 KB PDF).

We thank the U.S. Department of Agriculture Pacific Southwest Region and Pacific Southwest Station for financial support, the Sequoia National Forest for their administrative support, especially J. Whitfield, and dozens of field technicians for their assistance. J. Werren provided assistance with the figures and R. Schlexer provided editorial assistance. This manuscript was improved by the comments we received from the Subject Editor and anonymous reviewers.

Any use of trade, product, or firm names is for descriptive purposes only and does not imply endorsement by the U.S. Government.