The nine-banded Armadillo (Dasypus novemcinctus) is the only species of Armadillo in the United States and alters ecosystems by excavating extensive burrows used by many other wildlife species. Relatively little is known about its habitat use or population densities, particularly in developed areas, which may be key to facilitating its range expansion. We evaluated Armadillo occupancy and density in relation to anthropogenic and landcover variables in the Ozark Mountains of Arkansas along an urban to rural gradient. Armadillo detection probability was best predicted by temperature (positively) and precipitation (negatively). Contrary to expectations, occupancy probability of Armadillos was best predicted by slope (negatively) and elevation (positively) rather than any landcover or anthropogenic variables. Armadillo density varied considerably between sites (ranging from a mean of 4.88 – 46.20 Armadillos per km²) but was not associated with any environmental or anthropogenic variables.

Site Selection

Our study took place in Northwest Arkansas, USA, in the greater Fayetteville metropolitan area. We deployed trail cameras (Spypoint Force Dark (Spypoint Inc, Victoriaville, Quebec, Canada) and Browning Strikeforce XD cameras (Browning, Morgan, Utah, USA) over the course of two winter seasons, December 2020-March 2021, and November 2021-March 2022. We sampled 10 study sites in year one, and 12 study sites in year two. All study sites were located in the Ozark Mountains ecoregion in Northwest Arkansas. Sites were all Oak Hickory dominated hardwood forests at similar elevation (213.6 – 541 m). Devils Eyebrow and ONSC are public natural areas managed by the Arkansas Natural heritage Commission (ANHC). Devil’s Den and Hobbs are managed by the Arkansas state park system. Markham Woods (Markham), Ninestone Land Trust (Ninestone) and Forbes, are all privately owned, though Markham has a publicly accessible trail system throughout the property. Lake Sequoyah, Mt. Sequoyah Woods, Kessler Mountain, Lake Fayetteville, and Millsaps Mountain are all city parks and managed by the city of Fayetteville. Lastly, both Weddington and White Rock are natural areas within Ozark National Forest and managed by the U.S. Forest Service.  We sampled 5 sites in both years of the study including Devils Eyebrow, Markham Hill, Sequoyah Woods, Ozark Natural Science Center (ONSC), and Kessler Mountain. We chose our study sites to represent a gradient of human development, based primarily on Anthropogenic noise values (Buxton et al. 2017, Mennitt and Fristrup 2016). We chose open spaces that were large enough to accommodate camera trap research, as well as representing an array of anthropogenic noise values. Since anthropogenic noise is able to permeate into natural areas within the urban interface, introducing human disturbance that may not be detected by other layers such as impervious surface and housing unit density (Buxton et al. 2017), we used dB values for each site as an indicator of the level of urbanization.

Camera Placement

We sampled ten study sites in the first winter of the study. At each of the 10 study sites, we deployed anywhere between 5 and 15 cameras. Larger study areas received more cameras than smaller sites because all cameras were deployed a minimum of 150m between one another. We avoided placing cameras on roads, trails, and water sources to artificially bias wildlife detections. We also avoided placing cameras within 15m of trails to avoid detecting humans.

At each of the 12 study areas we surveyed in the second winter season, we deployed 12 to 30 cameras. At each study site, we used ArcGIS Pro (Esri Inc, Redlands, CA) to delineate the trail systems and then created a 150m buffer on each side of the trail. We then created random points within these buffered areas to decide where to deploy cameras. Each random point had to occur within the buffered areas and be a minimum of 150m from the next nearest camera point, thus the number of cameras at each site varied based upon site size. We placed all cameras within 50m of the random points to ensure that cameras were deployed on safe topography and with a clear field of view, though cameras were not set in locations that would have increased animal detections (game trails, water sources, burrows etc.). Cameras were rotated between sites after 5 or 10 week intervals to allow us to maximize camera locations with a limited number of trail cameras available to us. Sites with more than 25 cameras were active for 5 consecutive weeks while sites with fewer than 25 cameras were active for 10 consecutive weeks. We placed all cameras on trees or tripods 50cm above ground and at least 15m from trails and roads. We set cameras to take a burst of three photos when triggered. We used Timelapse 2.0 software (Greenberg et al. 2019) to extract metadata (date and time) associated with all animal detections. We manually identified all species occurring in photographs and counted the number of individuals present. Because density estimation requires the calculation of detection rates (number of Armadillo detections divided by the total sampling period), we wanted to reduce double counting individuals. Therefore, we grouped photographs of Armadillos into “episodes” of 5 minutes in length to reduce double counting individuals that repeatedly triggered cameras (DeGregorio et al. 2021, Meek et al. 2014). A 5 min threshold is relatively conservative with evidence that even 1-minute episodes adequately reduces double counting (Meek et al. 2014).

Landcover Covariates

To evaluate occupancy and density of Armadillos based on environmental and anthropogenic variables, we used ArcGIS Pro to extract variables from 500m buffers placed around each camera (Table 2). This spatial scale has been shown to hold biological meaning for Armadillos and similarly sized species (DeGregorio et al. 2021, Fidino et al. 2016, Gallo et al. 2017, Magle et al. 2016). At each camera, we extracted elevation, slope, and aspect from the base ArcGIS Pro map. We extracted maximum housing unit density (HUD) using the SILVIS housing layer (Radeloff et al. 2018, Table 2). We extracted anthropogenic noise from the layer created by Mennitt and Fristrup (2016, Buxton et al. 2017, Table 2) and used the “L50” anthropogenic sound level estimate, which was calculated by taking the difference between predicted environmental noise and the calculated noise level. Therefore, we assume that higher levels of L50 sound corresponded to higher human presence and activity (i.e. voices, vehicles, and other sources of anthropogenic noise; Mennitt and Fristrup 2016). We derived the area of developed open landcover, forest area, and distance to forest edge from the 2019 National Land Cover Database (NLDC, Dewitz 2021, Table 2). Developed open landcover refers to open spaces with less than 20% impervious surface such as residential lawns, cemeteries, golf courses, and parks and has been shown to be important for medium-sized mammals (Gallo et al. 2017, Poessel et al. 2012). Forest area was calculated by combing all forest types within the NLCD layer (deciduous forest, mixed forest, coniferous forest), and summarizing the total area (km2) within the 500m buffer. Distance to forest edge was derived by creating a 30m buffer on each side of all forest boundaries and calculating the distance from each camera to the nearest forest edge. We calculated distance to water by combining the waterbody and flowline features in the National Hydrogeography Dataset (U.S. Geological Survey) for the state of Arkansas to capture both permanent and ephemeral water sources that may be important to wildlife. We measured the distance to water and distance to forest edge using the geoprocessing tool “near” in ArcGIS Pro which calculates the Euclidean distance between a point and the nearest feature. We extracted Average Daily Traffic (ADT) from the Arkansas Department of Transportation database (Arkansas GIS Office). The maximum value for ADT was calculated using the Summarize Within tool in ArcGIS Pro.

We tested for correlation between all covariates using a Spearman correlation matrix and removed any variable with correlation greater than 0.6. Pairwise comparisons between distance to roads and HUD and between distance to forest edge and forest area were both correlated above 0.6; therefore, we dropped distance to roads and distance to forest edge from analyses as we predicted that HUD and forest area would have larger biological impacts on our focal species (Kretser et al. 2008).

Occupancy Analysis

In order to better understand habitat associations while accounting for imperfect detection of Armadillos, we used occupancy modeling (Mackenzie et al. 2002). We used a single-species, single-season occupancy model (Mackenzie et al. 2002) even though we had two years of survey data at 5 of the study sites. We chose to do this rather than using a multi-season dynamic occupancy model because most sites were not sampled during both years of the study. Even for sites that were sampled in both years, cameras were not placed in the same locations each year. We therefore combined all sampling into one single-season model and created unique site by year combinations as our sampling locations and we used year as a covariate for analysis to explore changes in occupancy associated with the year of study.

For each sampling location, we created a detection history with 7 day sampling periods, allowing presence/absence data to be recorded at each site for each week of the study. This allowed for 16 survey periods between 01 December 2020, and 11 March 2021 and 22 survey periods between 01 November 2021 and 24 March 2022. We treated each camera as a unique survey site, resulting in a total of 352 sites. Because not all cameras were deployed at the same time and for the same length of time, we used a staggered entry approach.

We used a multi-stage fitting approach in which we used Akaike’s Information Criterion (AIC) to select for the best detection covariate. We modeled the survey period (to allow detection to vary over time), year (to evaluate detection across the two years of the study), weekly mean precipitation (to evaluate if precipitation influenced Armadillo activity and thus detection), and weekly mean temperature (to evaluate if temperature influenced Armadillo activity and thus detection) as covariates for detection against null occupancy parameters and selected the top covariate model with lowest AIC score. The top-ranked detection covariate(s) was then used in all subsequent analyses of occupancy. We acquired temperature and precipitation data from the NOAA weather station closest to each site for each detection date.

For occupancy covariates, we used distance to the nearest water source, distance to the nearest road, elevation, slope, aspect, maximum ADT, maximum anthropogenic noise, developed open space, area of forest, and maximum housing unit density. We then evaluated all single variable models using an AIC approach with an apriori cutoff of 2 ∆AIC (Burnham and Anderson 2002).

Density Estimation

To calculate Armadillo density at each of the study sites, we used the Random Encounter Model (REM). The REM was developed to estimate density of unmarked animals through camera trap data (Rowcliffe et al. 2008). The three assumptions of the REM are: 1) that animals move randomly throughout their environment and thus cameras are not set on any features that might increase their detection probability (e.g., trails, roads, bait etc.), 2) detection episodes are of individual animals, and 3) that the study population is closed (Rowcliffe et al. 2008). We used the REM equation to calculate Armadillo density at each camera using Microsoft Excel (Microsoft corporation). 

