West African crocodylians are among the most threatened and least studied crocodylian species globally. Assessing population status and establishing a basis for population monitoring is the highest priority action for this region. Monitoring of crocodiles is influenced by many factors that affect detectability, including environmental variables and individual or population-level wariness. We investigated how these factors affect detectability and counts of the Critically Endangered Mecistops cataphractus and the newly recognized Crocodylus suchus. We implemented 195 repetitive surveys at 38 sites across Côte d’Ivoire between 2014 and 2019. We used an occupancy-based approach and a count-based GLMM analysis to determine the effect of environmental and anthropogenic variables on detection, and modeled crocodile wariness over repetitive surveys. Despite their rarity and level of threat, detection probability of both species was relatively high (0.75 for M. cataphractus and 0.81 for C. suchus), but a minimum of two surveys was required to infer absence of either species with 90% confidence. We found that detection of M. cataphractus was significantly negatively influenced by fishing net encounter rate, while high temperature for the previous 48h of the day of the survey increased C. suchus detection. Precipitation and aquatic vegetation had significant negative and positive influence, respectively, on M. cataphractus counts and showed the opposite effect for C. suchus counts. We also found that fishing encounter rate had a significant negative effect on C. suchus counts. Interestingly, survey repetition did not generally affect wariness for either species, though there was some indication that at least C. suchus was more wary by the fourth replicate. These results are informative for designing future survey and monitoring protocols for these threatened crocodylians in West Africa, and for other endangered crocodylians globally.

Crocodiles surveys_From 2014 to 2019, we conducted surveys of M. cataphractus and C. suchus in different habitat types throughout Côte d’Ivoire. We sampled 38 sites across the major ecoregions of the country, which are representative of those of West Africa: Guinean forest (50% of the country), the Sudano-Guinean zone (19% of the country), and the Sudanian region (31% of the country). We surveyed crocodiles predominantly during the dry season to increase detection rates (Fukuda et al., 2013).  We conducted standard nocturnal spotlight surveys (Chabreck 1966) from an inflatable, outboard-powered boat with 15 hp engine at a cruising speed of about 5.0 – 6.0 km/h, by inflatable kayak, and/or on foot.  A single observer conducted all surveys, who located crocodiles by their eyeshine using either a 78 lumen LED headlamp (80% of observations) or a 550 lumen LED spotlight (20% of observations), depending on the habitat , and approached individuals as close as possible to visually determine species and demographic class (i.e., by total length, TL). We classified crocodiles that submerged before species and total length could be determined as eyes only (EO).  We tracked all survey routes and took waypoints for each crocodile sighting using a handheld GPS. we surveyed each site on three (from 2014 to middle 2016) or five (end of 2016 to 2019) consecutive nights. We surveyed a minimum survey distance of 10 km at each site on each occasion.

Environmental variables_We examined the influence of 10 environmental and anthropogenic variables that were previously shown to have significant influence on crocodile detection probability and are relevant to our study species and habitat.  We measured six of these variables in the field before or during each survey: moon phase (0 – 4), wave height and wind speed (0 – 3), precipitation the day prior to the survey (0, 1), the amount of vegetation present along the shoreline and fishing net encounter rate.  We used a binary index of low or high vegetation where low vegetation denotes a visible shoreline with little to no overhanging vegetation and high aquatic vegetation denotes a shoreline completely covered by overhanging or aquatic vegetation (Gardner et al., 2016).  We counted the number of fishing nets seen on the survey as an index of the subsistence fishing threat (Shirley et al., 2009).  We assessed mean night air temperature and mean daily precipitation both on the day of the survey and for the previous 48 h from remote sensed data accessed through MODIS (Wan et al., 2015) and CHIRPS (Funk et al., 2015), respectively. Prior to further analysis, we standardized all continuous covariates and tested for multi-collinearity among independent variables using the VIF function in the R package car (John et al., 2020).  We found evidence of collinearity for wind speed with wave height, and mean night air temperature the day of the survey with mean night temperature for the previous 48h.  Wave height is often correlated with wind speed (Woodward & Marion, 1978) and generally not significant in the small river systems where we surveyed, so we removed wave height from subsequent models.  Likewise, we retained mean night temperature for the previous 48h over mean night air temperature the day of the survey because of its more significant individual effect in subsequent models (Couturier et al., 2013).  We conducted all subsequent analyses for each species independently.

Influence of environmental and anthropogenic variables on detection probability.¾We assessed the influence of environmental and anthropogenic variables on crocodile detection probability using both an occupancy framework and with linear mixed models.  For both model types we included all data from all surveys in all years, though treated missing data for repetitions four and five in years 2014 to mid 2016 differently.  We categorized missing values as NA in occupancy models, but used imputation methods to infer missing values in GLMM analyses (see below; Nakagawa & Freckleton, 2008).

