DOI: 10.5061/dryad.2bvq83bz2
title: Abundance-mediated species interactions between coyote, fisher, and marten in Northeastern US

description: Ecological theory posits that the strength of interspecific interactions
 is fundamentally underpinned by the population sizes of the involved
 species. Nonetheless, contemporary approaches for modelling species
 interactions predominantly centre around occupancy states. Here, we use
 simulations to illuminate the inadequacies of modelling species
 interactions solely as a function of occupancy, as is common practice in
 ecology. We demonstrate erroneous inference into species interactions due
 to bias in parameter estimates when considering species occupancy alone.
 To address this critical issue, we propose, develop, and demonstrate an
 occupancy-abundance model designed explicitly for modelling
 abundance-mediated species interactions involving two or more species.
 When modelling interactions as a function of abundance rather than
 occupancy, we uncover previously unidentified interactions. Through an
 empirical case study and comprehensive simulations, we demonstrate the
 importance of accounting for abundance when modelling species
 interactions, and we present a statistical framework equipped with MCMC
 samplers to achieve this paradigm shift in ecological research.
description: In the case study using empirical data, we examine the intraguild
 interactions between three carnivores, a top mesopredator in the system,
 the coyote, an intermediate mesopredator, the fisher<em>,
 </em>and a small carnivore, the American
 marten<em>.</em> There is a long history of examining
 intraguild interactions between fisher and marten through harvest (e.g.,
 Hardy, 1907; Krohn, Zielinski, & Boone, 1997). Recent
 harvest-based evidence was used to infer negative interactions between all
 three species, with fishers being limited through intraguild killing by
 coyotes, and martens being limited by both fisher and coyotes (Jensen
 & Humphries, 2019). Nonetheless, the three species co-occur over
 much of the marten’s limited range in New York State and recent analysis
 using Rota <em>et al.</em> (2016) co-occurrence models was
 inconsistent with previous hypotheses. This analysis found fisher
 occupancy was higher conditional on coyote presence, and marten occurred
 independently from both other species (Twining <em>et al</em>.
 In Press). Nonetheless, as explored in Simulation study I, a focus on
 occupancy states (and ignoring the abundance of species) to infer
 interactions may limit inference. As an example, we fit the
 occupancy-abundance model to the landscape scale camera trapping dataset
 on these species previously analysed using co-occurrence models in Twining
 <em>et al.</em> (In Press). Detection frequency data on
 coyotes and fisher, and detection-non detection data on American marten
 was obtained from camera traps deployed to monitor occurrence of the
 target species in northeastern New York State. Sampling was conducted from
 January-March 2016-2018 throughout the Adirondack and Tug Hill regions of
 northern New York State. This sampling used a stratified random sampling
 design to select 195 15 km<sup>2</sup> sample units across the
 region of interest. A standardized methodology was used across all
 surveys. At each site, a camera trap was deployed randomly within the 15
 km<sup>2 </sup>grid. Camera traps were secured to trees
 approximately 1.0-1.5 m above ground. A bait station was placed on a tree
 opposite the camera trap and secured to the tree using wire mesh. At all
 sites skunk-based call lures were applied. Cameras were deployed for 3
 weeks (21 days) at each location after which cameras were retrieved.
 Cameras and bait were checked halfway through 3 weeks of sampling with
 batteries and bait replaced and replenished as necessary<em>.
 </em>We used a weekly occasion length. For the
 detection/non-detection data only one detection of the subordinate species
 was possible per weekly period. For the count data, we allowed a single
 detection each 24-hr period over each week and summed the days with
 detections into weekly counts for the species. The same sample units were
 sampled in each of the three sampling years (except for 13 sites that were
 not sampled in 2017).