In the REM, the y represents the total detections of Armadillo at each camera. Total trap nights in hours (the measure of trapping effort) is represented by t. V refers to the day range (how far an individual travels in a 24-hour period) of the Armadillo. We derived a mean day range for Armadillo from day ranges reported in the literature (Table 3). Detection radius (r) and detection angle (θ) were measured at each camera in the field through walk tests. The walk tests involved walking directly towards each camera to calculate the detection radius and from each side at 5m from the camera to calculate the detection angle in degrees. Detection was determined by whether or not the detection light on the camera was triggered during the walk test. We then converted the detection angle to radians for density calculations (Caravaggi et al. 2015, Rowcliffe et al. 2008). We were not able to calculate the detection angle and radius at 14 of the cameras due to camera malfunction (no detection light during walk test), and so we used the average detection angle for the given camera model (Schaus et al. 2020).

We evaluated if Armadillo density correlated to anthropogenic or landcover variables. We used linear models in R, using the packages “lme4” and “AIC modavg”. We only included data from the second year of sampling in our density calculations as detection radius and angle were not collected during the first year of the study. We evaluated Armadillo density against HUD, anthropogenic noise, distance to water, forest area, development, and ADT. We modeled all single- and two-way combinations of these variables. However, we did not include ADT and HUD in the same models due to high correlation between these covariates. Thus, we evaluated 22 candidate models including the null and global models. In each model, we included site as a random effect. We considered models within 2 ∆AIC to be competitive (Burnham and Anderson 2002).