Within the occupancy framework, we used a single season occupancy model to estimate detection probabilities (p) (MacKenzie et al., 2006).  To do this, we created a detection history (0 = non detection, 1 = detection) for each site across all the survey repetitions.  For this analysis, we hypothesized that the populations were closed during the survey period, no heterogeneity in detection occurred, and the detection process was independent at each site (MacKenzie et al., 2002).  We used the method of “plausible combination” (Bromaghin et al., 2013) for model selection and covariate evaluation, which is increasingly recognized as a robust multi-stage strategy to assess the fit of single season occupancy models (Morin et al., 2020).  To derive detection probability, and better understand the influence of covariates on detection, we paired the most general sub-model for occupancy (ψ) with all candidate sub-models for detection probability (p) using the dredge function in the package MuMin (Barton, 2020) and returned the best model (i.e., ∆AIC threshold of 0) (Morin et al., 2020).  Ultimately, we assessed all combinations of the best detection covariates.  Because our focus was exclusively on detection probability, we held occupancy constant in the final analysis (i.e., (ψ.)p[covariate]) (Kroll et al., 2008; Moreira et al., 2016; Phumanee et al., 2020), a standard practice when focused on one component, or the other, in occupancy-based analyses (Cook et al., 2011; Jeffress et al., 2011; Wagner et al., 2019).  We ranked models using Akaike’s Information Criterion corrected for small sample size (AICc) and considered all models with ∆AICc ≤ 2 to be competitive models (Burnham & Anderson, 2002).  We considered a covariate significant if the 95% CI did not include zero (Bauder et al., 2017).  We conducted all analyses using the packages unmarked (Fiske et al., 2017) and AlCcmodavg (Mazerolle, 2015) in R v4.0.2 (R Development core team, 2020).

For GLMM analysis, as our surveys varied from three to five replicates, we replaced missing values (8.57% of the total dataset) in all sites with less than five replicates using a multiple imputation procedure (Nakagawa & Freckleton, 2008).  Specifically, we used multiple imputation to fill in 18 (of 210; 8.6%) missing values for each of crocodile encounter rate, moon phase, wind speed, amount of aquatic vegetation, fishing net encounter rate, and precipitation the day prior to the survey. We generated and combined 100 imputed datasets (Graham et al., 2007) using the R package mice (Buuren & Groothuis-Oudshoorn, 2011).  After imputation, we determined whether the MICE algorithm has converged by plotting parameters against the iteration number and found no definite trends, indicating good convergence in the dataset including imputed values (Buuren & Groothuis-Oudshoorn, 2011).  We modeled crocodile counts using the lmer function in the lme4 package (Bates et al., 2015).  We tested the same eight covariates included in the occupancy analysis as fixed effects with site as a random effect.  We fit 256 combinations for each species, including the null and global models, without interaction terms and ranked models by AICc.  We considered all models with ∆AICc ≤ 2 to be competitive models (Burnham & Anderson, 2002).  We obtained model-averaged coefficients using the Model.avg function in MuMin.  For each coefficient, we report associate 95% confidence intervals (CIs) and the coefficient estimate with shrinkage (also called “zero method”).  We considered a covariate significant if the 95% CI did not include zero (Bauder et al., 2017).

Effect of repetitive surveys on crocodile wariness.—We used proportion of “EO” and zero detections across repeated surveys as indexes to examine the effect of repetitive surveys on crocodile wariness.  We considered that a progressive reduction in the number of individuals seen or formally identified (i.e., increase in the number of “EO”) during surveys to be a reflection of an increase in wariness.  For each site where the presence of each crocodile species was confirmed, we determined the proportion of observations that were EO (EO/number of all observations) for each survey, where 0.0 represented no wariness and 1.0 (represented by either 100% EO observations or zero detections) represented complete wariness.  For sites where the two studied species were sympatric, we partitioned the EO observations by the ratio of the number of individuals actually attributable to either species.

Because our surveys varied from three to five replicates, we used multiple imputation as described above to estimate missing values for the following covariates: index of wariness, aquatic vegetation, fishing net encounter rate, and crocodile abundance for both species, resulting in 10.9% and 7.14% values derived from MI for M. cataphractus and C. suchus, respectively.  We assessed crocodile wariness as a function of the survey replicate, the most important covariates for each species identified in both occupancy and count-based GLMM analyses above, and crocodile encounter rate as an index of abundance.  We included encounter rate because we hypothesized that population abundance may represent unmeasured disturbance effects on individual wariness (i.e., for rare species, higher abundance sites likely have less impacted or threatened population histories, which may capture unmeasured/unmeasurable histories of disturbance or harassment of individuals).  The M. cataphractus model included survey repetition, abundance, fishing net encounter rate, moon phase, precipitation the day of the survey and for the previous 48h, and aquatic vegetation as fixed factors.  The model for C. suchus included survey repetition, abundance, and temperature for the previous 48h of the day of the survey and aquatic vegetation as fixed factors.  Both species models included site as a random effect.  We used generalized linear mixed effects modelling, fit by restricted maximum likelihood, using the package lme4 (Bates et al., 2015) and Satterthawaite’s approximation from package lmerTest (Kuznetsova et al., 2017) to assess each factor’s significance in R v4.0.2 (R Development core team, 2020).

There are misssing values on the dataset and I included the code used for imputed them.