The Luangwa River forms most of the eastern border of SLNP and the western border of Lupande GMA. As the largest perennial water source in the region, wildlife and human activity is centered along the river boundaries of SLNP and adjacent GMAs, particularly during the dry season (May–November). Large mammals move across the Luangwa River freely, particularly in the late dry season (Rosenblatt et al. 2014 ); thus, we considered our two study areas as segments of the same leopard population, which we termed the western study area (WSA) and eastern study area (ESA) (Fig. 2 ). The WSA (172 km 2 ) is located in SLNP on the western side of the Luangwa River ranging between the seasonal Katete and Luwi Rivers and includes areas within 6 km of the Luangwa River. The ESA (141 km 2 ) is located on the eastern bank of the Luangwa River, bounded by the seasonal Mwangazi and spring‐fed Chichele streams and also includes areas within 6 km of the Luangwa River. The ESA includes portions of the Lupande GMA and the Nsefu sector of SLNP, a small portion of the park situated on the eastern side of the Luangwa River. Other than differences in human use, the ESA and WSA were selected to be ecologically similar, with comparable compositions of edaphic grassland, deciduous riparian forest, miombo ( Brachystegia spp ) woodland, mopane ( Colophospermum mopane ) woodland and scrubland, dry deciduous forest and undifferentiated woodland (Astle et al. 1969 ; White 1983 ; Astle 1988 ). Our sampling on the ESA and WSA was designed to provide strong inferences in the following ways: (1) The two sites were selected to be similar for variables (other than those directly related to the level of protection) that would be expected to influence leopard density or demography (e.g., vegetation type and proximity to permanent water). (2) The two sites were of the same size, were spatially close, and were sampled over highly overlapping time periods. (3) Within primary sampling periods, the sampling design for both sites was identical.

Our study was conducted in two adjoining areas on the boundary of SLNP (S13.07958 E31.77407; 9050 km 2 ) and Lupande GMA (5660 km 2 ), allowing comparison of leopard densities across management regimes differing in the degree and type of human activity (Fig. 2 ). SLNP is strictly protected as an IUCN Category II Protected Area, with photographic safari tourism and law enforcement patrols, although some illegal wire‐snare and rifle poaching does occur. GMAs are IUCN Category VI areas intended as buffer zones allowing a variety of natural resource‐based uses (Chomba et al. 2011 ). Human settlements are permitted and increasing in the Lupande GMA (and most GMAs in Zambia), causing an array of conservation concerns (Lewis and Phiri 1998 ; Becker et al. 2013a ; Watson et al. 2013 , 2014 ). Legal trophy hunting of adult male leopards, other large carnivores and herbivores occurred in Lupande and other GMAs, except during a January 2013–April 2015 moratorium. Livestock densities were locally low, making human–carnivore conflict uncommon relative to other studies (e.g., Marker and Dickman 2005 ).

We established 25 and 26 unbaited camera‐trap sites (hereafter sites) in the ESA and WSA, respectively (Fig. 2 ). We selected sites by searching for leopard tracks within 100 m of each grid point (Silver et al. 2004 ). If we encountered no tracks, we selected the most active game trail within 100 m of the point. Cameras were attached to trees at a height and angle intended to maximize the likelihood of being triggered by leopards. One site was located more than 100 m from the grid point, because no trees were available within the 100 m radius. Vegetation varied between sites but all were in vegetation types used by leopards (Balme et al. 2007 ). At each site we set two Reconyx Hyperfire PC800 cameras (Reconyx, Inc., Holmen, WI) facing each other to photograph both sides of passing leopards, set to take five photographs in succession upon detection of movement. We visited sites on foot in small groups to minimize our potential impact on subsequent detections. We downloaded photographs when cameras were moved between locations (see below). We identified individuals using spot patterns and sexed them using genitalia and sexually dimorphic traits such as body and head size and the prominence of their neck dewlap (Balme et al. 2012 ). We did not assign ages due to the limitations of image quality and the difficulty of aging leopards accurately (Balme et al. 2012 ). We created capture histories for each individual denoting detections (1) and nondetections (0) on each day of camera trapping.

