Dataset archived in Dryad: https://doi.org/10.5061/dryad.9w0vt4bq9

This repository contains the R code and input data required to run the Spatially-Explicit Individual Based Model (SE-IBM) presented in the article "Resident mortality determines grey wolf range dynamics". The model was developed for the population of wolves in the Iberian Peninsula to investigate the demographic drivers of range dynamics, particularly how wolf mortality influences demographic processes affecting population growth and spatial expansion.

---

Repository contents

wolf_IBM.zip contains the following files and folders:

wolf_IBM.zip
│── README.md
│── runIBM.R
│── paramIBM.R
│── submodelsIBM.R
│── exceedNG_91
│── exceedNG_02
│── exceedNG_13
│
├── groupLoc/
│   ├── wolves91.RData
│   ├── wolves02.RData
│   ├── wolves13.RData
│   ├── wolves21.RData
│   ├── land91.RData
│   ├── land02.RData
│   ├── land13.RData
│   ├── land21.RData
│   ├── relInPop91.RData
│   ├── relInPop02.RData
│   ├── relInPop13.RData
│   └── relInPop21.RData
│
├── outputs/
│   ├── outputs91/
│   ├── outputs02/
│   └── outputs13/
│
└── maps_geotiff/
    ├── same_2030.tiff
    ├── same_2040.tiff
    ├── same_2050.tiff
    ├── red10_2030.tiff
    ├── red10_2040.tiff
    ├── red10_2050.tiff
    ├── red10r_2030.tiff
    ├── red10r_2040.tiff
    ├── red10r_2050.tiff
    ├── red10d_2030.tiff
    ├── red10d_2040.tiff
    └── red10d_2050.tiff

Recommended software

The repository can be used with the following free/open-source software:

• R (v4.2.1): to run the model code and reproduce simulations.
• RStudio (optional): recommended interface for running R scripts.
• QGIS (v3.22): to visualise GeoTIFF spatial outputs.
• Markdown editor (optional): recommended to inspect the README file with proper formatting.

---

Code and reproducibility

The code in this repository runs the SE-IBM used to simulate wolf population dynamics in the Iberian Peninsula. Specifically, it reproduces the simulations used to calibrate the model through Pattern-Oriented Modelling (POM) during past study periods (1991–2002, 2002–2013, and 2013–2021), generating the raw outputs from which demographic metrics reported in the associated article were derived.

The file runIBM.R runs 13,500 simulations for each past study period, corresponding to 135 parameterizations and 100 replicates per parameterization.

The file runIBM.R calls the following scripts and auxiliary files:

• paramIBM.R, which loads the initial wolf population (wolves91.RData, wolves02.RData, wolves13.RData, wolves21.RData), the landscape layers containing habitat quality and the initial location of group territories (land91.RData, land02.RData, land13.RData, land21.RData), and sets up model parameterizations, including loading relatedness tables (relInPop91.RData, relInPop02.RData, relInPop13.RData, relInPop21.RData).
• submodelsIBM.R, which contains the sub-models simulating survival, reproduction and dispersal processes.
• exceedNG_91, exceedNG_02, and exceedNG_13, which are used during calibration of past study periods to stop simulations that produced unrealistically large wolf populations.

Output files are generated automatically when runIBM.R is executed and saved in folders outputs/outputs91, outputs/outputs02, and outputs/outputs13, respectively.

The repository also includes the input files required to initialize the 2021 wolf population (wolves21.RData, land21.RData, and relInPop21.RData) for future projections. Users can project wolf dynamics for 2021–2050 by modifying paramIBM.R to run the parameterizations retained through POM under four mortality scenarios: no changes in mortality; a 10% reduction in overall mortality; a 10% reduction only for dispersers; and a 10% reduction only for residents (see Data S1 and Table S5 of the associated article). To save simulation outputs, users must create the folder outputs/outputs21 and modify the output path accordingly.

