INTRODUCTION
Anthropogenic disturbances such as logging and forest clearing constitute major threats to biodiversity (Michalski and Peres 2005). When a disturbance changes environmental conditions, the natural community can suffer profound, immediate, and cumulative effects (Allesina et al. 2006). However, as different species have different ecological attributes-their scale of movement, life history, longevity, and what constitutes its habitat-that determine their niche and their capacity to survive in a modified landscape. Thus, some become more abundant or expand their range, while others may drop and even become locally extinct (Boyle and Smith 2010, Chapman et al. 2010). The loss of a species or a change in its abundance, particularly for those that interact with many other species, can have a considerable effect on ecological processes (Lindenmayer and Fischer 2013). It is, therefore, important to consider multiple anthropogenic factors when assessing the synergistic effects of human activities, in order to understand how species, respond in fragmented landscapes (Benchimol and Peres 2013).
The west margin of the Magdalena River valley and east slope of the Cordillera Central have a history around 150 years of intensive transformation, and the natural land cover has been reduced by almost 70 %, transforming it into a highly fragmented landscape (IDEAM 2009). The floral association in the forest remnants have been impacted at different levels and different ways, affecting successional processes.
Four primates are found in this region: Red howler monkeys (Alouatta seniculus Linnaeus, 1766), brown spider monkeys (Ateles hybridus Geoffroy, 1829), white fronted capuchins (Cebus versicolor Pucheran, 1845) and white footed tamarins (Saguinus leucopus GUnther, 1877). A. seniculus is an atelid that weights around 7 kg that inhabits a broad range of forest covers from transitional vegetation until riparian and mature forest up to 3200 masl. It is a facultative leaf-eater monkey when the fruit is scarce (Link et al. c2021c). A. hybridus is an atelid primate too around that weights 7.5 kg. It is dependent of mature forest, but some groups could be found in fragmented forests. It is highly fru-givorous and feed largely on ripe fleshy fruits, which comprise approximately 80 % of its diet (Link et al. c2020b). C. versicolor is a cebid primate that weights around 2.5 kg that inhabits tropical lowland forests and seasonally flooded forests. It has an omnivorous diet, eating predominantly ripe fruits and invertebrates, but also other vegetative parts of plants, vertebrate prey and eggs (Link et al. c2021d). S. leucopus is a small (460 g) callitrichid that inhabits from transitional to mature, mainly humid, and subhumid forests. Its diet is based on fruits and is supplemented with invertebrates, bark, and flowers (Roncancio et al. 2012). and S. leucopus are endemic to Colombia. C. versicolor is Endangered (Link et al. c2021d) and S. leucopus is Vulnerable (Link et al. c2021a). A. hybridus is categorized as Critically Endangered (Link et al. c2020b). A. seniculus is categorized as Least Concern even though its populations are decreasing (Link et al. c2021c). The main threat to these taxa is the reduction and fragmentation of their habitat resulting from the combined effects of selective extraction, hunting and trade, and the lack of no-take protected areas as National Natural Parks (Londoño and Gómez-Posada 2010, Roncancio et al. 2011, Defler et al. 2013).
The responses of these primates to disturbances, and their specific tolerance thresholds to the resulting landscape configuration are, however, not well understood (French and Smith 2005, Ewers and Didham 2007) due to the lack, among other things, of a landscape model grounded on a robust mathematical approach (Reckhow 1990, Wade 2001, Fahrig 2003, Arroyo-Rodríguez et al. 2013). In this study, I estimated the relationship of different traits such as size, shape, isolation distance and land cover composition, plant structure, and plant diversity, with the population density of these four primates. I also estimated the relationship between the population densities of these primates under a theoretical framework of compensation density (Peres and Dolman 2000). I assumed that A. hybridus is a better competitor than other species due to its large size and specialized frugivorous diet (Defler 2010). My hypothesis was that the changes in the population density of primates could be explained by changes in size, shape, isolation distance, composition of plant covers and plant structure and diversity, although, differentially in type and intensity and affected by the interspecific interactions of the assembly.
MATERIALS AND METHODS
Study area
I conducted fieldwork in 20 localities (local landscape), randomly selected with respect to the primate presence. In other words, the presence of the primates was not a criterion to select the localities (WGS 84 - decimal degrees: Amaní -74.88611W, 5.56083N; Berlin -74.57354W, 6.52191N; Cerro Las Campanas -74.45416W, 8.23381N; Cerros de Corcovado -74.46314W, 8.46191N; Charca de Guarinocito -74.73472W, 5.33972N; EL Pajuil -74.69466W, 5.88033N; Falan -74.97000W, 5.11889N; Fresno -75.01475W, 5.20222N; Guamo 2 -74.69180W, 5.86057N; Guamo1 -74.70102W, 5.86893N; La Florida -74.19829W, 8.57968N; La Pedrera -74.95722, 5.57611; La Primavera -74.96667W, 5.51667N; Las Margaritas -74.59761W, 6.56341N; Mariquita -74.94606, 5.22742N; Playas -74.92331W, 6.26989N; Porce III -75.18909W, 6.83280N; Río Manso -74.77361W, 5.67556N; Venecia -74.83861W, 5.64278N; Yalí -74.87281, 6.67851N). These areas are located in the middle Magdalena River valley in Colombia, distributed from the north of the department of Tolima to the south of the department of Bolivar, extending through the eastern parts of the departments of Caldas and Antioquia (Fig. 1). The altitude of the locations ranged from 190 to almost 1400 m (Table 1), with mean temperatures of 19 to 27 °C and mean annual precipitation ranging from 1800 to 4200 mm (Hijmans et al. 2005). The region is approximately 48 000 km2, comprising crops (23 %), pasture (22 %) and natural vegetation (32 %) (the remaining 23 % is for other land use). Natural vegetation includes dense forest fragments of less than 0.1 to 300 and 6000 km2 (mean = 60 km2 ±1640), with shapes that tend to be irregular (Mean Shape Index = 1.8) and medium connectivity (interspersion juxtaposition index = 52, range: 0 < IJI ≤100) (Statistics estimated for author from IDEAM 2009).
Primate population density
Sampling was done between 2005 and 2011. To estimate primate population densities, I used distance sampling along linear transects (Buckland et al. 2001) and censuses in two localities where I could assume that complete counts were possible. With this method, transects can be resampled to increase the sample size to calculate encounter rate, and individual recognition is not necessary. The number and length of transects varied depending on the extension, shape and topography of each site (Table 2). Transects were distributed throughout each site to obtain a representative density estimation. In two localities, Guarinocito and Guamo 1, population density was estimated using censuses, because these localities are fragments without structural connectivity of forest, of 0.21 and 0.1 km2, respectively, and are elongated. I therefore, assumed that it was possible to evaluate the entire area of interest and count all its primates individually. Population densities of A. seniculus, A. hybridus and C. versicolor in Florida, A. seniculus in Pedrera and Rio Manso, and A. hybridus and C. versicolor in Margaritas, were estimated using the detection function obtained in similar localities of this study, because the number of detections were too low (<10) in these localities to enable a robust analysis in Distance. Some population densities were taken from secondary information (Santamaria et al. 2007, Rojas-V. 2009). I assumed that the results are comparable given that they used de same methodology.
Local features and landscape and class metrics of each site
Initially, I defined a local landscape for each locality through a circle that was the smallest perimeter containing the transects used to estimate primate population densities (Fig. 2). For Guarinocito and Guamo 1, the entire fragments were contained in the smallest perimeter. This circle was overlaid with a raster of land cover (IDEAM 2009), a raster of mean annual temperature and mean annual precipitation from BIOCLIM (Hijmans et al. 2005), and with a Digital Elevation Model (USGS 2000). For each locality I measured three local features, two landscape metrics and twelve class metrics (in the class where the population density was estimated). Local features included: mean altitude, mean annual temperature, and mean annual precipitation. Landscape metrics included: Shannon index (SI) and a weighted class index (WCI). I developed the WCI to integrate an index that allowed to measure the land cover composition and structure, not only the land cover structure as SI (Fahrig et al. 2011). To estimate the WCI, I qualified class types as follows: Dense forest = 8, Open forest = 7, Fragmented forest = 6, Scrubland = 5, Grassland (natural) = 4, Crops = 3, Pasture = 2, and wetland, open soil and urbanized area = 1. Then, I multiplied the value of each class type by its relative size in the local landscape. Finally, I added the results of the multiplication to obtain the index values. In this way, the value for a local landscape with predominantly dense forest, was close to eight, while for one with more pasture and crops, the values were between two and three. Class metrics, for the class in the local landscape where the population density estimation was done, were as follows: mean core area (MCA), total core area (TCA), total core area index (TCAI), core area density (CAD), interspersion juxtaposition index (IJI), area weighted mean shape index (AWMSI), mean shape index (MSI), mean patch fractal dimension (MPFD), area weighted mean patch fractal dimension (AWMPFD), total edge (TE), mean patch size (MPS) and core area (CA). All analyses were done using ArcGis 10.5, Patch Grid tool.
Plant structure and diversity
In the forest where the primate densities were estimated of each locality, six or seven rectangular plots of 0.20 ha (50 m x 4 m) were randomly placed. I stopped sampling when the standard deviation of the diameter breast height (DBH) was stable. I took samples for species identification, and determinations were carried out in the herbaria of the Antioquia -HUA-, and Caldas -FAUC- universities. Mean DBH and height, tree density, and the Inverse Simpson Index (ISI) were estimated for each locality as explanatory variables.
Data Analyses
Population density. I analysed data using Distance 7.1 software (Thomas et al. 2010). Analyses were performed separately for each locality. Detection functions were selected according to the fit between the frequency distribution of detection distances and theoretical models (Key functions and adjustment terms). The models tested were: Half-normal (Cosine, Hermite polynomial), Uniform (Cosine, simple polynomial), and Hazard rate (Cosine, simple polynomial). I chose the best fitting model according to the Akaike Information Criterion. Population density variance (CV) was empirically calculated as the sum of the sampling variances of encounter rate, estimates of detection probability, and group size. I used the expected group size estimated from the regression between group size and detection probability to estimate the density (Buckland et al. 2001). In Guarinocito and Guamo 1 where I used censuses, the density was estimated divided the number of individuals over the size of the fragments. Distance sampling is one of the most often used methods to obtain precise and unbiased estimations with a given sampling effort, so I assumed that the results obtained with distance sampling and the censuses are comparable (Cassey and McArdle 1999, Norvell et al. 2003).
Relationship of the explanatory variables in primate density. Because the vegetation was sampled in only fifteen out of the 20 localities, I analysed the relationship of the density of primates and the local features, landscape, and class metrics of the localities independent of the plant structure and diversity. First, I performed analyses of association between the explanatory variables in order to eliminate collinearity (Autocorrelation) (Cush-man et al. 2008). I used the Spearman coefficient of correlation to evaluate this association (data did not distribute normally, Shapiro-Wilk test: W = 0.90, P = 0.04). When two variables were correlated (Spearman rank correlation coefficient: P-rs < 0.05), one of them was eliminated (Sokal and Rohlf 1995), considering the precision of the measurement. However, the main criterion was mathematical. Finally, four out of the 17 local features, landscape and class metrics were considered (altitude, interspersion juxtaposition index, mean shape index and weighted cover index) and, likewise, only the density of the Variegated spider monkey out the four primates and the four plant structure and diversity variables were considered as explanatory variables.
Given that the population density is a variable expressed in area, I assumed that it has Poisson probability distribution and, thus, the relationships were evaluated using Poisson regressions. This was done using a Bayesian approach and using an uninformative prior distribution to the precision of the explanatory variable effects, and an uninformative prior distribution to the intercept (alpha) and slopes (beta). I assumed that the posterior distributions of the intercept and the slope of each variable had a normal distribution (McCarthy 2007, Pfeiffer et al. 2008). I did not consider the interaction (multiplicative effects) between the explanatory variables. In order to select the best model, the Deviance Information Criterion (DIC) was used. Initially, I ran the model with all of the possible explanatory variables and estimated the DIC. Later, I extracted one variable, ran the model and estimated de DIC again. If the DIC was higher, I returned the variable and extracted another. If the DIC resulted lower, I removed definitely that variable and extracted another, until I got the model with the lowest DIC. Estimation of the intercept and slopes was done using Markov chains with 100 000 iterations, and taking account from the 10 001 iteration to the final estimation. The DIC was estimated with 100 000 additional iterations. The analyses were performed using OpenBugs 3.2.2 software (Lunn et al. 2000).
RESULTS
Population densities
Saguinus leucopus was found in 20 localities with densities of 6.0 to 149 ind/km2. A. hybridus was found in seven localities with densities of 0.62 to 58 ind/km2. C. versicolor was found in six localities with densities of 1.04 to 103 ind/km2, and A. seniculus was found in nine localities with densities of 0.21 to 75 ind/km2 (Table 3).
Local features, landscape, and class metrics
The localities were at elevations between 191 m and 1385 m. Corcovado, Campanas and Margaritas had the largest area of dense forest while Guarinocito and Fresno were the localities with largest areas of crops and pastures (Table 1; WCI). The localities with the more aggregated forests were Fresno, Río Manso, Guamo 1, Guamo 2, and Corcovado (One fragment with all pixels together), while Margaritas had the most disaggregated forest patches (Table 1; IJI). The localities with the most regular shaped forests were Guamo 1, Guamo 2 and Pajuil, while the forests of Pedrera and Primavera were the most elongated (Table 1; MSI).
Regarding plant structure and diversity, Guarinocito showed the highest mean DBH following by Corcovado and Río Manso. Playas had the smallest mean DBH. Río Manso, Amaní and Berlín had the tallest trees. The highest tree densities were found in Primavera and Berlín, and the lowest in Porce III (Table 4). The localities with the highest tree diversity were Río Manso and Berlin, while Porce III, Playas and Guarinocito showed the lowest tree diversity values.
Relationship of the explanatory variables in primate densities
Higher population densities of S. leucopus were recorded in localities with elongated, denser, and taller forests, in localities at lower elevations and with higher structural connectivity. These results also show that S. leucopus densities decrease when the density of A. hybridus increase (Table 5). Population densities of C. versicolor and A. hybridus tended to be higher when there was less connectivity. Regarding with plant structure and diversity, these results suggest that C. versicolor densities are positively affected only by tree DBH, while A. hybridus populations tend to increase in forests with higher DBH values, lower tree-height, and higher plant diversity (Table 5).
Population densities of A. seniculus were higher with denser forests, and lower in more elongated forests. Contrary to the pattern found for S. leucopus and C. versicolor, the population densities of A. seniculus increase when the density of A. hybridus increase and decrease when the lower DBH and tree density increase (Table 5).
DISCUSSION
These results show that the population densities of the four studied species are very variable among different localities, from local extinction, for example to A. hybridus in all localities in Caldas, to some of the highest population densities estimated for the taxon or congeneric species (Crockett et al. 1987, Chapman and Balcomb 1998, Gómez-Posada et al. 2007, 2009, Roncancio and Gómez-Posada 2009, Gómez-Posada et al. 2010, Williams-Guillén et al. 2013, De Luna and Link 2018). This indicates that different zones in the Magdalena River valley have been affected in different ways. Differences in the local occupation between species indicate that S. leucopus in this landscape could resists better the degradation, reduction and fragmentation of its habitat and synergic threats such as trade, given that, in spite of the large variation in its densities, it was found in all 20 localities. The other three species, in contrast, show less tolerance to changes in their habitats, and face a local extinction probability close to 60 % even though there are remnant forests, based in proportion of the occupied localities to each primate species.
Results obtained for C. versicolor and A. hybridus are consistent with the patterns found for Cebus capucinus (Linnaeus, 1758) and Ateles geoffroyi (Kuhl, 1820) in Nicaragua and with Cebus libidinosus (Spix, 1823) in Bolivia (Pyritz et al. 2010) where they have suffered local extinctions. Usually, A. seniculus and other howlers are recognized as being able to tolerate relatively higher levels of habitat fragmentation, as in many areas it is the only primate taxon to persist and in some it can be found at very high densities (Crockett et al. 1987, Chapman and Balcomb 1998, Gómez-Posada et al. 2007, 2009, Roncancio and Gómez-Posada 2009, Gómez-Posada et al. 2010, Williams-Guillén et al. 2013, De Luna and Link 2018). Red howlers have been locally extirpated in a number of localities in the Magdalena valley, and their densities are very variable in the localities where they are still present. A study of Alouatta guariba (Humboldt, 1812) carried out in southern Brazil showed a similar situation to that found for red howler populations in this study (Silva and Bicca-Marques 2013). The population depletion of howler monkeys could become driven by synergic threats and the degradation of the ecological conditions in the forest fragments.
Different features, landscape and class metrics are related differently the population densities of the four studied primate species, even though forest shape had the most significant effect on population density for all species. This situation was also found in Mexico, where the occupancy of all primates was negatively affected by shape irregularity (Arroyo-Rodríguez et al. 2008). In this study the effect of the shape was not the same for the four primates. While the population density of A. seniculus tended to decrease in elongated forests, the densities of the others three primates tended to increase (Lenz et al. 2014). High densities, however, do not always mean high abundance or better conditions for the population. High densities in forest remnants have been explained as a crowding effect (Chiarello and de Melo 2001, Peres 2001, Baranga 2004, Martins 2005, Rode et al. 2006). If the crowding is the first effect that suffer these primates, probably, the howler monkeys become affected later than the other species given that it could use some types of forest land cover that the other species do not use and in this way the population is diluted in a more extended area (Crockett 1998, Gómez-Posada and Londoño 2012). On the other hand, depending on the moment of the fragmentation process and the regeneration dynamics, early fragmentation, reversal isolation by regeneration in transitional crops, the smaller or more mobile species could or must use the elongated fragments to satisfy their ecological requirements, increasing their relative abundance, while the howler just need a small home range in the core areas (Matthews 2009, Spehar et al. 2010, Palma et al. 2011, Alba-Mejia et al. 2013). In this sense the fragment shape is an important factor that determines population size in fragmented landscapes, thus, habitat restoration efforts may be more effective if they focus on saving the forest fragments with less elongated shapes and keep a minimal width in the corridors (Ewers and Didham 2007).
Forest cover (WCI) was the second most important factor, when WCI was higher the population densities of all four species tend to be higher. Similarly, forest land cover was one of the main predictors of species persistence in isolated forests in a continental-scale analysis of neotropical primates (Benchimol and Peres 2013) and in two primate-assemblages in French Guiana and Ecuadorian Amazonia (Youlatos 2004). Forest land cover could be a determinant factor in the survival of these species (Sorensen and Fedigan 2000). Additionally, the general structure of the land cover landscape, associated with the matrix has been recognized as a factor that determines the presence of some primates in a fragment (da Silva et al. 2015).
Regarding plant structure, DBH was the main determining factor but, while the howler densities tend to be lower in forest with higher densities of large-trunked trees, densities of spider monkeys and capuchins tends to higher. In contrast, Alouatta palliata (Gray, 1849) has been shown to prefer fragments with a greater density of large trees (DBH > 60 cm) (Arroyo-Rodríguez et al. 2008). The equitability in the plant diversity was important for A. hybridus only, highlighting the high ecological requirements to this species and determine the necessity of focal intervention as specific plant enrichment to habitat restoration (Boyle and Smith 2010).
There is broad evidence that some types of agro-ecosystems in a fragmented landscape could be useful as temporary or additional habitat for different taxa (Estrada et al. 2006, Guzmán et al. 2016). However, I found that the primates in the Magdalena River valley are concentrated in zones with considerable areas of mature forest (WCI), including S. leucopus, even though tamarins and marmosets have preferences for regenerating forest (Snowdon and Soini 1988, Rylands 1996, Garber 1998, Defler 2010). This pattern is probably related to the ecological stability and resilience provided by mature forests. Thus, the biodiversity conservation at a landscape scale requires considerations regarding a minimum proportion of dense forest in the mosaic.
The densities of A. hybridus were significantly related to those of the other three primates, but, while howlers tend to be found in higher densities in forests with higher densities of Ateles, the other two species tend to be less abundant. The reason may be related to hunting pressure Howler monkeys come second to spider monkeys as preferred game by hunters (Freese et al. 1982, Rodríguez-Luna et al. 1996, Sorensen and Fedigan 2000, da Silva et al. 2005). So, it could be possible that the populations of A. seniculus will not start to decline by hunt until the density of A. hybridus is impacted by overhunting.
Ateles hybridus, by size, could be a better competitor than Cebus versicolor and Saguinus leucopus. The three species are overlapping in some niche dimensions (ecological requirements), so, by a compensation density process (MacArthur et al. 1972), the abundance of A. hybridus could regulate the density of the smaller primates. As a result, when the density of the bigger primates decreases, the density of the smaller primates increases (Peres and Dolman 2000, Hilário and Ferrari 2015), as has been observed in fragments where there were no spider monkeys, howlers or capuchins, and the density of the tamarins was of over 100 ind/km2. Thus, the high densities of S. leucopus reflect a very disrupted assemblage of primates, and I must not interpret them as a good indicator for the whole system.
I therefore suggest the use of A. hybridus as a landscape species to play the role as umbrella, keystone and maybe flagship surrogate to guide the priorities in the management of the humid and sub-humid forests of the Magdalena River valley (Sanderson et al. 2002). A. hybridus was the most affected species and by more factors (local vulnerability) (Rimbach et al. 2013) and it is a keystone species in this system (Christianou and Ebenman 2005). Thus, consolidating sustainable conservation scenarios for this primate, will satisfy the requirements for the other species, as a good representative of the ecological integrity of the whole system in a landscape approach.
The landscape approach to the conservation management usually demands the mosaic design that determines the minimal proportion and type of the core areas, their shape and aggregation and landscape composition. The landscape species approach based in the requirements or in the landscape structure effects allow us to plan and monitor the management actions in more cost-effective way.
However, landscape structure needs to be measured at the appropriate scale (Arroyo-Rodriguez and Fahrig 2014). The majority of the studies, including this, quantified landscape predictors within a single spatial scale, potentially missing significant primate-landscape responses (Galán-Acedo et al. 2019). Thus, studies on landscape effects on primates need to use a multiscale approach to ensure that landscape-species associations are correctly evaluated (Arroyo-Rodriguez and Fahrig 2014, Jackson and Fahrig 2015, Galán-Acedo et al. 2019). This study was done as patch-landscape analysis, and the scales were adjusted to each local landscape for a relative criteria based in the inferential area to primate density estimation. I think that standard measures of scales or direct scales are not applicable since the size of the inferential areas was highly variable, but a multiscale stratified analyses, with enough sampling sites, could be addressed in regions with similar land use patterns such as Caldas western region, Antioquia southeast and Serranía de San Lucas, that allow the use of standard scales, based on ecological criteria, to determine the optimal scale for each species to each region (Jackson and Fahrig 2012, 2015).