We used a systematic camera grid to photograph leopards within each study area and used closed robust design capture–recapture models to estimate population size and annual survival rates (Pollock 1982 ; Kendall et al. 1995 ). In both study areas, we placed cameras using a square grid that was random in its origin and orientation. Spacing for this grid followed established procedures for large felids to meet the assumptions of closed mark–recapture models (Otis et al. 1978 ; Karanth and Nichols 1998 ; Balme et al. 2009a ). We based grid cell size on the smallest home range estimate available (14 km 2 ) for an adult female leopard in Zambia's Luambe National Park (approximately 60 km from our study; Ray 2011 ), and spaced trap sites 2.5 km from each other (Fig. 2 ). This spacing was intended to place multiple trap sites within the home range of each individual (Karanth and Nichols 1998 ).

Robust design model selection

We used an extended robust design model to estimate population size (N) for each study area, annual survival (S), detection and redetection probabilities (p and c respectively), and rates of temporary emigration (γ″, the probability of an individual temporarily moving off of the study area and becoming unavailable for capture and γ′, the probability of an individual remaining outside of the study area and thus remaining unavailable for capture). We hypothesized that density would vary by study area and annual survival rates would vary by gender, study area, and time, with a potential interaction between gender and study area. We compared models with random, Markovian, or no temporary emigration (Kendall et al. 1997). Finally, we tested for effects on detection probability of study area and season (see below) and whether detection probability differed from redetection probability to evaluate whether our activity at trap sites was impacted leopard behavior.

For each primary sampling occasion on each study area, we estimated population size during an 87‐day period, a period we selected a priori to satisfy the assumption of population closure based on other large felid camera‐trap studies (Karanth and Nichols 1998). Each of these primary sampling occasions fell entirely within the cold dry season (CD, May–August) or the hot dry season (HD, September–November), so we refer to primary sampling occasions as seasons hereafter. Wet season data collection was not possible because portions of both study areas were inaccessible. Because we had too few cameras to survey all locations simultaneously, we deployed camera traps in a random rotation across four “sections” within each study area (Karanth 1995). Each section consisted of 6–7 sites and was sampled for 21 days (21 days/section X four sections + 3 days to redeploy cameras = 87 days). We created encounter histories for each individual pooled across sections for a 21‐day period for each study area in each season. We broke each primary occasion into three 7‐day secondary occasions a priori to maximize detection probability.

In total, five seasons were surveyed over 3 years in the WSA (CD 2012, CD 2013, HD 2013, CD 2014, and HD 2014) and four seasons were surveyed over 2 years in the ESA (CD 2013, HD 2013, CD 2014, and HD 2014). With the rotation design to maximize spatial replication, there were staggered periods of 69–70 days between CD and HD samplings. We considered these periods as open to contribute to survival estimates, in addition to the open periods during wet seasons between HD and CD sampling (CD 2012 – CD 2013 = 349 days, HD2013 ‐ CD2014 = 262 days). One ESA section was not sampled in the HD season of 2013 due to early rains, so our estimates for this season are based on 76% of the camera‐days used for estimate in the other seasons. Survival rates are expressed as annual rates, exponentiating as needed to account for the time over which survival was estimated.