The folder maps_geotiff contains GeoTIFF rasters corresponding to projected distributions of wolf groups in the Iberian Peninsula for 2021–2050, obtained from the SE-IBM. These files correspond to the projection scenarios presented in Figure 3 of the associated article and are the only simulation outputs included in the repository; all other outputs can be generated using the code provided here.

Detailed comments describing model processes and workflow are included within each script.

---

Description of the data and file structure

Initial wolf population layers

wolves91.RData, wolves02.RData, wolves13.RData, wolves21.RData

These files contain the initial wolf population used to launch the SE-IBM for the study periods 1991–2002, 2002–2013, 2013–2021, and 2021–2050, respectively. Each file is stored as a SpatialPointsDataFrame object and contains the spatial location of individuals (in map units) and their characteristics. Individuals belonging to the same wolf group share the same coordinates, representing the group location.

The SpatialPointsDataFrame attribute table contains the following individual-level variables:

Variables and definitions:

• who: individual identity. Each wolf has a unique numerical identity.
• disp: status of the individual as resident or disperser. Categorical variable with three levels: 0 corresponds to group members (“residents”), 1 to dispersers on the move (“seekers”), and 2 to single dispersing females that restricted her movements to a suitable territory but have not yet formed a group (“receivers”). Dispersing males are always classified as seekers.
• packID: territory identity for residents and receivers. Each group and receiver hold a territory with a unique numerical identity.
• dominant: dominance status. Categorical variable with two levels: dominant (1) or not (0). Despite the dominance status may be meaningless for dispersers, we assigned value 1 to receivers for coding purposes.
• sex: female (F) or male (M).
• age: age in calendar years. 1st calendar year corresponds to pups (2–11 months), 2nd calendar year to yearlings (12–23 months), 3rd calendar year to 2-year-old wolves (24–35 months), etc. Individuals in their 1st, 2nd and 3rd or upper calendar year are referred to as pups, yearlings and adults respectively.
• heading: direction during last movement for seekers. Eight possible directions in degrees from 0º to 360º (0º is heading North) and every 45º.
• firstB: breeding status for dominant pairs. Categorical variable with two levels: the dominant pair already produced pups together (meaning that pups emerged from the den; 1) or not (0).
• dterr: number of days as resident for residents in groups.
• ddom: number of days as dominant for dominants in groups.
• dpaired: number of days paired for dominants in groups.
• daughter: information for dominant females. Categorical variable with two levels: the dominant female is the daughter of the previous dominant female (1) or not (0). We assign 0 to receivers.
• MotherID: mother’s identity.
• FatherID: father’s identity.

---

Landscape layers

land91.RData, land02.RData, land13.RData, land21.RData

These files contain the landscape layers used to initialize the SE-IBM for the study periods 1991–2002, 2002–2013, 2013–2021, and 2021–2050, respectively. Each file is stored as a RasterBrick object composed of two raster layers representing habitat quality and the initial location of group territories in the Iberian Peninsula.

The study area is divided into grid cells of 100 km² resolution. Each raster contains two layers:

• layer.1: Habitat layer. We used the map of habitat quality from Grilo et al. (2019). Cell values represent the predicted probability of wolf occurrence based on environmental variables related to topography, water availability, vegetation, prey availability, and human disturbance. The original habitat layer (ETRS89 / UTM zone 30N, EPSG:25830) was transformed to map units. This is the spatial reference for locating individuals, territories and movements. We assigned zero value to the sea cells bordering the Iberian Peninsula. We did not provide information on habitat quality to the cells corresponding to France.
• layer.2: Territory layer. Cells occupied by a group or a receiver hold the identity number of the territory (packID). The rest of cells within the Iberian Peninsula have value 0 meaning that such areas are not occupied by groups or receivers. We do not assign information to cells bordering the Iberian Peninsula.

---

Relatedness tables

relInPop91.RData, relInPop02.RData, relInPop13.RData, relInPop21.RData