description: # Abundance-mediated species interactions # Summary These are MCMC
 samplers, data simulators, and processing and run scripts for a
 occupancy-abundance model for abundance-mediated species interactions with
 a full example case study. This model framework models
 detection/non-detection data to estimated occupancy, abundance, and
 interactions between the species of interest. This model enables the user
 to apply summarized detection/non-detection data collected over repeat
 surveys to model interactions as a function of abundance. These samplers
 are presented in Twining et al. 2024, and are based on adaptations of the
 [Waddle et al. (2010)](https://www.jstor.org/stable/25680391) formulation
 for modelling species interactions within an occupancy model, but instead
 of modelling the state model of a subordinate species a function of the
 occupancy states, it is modelled as a function of abundance (N). We
 provide a range of MCMC samplers for different ecological scenarios
 between two or more species including disease- and predator- mediated
 competition, intraguild predation, mesopredator release, and tri-trophic
 cascades. # **The working directory** Below you fill find descriptions of
 each folder in this repository and files contained within them. # Models
 and data simulators
 (./occupancy\_abundance\_model/models\_and\_simulations) This folder has
 seven R scripts. Within each script are MCMC samplers and a data simulator
 for different iterations of the occupancy-abundance model. Model code and
 implentation is conducted in the nimble package [(de Valpine et al.
 2022).](https://cran.r-project.org/web/packages/nimble/index.html) **1.
 two_species_occupancy_abundance_data_simulator_and_model.R** This script
 contains a data simulator and model for simulating a two species occupancy
 abundance model with spatially varying interaction terms. The code is
 commented throughout to describe each part of the simulator and script.
 **2.three_species_occupancy_abundance_abundance_model_and_simulator_pinemarten_greysquirrel_redsquirrel_example.R** This script contains a data simulator and model for simulating a three species occupancy abundance model - modeling the occupancy of the dominant species on abundance of intermediate, and abundance of subordinate, with additional interactions between abundance of intermediae with the subordinate. The code is commented throughout to describe each part of the simulator and script. **3.three_species_abundance_abundance_abundance_model_and_simulator_otter_urchin_kelp_example.R** This script contains a data simulator and model for simulating a three species abundance model - modeling the abundance of the dominant species on abundance of intermediate, and abundnace of subordinate, with additional interactions between abundance of intermediae with the subordinate. The code is commented throughout to describe each part of the simulator and script. **4. three_species_abundance_occupancy_occupancy_model_and_simulator_deer_liverfluke_moose_example.R** This script contains a data simulator and model for simulating a three species abundance model - modeling the abundance of the dominant species on occupancy of intermediate, and occupancy of subordinate, with additional interactions between occupancy of intermediae with the subordinate. The code is commented throughout to describe each part of the simulator and script. **5. occupancy_vs_abundance_mediated_interactions_model_and_simulator.R** This script contains a data simulator and model for simulating a two species occupancy abundance model with and fits models with interactions as a function of occupancy and as a function of abundance. The code is commented throughout to describe each part of the simulator and script. **6. occupancy_abundance_model_three_species_simulator_and_model_variable_detection.R** This script contains a data simulator and model for simulating a three species abundance model - modeling the abundance of the dominant species on abundance of intermediate, and occupancy of subordinate, with additional interactions between abundance of intermediae with the subordinate. The code is commented throughout to describe each part of the simulator and script. **7. two_species_occupancy_abundance_simulator_spatially_varying_interaction_terms** This script contains a data simulator and model for simulating a two species occupancy abundance model with spatially varying interaction terms. The code is commented throughout to describe each part of the simulator and script. # The coyote-fisher-marten case study folder (./occupancy\_abundance\_model/case\_study) This folder has nine files. **1. alNY_2013-2021_6kmbuffer_allsitecovs.csv** This file contains all of the summarized spatial covariate data at a 6km scale used in the analysis. Each row is a different 6km 2 pixel in New York State, each column is a covariate, and each cell is a value. NAs indicate that data was not avaialble for that grid cell.  The covariates used in the analysis. | **Covariate** | **Description** | **Source** | | :--------------- | :------------------------------------------------------------------------------------------- | :----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- | | Deciduous | Proportion of a 6 km2 buffer around the detector made up of deciduous forest | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Coniferous | Proportion of a 6 km2 buffer around the detector made up of coniferous forest | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Mixed | Proportion of a 6 km2 buffer around the detector made up of mixed forest | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Pasture | Proportion of a 6 km2 buffer around the detector made up of pasture | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Cultivated.Crops | Proportion of a 6 km2 buffer around the detector made up of cultivated crops | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | road\_density | Mean number of km of road per km2 in each 6 km2 buffer | calculated from primary and secondary roads raster provided by the NYSDEC, hosted on the github | | snow\_depth | Mean daily snow depth(m) of the 6 km2 buffer around each detector across the sampling period | National Operational Hydrologic Remote Sensing Centre, 2004. Snow data assimilation system (SNODAS) products ([https://doi.org/10.7265/N5TB14TC](https://doi.org/10.7265/N5TB14TC)). | | forest\_edge | Edge density of combined class of all forest in each 6 km2 buffer around detectors | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | GPP | The amount of carbon captured by plants (kg C MJ-1) in at the detectors | MODIS Land Satellite, 2017 ([https://lpdaac.usgs.gov/products/mod17a2hv006/](https://lpdaac.usgs.gov/products/mod17a2hv006/)) | | Camera | The camera trap model used at each site (Bushnell, Recoynx, Browning) | Fieldwork datasheets | **2.allNY__2013-2021_15kmbuffer_allsitecovs.csv** This file contains all of the summarized spatial covariate data at a 15km scale used in the analysis. Each row is a different 15km 2 pixel in New York State, each column is a covariate, and each cell is a value. NAs indicate that data was not avaialble for that grid cell.  The covariates used in the analysis. | **Covariate** | **Description** | **Source** | | :--------------- | :-------------------------------------------------------------------------------------------- | :----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- | | Deciduous | Proportion of a 15 km2 buffer around the detector made up of deciduous forest | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Coniferous | Proportion of a 15 km2 buffer around the detector made up of coniferous forest | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Mixed | Proportion of a 15 km2 buffer around the detector made up of mixed forest | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Pasture | Proportion of a 15 km2 buffer around the detector made up of pasture | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | Cultivated.Crops | Proportion of a 15 km2 buffer around the detector made up of cultivated crops | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | road\_density | Mean number of km of road per km2 in each 15 km2 buffer | calculated from primary and secondary roads raster provided by the NYSDEC, hosted on the githu | | snow\_depth | Mean daily snow depth(m) of the 15 km2 buffer around each detector across the sampling period | National Operational Hydrologic Remote Sensing Centre, 2004. Snow data assimilation system (SNODAS) products ([https://doi.org/10.7265/N5TB14TC](https://doi.org/10.7265/N5TB14TC)). | | forest\_edge | Edge density of combined class of all forest in each 15 km2 buffer around detectors | NLCD, 2019 ([https://www.mrlc.gov/data/nlcd-2019-land-cover-conus](https://www.mrlc.gov/data/nlcd-2019-land-cover-conus)) | | GPP | The amount of carbon captured by plants (kg C MJ-1) in at the detectors | MODIS Land Satellite, 2017 ([https://lpdaac.usgs.gov/products/mod17a2hv006/](https://lpdaac.usgs.gov/products/mod17a2hv006/)) | | Deer | The probability of occupancy (ψ) of white-tailed deer in a 15 km2 buffer around each detector | Calculated from detection/non-detection data and covariate hosted in the repository (see below) | | Camera | The camera trap model used at each site (Bushnell, Recoynx, Browning) | Fieldwork datasheets | **3. 7daycountDetectionsmartenNZ.csv** This file contains the count data for American marten between the years of 2016-2018 in New York state. Each row is a site, each column is an occasion, each cell is a summary of 7 daily detection/non-detection record. NAs indicate that the site was not sampled on that occasion. .  **4. allNY_2013_2021_7dayocc_coyote_counts.csv** This file contains the count data for coyote between the years of 2013-2021 in New York state. Each row is a site, each column is an occasion, each cell is a detection/non-detection record. NAs indicate that the site was not sampled on that occasion.  **5. allNY_2013-2021_7dayocc_fisher_counts.csv** This file contains the count data for fisher between the years of 2013-2021 in New York state. Each row is a site, each column is an occasion, each cell is a detection/non-detection record. NAs indicate that the site was not sampled on that occasion.  **6. 7dayoccDetectionswhitetaileddeerNZ.csv** This file contains the detection/non-detection data for white-tailed deer between the years of 2016-2018 in the northern zone of New York state. Each row is a site, each column is an occasion, each cell is a detection/non-detection record. NAs indicate that the site was not sampled on that occasion.  **7. juliandays_allNY_2013-2021.csv** This file contains the sampling dates (in ordinal format) for each sampling occasion for the winter surveys from 2013-2021. Each row is a site, each column is an occasion. NAs indicate that the site was not sampled on that occasion.  **8. 'coyote_fisher_marten_abu_ocu_binomialobsmodel_nimble_model.R'** This is the nimble occupancy-abundance model that is fit to the coyote-fisher-marten data files above. The code is commented out to describe each part of the model. **9. nimble_occu_abu_coyote_fisher_marten_binomialversion_processingandruncode.R** This is the data loading, formatting, and run script for the case study three species occupancy-abundance model. The code is commented out to describe each stage of the process. **Updates log** This repository was updated during revisions of associated manuscript due to amendments to the analysis (changes to code to use 7-day summarized detection/non-detection data data vs. weekly detection/non-detection data. The data used in the revised analysis was added).