The robust design model includes several parameters (p, c, γ′, and γ”) that must be estimated to provide unbiased estimates of population size and survival rates, but were not of direct interest in this study. To focus on the parameters of interest, we evaluated models in two steps, using Akaike's information criteria corrected for small sample size and overdispersion (QAICc) in the RMark package of R (Laake 2013). To correct for overdispersion, we estimated a median value by collapsing the secondary sessions within each season and fitting a time‐varying Cormack–Jolly–Seber model to the data, as suggested by J. Laake (pers. commun.). All confidence intervals were then corrected for overdispersion in survival rates using this median value as a variance inflation factor. In the first step of model selection, we identified the best model of annual survival out of 10 candidate models (Table 1) with a single estimate of detection probability and no temporary emigration (γ″ = 0 and γ′ = 1). The model(s) receiving the majority of the QAICc weight was selected as the most likely parameterization for annual survival. In the second step, we used the best model(s) for annual survival to test our hypotheses for each of the remaining parameters, resulting in a set of 72 candidate models. In fitting these 72 candidate models, we eliminated any that showed signs of overparametrization. We identified the top models of remaining candidate models (n = 32; Table A1) using QAICc and used model averaging (with the collect.modelsl() and model.average() functions of the RMark package) to estimate parameter values across sex, study area, and time. We calculated overdispersion‐corrected 95% confidence intervals for parameters in each study area using model‐averaged seasonal estimates with pooled variances.

Table 1. Model selection results using QAICc to determine the best‐supported robust design model of survival (S): In the text, this is step one of model selection. Models varied only by their parameterization of S. In all models, there was no temporary emigration (γ″=0, γ′=1) and a single detection probability (p(.)), and population size was estimated by season and study area (N). From these results, we used the top three parameterizations of S for the second stage of model selection Model Parameters Delta QAICc QAICc weight S(.),γ″(0), γ′(1),p(.),N(season+area) 8 0.00 0.50 S(sex), γ″(0), γ′(1),p(.),N(season+area) 9 1.61 0.22 S(area), γ″(0), γ′(1),p(.),N(season+area) 9 2.37 0.15 S(sex+area), γ″(0), γ′(1),p(.),N(season+area) 10 3.99 0.07 S(area*sex), γ″(0), γ′(1),p(.),N(season+area) 11 6.44 0.02 S(time), γ″(0), γ′(1),p(.),N(season+area) 11 6.94 0.02 S(time+sex), γ″(0), γ′(1),p(.),N(season+area) 12 8.56 0.01 S(area+time), γ″(0) γ′(1),p(.),N(season+area) 12 9.34 0.00 S(time+sex+area), γ″(0), γ′(1),p(.),N(season+area) 13 11.00 0.00 S(time+area*sex), γ″(0), γ′(1),p(.),N(season+area) 14 13.62 0.00

We estimated density for each study area by dividing population estimates by the area surveyed. We estimated this area by calculating the mean maximum distance moved (MMDM; Stickel 1954; Wilson and Anderson 1985) across all individuals from both study areas and buffering each trap site by half of the MMDM distance (HMMDM). Balme et al. (2009a) found that without telemetry data, using HMMDM and buffering each camera‐trap site (Silver et al. 2004) was the least biased estimator for leopard density when compared to independent estimates of density from intensive telemetry data. Some recent research questions whether MMDM or HMMDM compares more closely to density estimates derived from telemetry data and spatially explicit capture–recapture models (SECR; Efford 2004) in large felids, and that HMMDM may overestimate density estimates by underestimating space use of individuals (e.g., Tobler and Powell 2013). We chose to follow Balme et al. (2009a)'s leopard‐specific recommendation for density estimates, but we provide density estimates based on average MMDM measures in the Appendix (Table A2). A strength of our study design is that the choice of HMMDM vs. MMDM has no effect on differences in density between the ESA and WSA (the primary interest of this study). In this study, the choice is relevant only for comparisons of our density estimates to those from other studies. We did not implement spatially explicit capture–recapture models (SECR; Efford 2004) because current implementations require assumptions about space‐use distributions that are not likely to be met by leopards, and our sampling was carefully designed to limit differences in the area sampled for the two sites. By sampling in immediately adjacent areas with identical sampling grids and methods, we minimized problems related to estimation of sampling area that can arise during conversion of population size to population density, thus avoiding the primary problem that SECR attempts to address.