These files contain auxiliary tables used to identify related wolves in the initial population for the study periods 1991–2002, 2002–2013, 2013–2021, and 2021–2050, respectively. Each file is stored as a data frame object and contains the identity of each wolf individual present in the initial population (who) and the identity of the wolf group to which it belonged (packID). In the model, individuals from the initial population are considered related if they belonged to the same wolf group.

These files are used by the relatedness function in the model code. For pairs of individuals from the initial population, relatedness is assigned as 1 if both individuals belonged to the same group and 0 otherwise. For individuals born during simulations, relatedness is additionally evaluated using parental identities stored in the population object.

---

Calibration files

exceedNG_91, exceedNG_02, exceedNG_13

These files contain auxiliary objects used during model calibration to track parameterizations that produced unrealistically large wolf populations during simulations for the study periods 1991–2002, 2002–2013, and 2013–2021, respectively. Each file is stored as a data frame. For each parameterization (sample.n), the files record the number of simulation replicates (numExceed) that exceeded a predefined threshold in the number of wolf groups at any simulated day. We set this threshold to twice the maximum number of wolf groups documented in the Iberian Peninsula during 1991–2021 (i.e., twice the maximum documented value of 356 groups).

These files are used internally by runIBM.R to stop running parameterizations that consistently generated unrealistic population sizes, thereby improving computational efficiency during model calibration.

---

Output folders

outputs/outputs91, outputs/outputs02, outputs/outputs13

This folder contains empty subfolders (outputs91, outputs02, and outputs13), included to preserve the original workflow of the model code and used to store outputs generated when running the SE-IBM for the study periods 1991–2002, 2002–2013, and 2013–2021, respectively.

Executing runIBM.R generates the following outputs for each simulation:

• popSimRep: a list of yearly SpatialPointsDataFrame objects storing the state of the wolf population at each simulated year (March 10th), including individual-level demographic information and wolf locations. Individuals belonging to the same wolf group share the same coordinates and group identity (packID).

• mapSimRep: a list of yearly RasterLayer objects representing wolf group territories at each simulated year (March 10th). Raster cell values correspond to group identities (packID), allowing linkage with popSimRep.

• packDyn: a data.frame containing annual parameters related to changes over time in the number of groups and individuals across different demographic stages (e.g., pups, residents, dispersers, dominant individuals).

• dispTable: a data.frame containing distance, duration, outcome, and other characteristics associated with each dispersal event.

These objects were used to derive the population dynamic metrics reported in the associated article. Of the variables available in packDyn and dispTable, we only used those described in Table S3 and Table S4, respectively. The remaining variables are available for further studies on population and dispersal dynamics. To obtain the metrics, we excluded the first three years of simulation data, which could be biased by the model’s initial conditions (see Supplementary Methods).

---

Spatial outputs

maps_geotiff

The folder maps_geotiff contains GeoTIFF rasters corresponding to projected distributions of wolf groups in the Iberian Peninsula for 2021–2050 obtained from the SE-IBM. Each raster cell stores the proportion of simulations (out of 300) predicting the presence of a wolf group under a given mortality scenario and projection year. Raster values range from 0 to 1, where 0 indicates that no simulations predicted wolf group presence and 1 indicates that all simulations predicted presence. The coordinate reference system (CRS) is ETRS89 / UTM zone 30N (EPSG:25830).

Files follow the naming structure: [scenario]_[year].tiff

Mortality scenarios:

• same → no changes in mortality
• red10 → 10% reduction in overall mortality
• red10r → 10% reduction in resident mortality
• red10d → 10% reduction in disperser mortality

Projection years:

• 2030 → 10 years after the 2021 baseline
• 2040 → 20 years after the 2021 baseline
• 2050 → 30 years after the 2021 baseline

These are the only simulation outputs included in the repository and were derived from mapSimRep object.

---

Version note: This second version includes minor edits to the README only.