Introduction

Agricultural expansion is the leading cause of deforestation and biodiversity loss across the tropics (Curtis et al. 2018). The demand for commodities and the need to sustain an increasing human population means further agricultural expansion will be unavoidable (Laurance et al. 2014). Oil palm (Elaeis guineensis) is a prominent driver of tropical deforestation, with 21.5–23.4 million ha already under production globally (Descals et al. 2020). Biodiversity is heavily impoverished in intensively managed oil palm plantations, which typically support < 15% of the species found in native forests (Fitzherbert et al. 2008; Meijaard et al. 2018). However, there is growing evidence that by retaining forest remnants within plantations and managing them effectively, some wildlife species can survive in these human-modified landscapes (Deere et al. 2018, 2019; Mitchell et al. 2018). There is thus an increasing emphasis on including human-modified landscapes into the broader agenda of wildlife conservation in tropical countries (Ng et al. 2020).

Agricultural expansion contributes to the partitioning of remaining habitats into smaller, more isolated patches (Haddad et al. 2015). This fragmentation can limit individual movements over a landscape, thereby restricting the exchange of genetic information within and among populations, and inhibiting range shifts in response to environmental change (Årevall et al. 2018; Lino et al. 2019). Enhancing the connectivity value of human-modified landscapes is therefore a central theme in conservation by, for example, facilitating species movement between forest patches, which increases the population viability. Linear corridors and ‘stepping-stone’ patches of natural habitat are key ways by which such connectivity can be achieved (Keitt et al. 1997; Baum et al. 2004; Saura et al. 2014; Carroll et al. 2015).

Emerging environmental sustainability standards are beginning to recognise the importance of maintaining functional connectivity in agricultural landscapes. The Roundtable on Sustainable Palm Oil (RSPO), for example, is an international sustainable certification standard, which aims to alleviate both environment degradation and social impacts associated with oil palm production (https://rspo.org/). Recent uptake of zero-deforestation commitments by the RSPO seeks to align environmental sustainability and development goals by decoupling deforestation and agricultural expansion to ensure ecologically functional forest mosaics are retained during the development of new plantations (Meijaard et al. 2018; Deere et al. 2019). The High Conservation Value (HCV) concept has emerged as the principal methodology for palm oil producers to honour these commitments, and the High Carbon Stock (HCS) approach is integrated into the RSPO standard to assist in the quantification of HCV (http://highcarbonstock.org/leading-palm-oil-certification-system-adopts-no-deforestation-requirements/). HCS is a transparent land-use planning tool that aims to direct agricultural development towards degraded land of limited conservation value, while prioritising conservation set-asides based on ecological conditions, such as forest quality, fragment size and connectivity. While previous research on certification standards has emphasised the importance of forest fragments and riparian margins for biodiversity (Deere et al. 2018, 2019), the extent to which these habitat remnants promote functional connectivity at local scales remains poorly understood (Scriven et al. 2019).

Acquiring sufficient empirical data to assess connectivity is often prohibitively expensive, time consuming, and may raise welfare concerns if animals have to be habituated to the presence of observers or captured and restrained for individual identification or for tracking devices to be fitted (Gutema 2015). Therefore, ecological modelling remains an important tool for investigating landscape connectivity and the potential for animal movements across landscapes (Kool et al. 2013). Advances in spatially-explicit population modelling and the incorporation of stochastic environmental and biological processes may provide more realistic model outcomes than correlative approaches, particularly when applied to highly complex landscapes (Kearney and Porter 2009; Urban et al. 2016). Mechanistic dispersal models incorporate stochastic movement, whereby simulated individuals make probabilistic decisions governed by movement rules and, as a result, are likely to more closely align with reality (Palmer et al. 2011; Aben et al. 2016). If coupled with spatially explicit demographic models, the long-term effects of management options on population viability and species movement can be assessed concurrently (Bocedi et al. 2014; Cabral et al. 2017). This allows detailed investigations into management scenarios, which can provide a powerful tool to inform effective land-use planning and to direct research.

Here, we apply a spatially explicit individual-based model to test the conservation impact of several potential management scenarios for a highly degraded landscape in Sabah, Malaysian Borneo—a major palm oil producing region (Meijaard et al. 2018). Retaining forest fragments in the landscape is particularly important when new development takes place: i.e. through informed land-use planning such as the HCV and HCS approaches. Around 25% of the land in Sabah is planted with oil palm, with only a small proportion of the previous forest cover remaining in these plantations. However, every 25–30 years palms need to be removed and replanted, providing an opportunity to incorporate and restore additional forest fragments within existing farmland. Thus, the effectiveness of sustainable certification standards will be critical to ensuring positive conservation outcomes under future development. Our appraisal focuses on the Bornean orangutan (Pongo pygmaeus), a large-bodied flagship species characterised by a slow life-history and low population densities, attributes that make the species particularly vulnerable to the effects of habitat fragmentation (Marshall et al. 2009). On Borneo, orangutans have already lost substantial habitat, and research suggests further reduction of up to 57,000 km2, equal to a 20% decline (Struebig et al. 2015; Voigt et al. 2018) is possible by 2050 under a business-as-usual scenario. Recent surveys provide evidence that orangutans can persist at low densities in fragmented landscapes and oil palm estates where forest remnants are retained (Ancrenaz et al. 2015; Spehar and Rayadin 2017; Seaman et al. 2019). Initial estimates suggest as many as 10,000 orangutans, or between 10 and 15% of the remaining Bornean population, may already occur in oil palm plantations (Meijaard et al. 2017). As most of the remaining orangutan range is found outside of protected areas in Indonesia, this number will almost certainly increase with future forest conversion (Wich et al. 2012). It is therefore essential to understand how the processes of habitat modification affect orangutan population dynamics and connectivity across disturbed landscapes to provide some conservation options for the species. Hence, our study aimed to investigate the potential of prospective management scenarios to create agricultural landscapes that can also support orangutan populations and promote functional connectivity between large protected forests.

Methods

Study system

We modelled orangutan movement across the southern district of Tawau, an area of ca. 208,000 ha (51.8 km × 40.3 km) encompassing the Stability of Altered Forest Ecosystems project (SAFE; www.safeproject.net) and surrounding oil palm estates. The landscape is a heterogeneous mosaic of degraded forest remnants embedded within oil palm plantations, and to a lesser extent, non-native softwoods (Fig. 1, panel I). Palms range in age from newly established plantations to mature stands of > 15 years (Mitchell et al. 2018). The landscape also harbours several large areas of intensively logged lowland and hill dipterocarp forest, including protected forest in Mt. Louisa Forest Reserve (contiguous with ca. 1 million ha of forest, of which 18,160 ha is situated in the study system), Tawau Hills (28,000 ha; of which 3890 ha is situated in the study system), and four smaller Virgin Jungle Reserves (525–2000 ha). The study area contains multiple settlements and an extensive network of permanent roads since much of the landscape is actively managed as plantation. Orangutans have been documented in all forest types across the landscape (Ancrenaz et al. 2004a; Bernard et al. 2016; Seaman et al. 2019).

Fig. 1
figure 1

Model framework and study system. I The study system is a highly fragmented landscape in the north of Borneo (location in inset). II Modelling framework: a model inputs used by RangeShifter 2.0. b The stage-structured demographic model of nine stages ni, where σ signifies survival probability of each stage, y denotes the probability of an individual developing to the next stage, φ signifies individual fecundity and d indicates where density dependence in development is applied (Modified from Neubert and Caswell 2000). c The dispersal model, both emigration and settlement probabilities (P) are dependent on the density relative to the carrying capacity of the habitat patch (K). The transfer process uses a stochastic movement simulation implemented in RangeShifter 2.0. d Model outputs are both temporally and spatially explicit

Landscape scenarios

We modelled four land-use scenarios that could be reasonably expected given current environmental policy and conservation approaches in Sabah and typical to other human-modified tropical landscapes: ‘Land Sparing’, ‘Uncertified’, ‘Certified’, and ‘Conservation Enhanced’. ‘Land Sparing’ assumes the conversion to agriculture of all land except that which is strictly protected (Class I Forest Reserve and Class IV Virgin Jungle Reserve), thus representing a worst-case scenario (Fig. 2). We simulated an ‘Uncertified’ Landscape, by applying the minimum environmental policies currently in place for Sabah, whereby all areas were converted except the protected areas, riparian buffers 20 m each side of permanent rivers and areas above 25° slope. A ‘Certified’ Landscape was created by following the High Carbon Stock (HCS) approach decision tree to prioritise land for conservation (Rosoman et al. 2017), through which forests were classified into strata defined by thresholds of carbon density, and forest patches prioritised based on their core area after applying a 100 m internal buffer into High (> 100 ha), Medium (10–100 ha), of Low priority (< 10 ha). Forest patches were further prioritised for conservation or development by patch proximity, forest quality and risk from anthropogenic activities. Lastly, we constructed a ‘Conservation Enhanced’ landscape, using recommendations from published literature. Here, we increased the riparian buffer width to 45 m each side of the river, as this is recently recommended to improve biodiversity outcomes (Gray et al. 2014; Mitchell et al. 2018), and increased the core area of the HCS medium priority patches to 100 ha but removed low priority patches in the ‘give and take’ process defined by the HCS protocol.

Fig. 2
figure 2

Landscape scenarios based on plausible management options in southern Tawau district, Sabah

Modelling framework

We applied a modified version of RangeShifter 2.0 (Bocedi et al. 2020), a freely available individual-based modelling platform, to model orangutan population persistence and connectivity across our study landscape (Fig. 1, panel II) and landscape scenarios. RangeShifter simultaneously models population dynamics and landscape connectivity, by integrating spatially explicit demographic and dispersal models (Bocedi et al. 2014).

Patch allocation and orangutan demography

The landscape was defined on a gridded system at a resolution of 30 m × 30 m. We used a patch-based approach to model population demography; whereby adjacent cells of suitable habitat were aggregated into discrete patches. The equilibrium density of each patch was based on habitat type, using existing orangutan density estimates from the same site (Table 1) (Seaman et al. 2019). In large forest areas we reduced density estimates by half to represent the female density. However, for small remnant forest patches that are likely to be occupied by resident females (Ancrenaz et al. 2021) we maintained density estimates to reflect this. Female orangutan home ranges are difficult to determine, and on Borneo estimates vary from 40 to 600 ha (Singleton et al. 2009). Although density estimates varied widely across the landscape, the mean density in remnant forest in oil palm was 0.82 individuals/km2, equating to a minimum of 122 ha to support a single orangutan. We therefore considered 122 ha of suitable habitat to be the minimum patch size (which would therefore sustain a single adult female) (Seaman et al. 2019). We considered suitable habitat to be any area of natural forest defined following HCS protocols as cells with 35 t C ha−1, produced using LiDAR data from the Carnegie Airborne Observatory (Asner et al. 2018). Mean carbon values extracted from transects in remnant forest in oil palm estates where orangutan nests have been observed range from 3.04 to 106.03 t C ha−1, and therefore, forest delineated by this carbon threshold is known to be capable of supporting orangutans. Habitat patches large enough to support one or more orangutan (i.e. > 122 ha) were numbered with a unique numeric identifier (detailed in Supplementary Information SI.2).

Table 1 Parameter values included within the model

We developed an overlapping stage-structured demographic model limited to females, and comprising nine age-related stages (Table 1). There is a slight male-bias sex ratio in orangutans of 55%, so we set survival probability of 0.45 at the neonate stage. Subsequent stage survival probabilities were derived from the 2019 Bornean Orangutan Population and Habitat Viability Analysis (PHVA; Utami-Atmoko et al. 2017). We added density dependence in development between the young adult and adult stages, so that females becoming reproductive would be delayed if their resident patch was at, or close to, population equilibrium density. We modelled a yearly reproductive season. However, as orangutans have a long interbirth interval ranging from 6 to 9 years (van Noordwijk et al. 2018), we set the annual fecundity at 0.167 (equating to a 6-year inter-birth period). Density dependence was also incorporated in fecundity (with a bespoke function added to match Utami-Atmoko et al. 2017), with the interbirth period increasing as patches became close to equilibrium density (Supplementary Information SI.1).

Dispersal

Females are highly philopatric (van Noordwijk et al. 2012); however, little is known about how female orangutans disperse in oil palm landscapes. Since females have been observed in forest patches in human-modified landscapes many years after conversion (Ancrenaz et al. 2015; Spehar and Rayadin 2017; Seaman et al. 2019), it is likely these animals continue to display a high level of home range fidelity even within these highly disturbed habitats. To account for this, we used a strongly density-dependent emigration probability, meaning orangutans remained within their natal patch until the patch reached its equilibrium density, at which point there was a conservative 0.2 maximum emigration probability (Fig. SI.2a). Settlement probability was also density-dependent. As female orangutans tend to stay close to their mothers’ home range (Goossens et al. 2006), we set a probability of 1 of an orangutan immediately settling in a new patch, unless that patch was near equilibrium, at which point there was a shallow decline in settlement probability (Fig. SI.2b). The emigration and settlement values we selected resulted in a pattern where females only emigrate if absolutely necessary (i.e. when there is no possibility of reproducing in the natal patch) and would likely settle at the first available opportunity, which is supported by our current understanding of orangutan behaviour (Ashbury et al. 2020).

We modelled dispersal movement with the stochastic movement simulator (SMS), which simulates stepwise nearest-neighbour movements informed by a cost surface (Palmer et al. 2011). Additionally, individual movements depend on three parameters: perceptual range (the distance at which the individual can evaluate its surroundings), directional persistence (DP, the tendency of an individual to move in a straight line) and memory size (the number of previous movement steps used to calculate the directional persistence).

As orangutans are arboreal, they are likely to have a large perceptual range. We therefore set the parameter to 25 cells (750 m). Simulated trajectories resulting from a range of values were visualised and plausible combinations were selected by expert judgement (the approach is consistent with that used for determining the cost values). We found plausible paths for intermediate values of directional persistence (between 2.0 and 3.0 when memory was set to 10: detailed in Supplementary Information SI.4). We thus used directional persistence = 2.5 and memory = 10 as our baseline values but also undertook a sensitivity analysis to investigate the impact of changing these parameters.

The study area is characterised by little hunting or conflict killings of orangutans and has a large number of small patches that orangutans can potentially utilize for resources. We therefore assumed mortality rates in the landscape to be low (per step mortality 0.001) and set the maximum number of steps per year to 3000 (equating to a path length of approx. 108 km), which seems reasonable based on daily path length from wild orangutans multiplied over the same period (Singleton et al. 2009) and visual inspection of dispersal trajectories. We also set the total maximum number of steps an individual may make during dispersal to 12,000 (approx. 435 km), and any individual exceeding that limit would die. Since the area is part of a much larger multi-use landscape and dispersing animals are not constrained to the area, we applied a reflective study boundary as orangutans may also enter the landscape from outside the study system.

Cost surface

Orangutans are primarily arboreal but will also travel terrestrially where necessary (Ancrenaz et al. 2014). However, this form of terrestrial locomotion incorporates additional energetic costs (Thorpe et al. 2007). To capture this cost, we produced a composite cost surface layer informed by both expert-derived landscape resistance and a time cost model (Frakes et al. 2015) (detailed in Supplementary Information SI.2). The resulting cost surface is time travel distance in seconds, weighted by the resistance of the landcover type.

Model initialisation and metrics

We initialised the model on each landscape scenario so that all suitable patches were occupied at female carrying capacity (based on habitat type), with the population distributed with a right-skew towards older ages classes (20% of the population among Adolescent or below and 80% among Young adult and above Ancernaz et al. 2004b). For each scenario, we modelled population dynamics over 100 years for 100 iterations. To evaluate the effects of each scenario, we derived several metrics, with respect to population dynamics and movement, recorded after the population had reached equilibrium (see Fig. SI.3). Mean Population size was calculated for the whole landscape and agricultural areas only (i.e. excluding all protected areas and VJRs) at 100 years. Dispersal Distance was derived from the centre of the cell from which the individual initially dispersed to the centre of the cell at which it settled in a non-natal patch. We determined Dispersal Success as the proportion of dispersing individuals that either settled in a non-natal patch, returned and settled in their natal patch, or died during transfer. Relative Dispersal Success was derived by comparing the annual number of individuals that either returned and settled in their natal patch, settled in a non-natal patch or died during transfer, compared to the ‘Uncertified’ scenario. We also created network maps, by plotting links between patches where individuals had successfully dispersed and settled.

Sensitivity analysis

As our movement parameters were largely based on anecdotal evidence and expert opinion, we undertook a sensitivity analysis to evaluate the robustness of our model to permutations in parameter values (detailed in Supplementary Information SI.5). To assess the model sensitivity to permutation in parameter values, we compared Mean Population size at 100 years between landscape scenarios and percentage change in Dispersal Success relative to the baseline scenarios, under the range of parameter values.

Results

Population size in each landscape configuration

We modelled orangutan population dynamics over 100 years in four different landscape scenarios. At a mean of 645 individuals in year 100 [± standard error (SE) 2.27 individuals], the ‘Conservation Enhanced’ configuration supported the largest simulated orangutan population across the whole landscape. The ‘Certified’ Landscape supported a comparable Mean Population of 612 individuals (± 2.12). The smallest population sizes were predicted for the ‘Uncertified’ Landscape (450 ± 1.70 individuals) and the ‘Land Sparing’ scenario (445 ± 1.62 individuals). When we excluded protected areas and assessed the population in the agricultural landscape only, the largest population was again predicted for the ‘Conservation Enhanced’ scenario (214 ± 1.23 individuals), with only a slight reduction estimated for the ‘Certified’ landscape (181 ± 1.18 individuals). The ‘Uncertified’ landscape had a final estimated Mean Population of 20 individuals (± 0.33), and with no suitable habitat within the production landscape the ‘Land Sparing’ scenario had an estimated Mean Population of < 1 individual (± 0.04; Fig. 3).

Fig. 3
figure 3

Raincloud plot showing Population size at 100 years over 100 iterations for each of the four management scenarios. Right side of each plot shows the probability density, with the left showing the raw data and boxes denoting the median, first and third quartiles, and whiskers indicating 95% confidence intervals. a Population for the whole landscape including the protected areas and remnant forest patches in the production landscape, and b population size for patches in the production landscape only

Dispersal Distance and Success

In terms of absolute numbers, when considering all patches in the landscape the mean number of annual successful dispersers settling into non-natal patches averaged 2.83 (± 0.02, Fig. 4) in the ‘Certified’ scenario and 3.16 (± 0.03) in the ‘Conservation Enhanced’ scenario, which was an increase relative to the ‘Uncertified’ landscape (1.08 and 1.33 individuals respectively). If limited to only individuals dispersing from one of the eight protected areas, annually an average of 2.03 (± 0.02) individuals settled in non-natal patches in the ‘Certified’ scenario and 2.16 (± 0.02) in the ‘Conservation Enhanced’ scenario. The mean number of dispersing orangutans that went back and settled in natal patches annually was comparable across the ‘Certified’ (2.07 ± 0.02), ‘Conservation Enhanced’ (2.00 ± 0.02) and ‘Uncertified’ (1.93 ± 0.02) landscape scenarios when considering all patches. When constricted to the protected areas only, fewer individuals settled in their natal patches in the ‘Certified’ (1.68 ± 0.02) and ‘Conservation Enhanced’ (1.54 ± 0.01) scenarios, compared to the ‘Uncertified’ (1.89 ± 0.02) landscape. The ‘Land Sparing’ scenario which only included the eight protected areas, experienced the lowest annual number of individuals settling into non-natal patches (1.05 ± 0.01) and the highest number of individuals returning and settling in their natal patches annually (2.28 ± 0.02).

Fig. 4
figure 4

Raincloud plot showing of mean annual Dispersal Success. We define successful dispersers as individuals that either settled in a non-natal patch or returned emigrated but then returned to their natal patch. Unsuccessful dispersers were individuals that died during dispersal

Dispersal mortality rate was 46% in the ‘Land Sparing’ scenario, with a similar rate recorded in the ‘Uncertified’ landscape (44%, Fig. 5). Mortality rate during dispersal was lower in the ‘Certified’ and ‘Conservation Enhanced’ landscapes (36% and 35% respectively). In the ‘Land Sparing’ scenario, only 13% of dispersing individuals settled in non-natal patches. This increased slightly in the ‘Uncertified’ landscape to 18%, and again increased further to 30% in the ‘Certified’ and 32% in the ‘Conservation Enhanced’ scenarios. In the ‘Land Sparing’ and ‘Uncertified’ scenarios a slightly larger proportion of dispersing individuals returned and settled in their natal patches (29% and 25% respectively), compared to the ‘Certified’ (22%) and ‘Uncertified’ (21%) landscapes (see Fig. 5). When only considering successful dispersers (those that either returned to their natal patch or settled in a non-natal patch), the difference between scenarios was greater, with most individuals in the ‘Land Sparing’ and ‘Uncertified’ scenarios returning to their natal patches (68% and 57% respectively). In the ‘Certified’ and ‘Conservation Enhanced’ scenarios over half the successful dispersers settled in non-natal patches (58% and 61% respectively) compared to 41% in the ‘Uncertified’ landscape, and just 32% in the ‘Land Sparing’ scenario.

Fig. 5
figure 5

Dispersal network maps for the four management scenarios. Lines on maps indicate where successful dispersal events between suitable habitat patches occurred, scaled in size and colour for mean number of individuals per year over a 40-year period. Each point denotes the Mean Population size for individual patches at the end of each simulation (n = 100). The probability density plots show the distribution of Dispersal Distances for each management scenario, with dotted red lines denoting the median distance travelled during dispersal. Dispersing individuals are characterized into three dispersal statuses, settled in their natal patch, settled in a non-natal patch, or died during transfer. The proportion of dispersing individuals within each dispersal status is represented by stacked bar charts. Red patches indicate patches which produced no successful emigrants

Across all four scenarios, Dispersal Distance ranged between 0.03 and 58 km, with the distribution heavily skewed towards shorter distances (Fig. 5). Under the ‘Certified’ and ‘Conservation Enhanced’ scenarios median Dispersal Distances were 5.96 km and 5.82 km respectively. The largest median Dispersal Distance was estimated for the ‘Uncertified’ scenario (7.35 km), with the shortest median distance modelled in the ‘Land Sparing’ landscape (4.15 km). Across scenarios, there were successful emigration events from all habitat patches, with the exception of 3 patches embedded in the plantation matrix in the Uncertified landscape, indicating possible population sinks (Fig. 5).

The higher number of suitable patches in the ‘Certified’ and ‘Conservation Enhanced’ landscapes (Fig. 2) provided many more potential connections, with multiple linkages between patches (mean number of connections across all model iterations 37.6 ± 0.42 and 49.2 ± 0.52 respectively), compared to the ‘Uncertified’ and ‘Land Sparing’ scenarios (11.3 ± 0.19 and 1.62 ± 0.65 respectively). Across all scenarios the largest number of potential connections summed over all 100 iterations were for the ‘Conservation Enhanced’ scenario with 348 and the ‘Certified’ Landscape with a total of 248 connections (Fig. 5). The development focused scenarios again had substantially fewer connections: 62 for the ‘Uncertified’ landscape and just 13 in the ‘Land Sparing’ scenario (Fig. 5). The mean number of connections per patch also differed between scenarios, with 1.77 (± 0.02) for the ‘Conservation Enhanced’ scenario, 2.21 (± 0.02) in the ‘Certified’ Landscape, 1.15 (± 0.02) in the ‘Uncertified’ landscape and with the fewest 0.44 (± 0.02) in the ‘Land Sparing’ scenario.

Sensitivity analysis

With the exception of the ‘Land Sparing’ and ‘Uncertified’ scenarios, the model outcomes did not change the ranking of the scenarios for most parameter combinations. Variations (± 99 individuals) in population size did change the ranking of ‘Land Sparing’ and ‘Uncertified’ scenarios; however, their baseline population size did not differ substantially (‘Land Sparing’ 445 ind. and ‘Uncertified’ 450 ind.: see Fig. 6). The model was sensitive to Per Step Mortality (the probability of an individual dying at each step taken during dispersal) with final population estimates changing the ranking of several of the scenarios. When Per Step Mortality was reduced to 0.0001, the Mean Population estimate for the ‘Uncertified’ landscape increased by 99 individuals (18%), increasing its ranking above the ‘Certified’ landscape (Fig. 6). Similarly, reducing Per Step Mortality increased the population estimate of the ‘Certified’ landscape by 100 individuals (14%) elevating it above the baseline estimate for the ‘Conservation Enhanced’ scenario by 67 individuals (Fig. 5). When Per Step Mortality was increased to 0.01, the final population estimates for the ‘Conservation Enhanced’ landscape reduced by 12%, which was 37 fewer individuals than the baseline scenario for the ‘Certified’ scenario (Supplementary Information Table SI.5.1).

Fig. 6
figure 6

Outcomes of sensitivity analyses showing the Mean Population size at 100 years across 100 iterations of each model parametrisation. Model parameters on the x axis, including pairs of extreme higher (+) and lower (−) parameter values. Solid black lines show the Mean Population size for each scenario, with the shaded area indicating the standard error

Dispersal Success also appeared sensitive to a reduction in per step mortality, with > 100% increase in individuals settling in non-natal patches in all scenarios with the exception of the ‘Land sparing’ scenario which saw a large increase (186%) in individuals returning to their natal patches. A similar pattern was seen when per step mortality was increased but to the lesser degree of change < 50% change in the number of individuals (Supplementary Information Table SI.5.2) (Fig. 7).

Fig. 7
figure 7

Sensitivity analysis results, showing percentage change Dispersal Success mean number of individuals per year, compared to baseline models. Model parameters on the x axis, including pairs of extreme higher (+) and lower (−) parameter values

Discussion

Our models pertaining to animal movement and population ecology in a highly fragmented landscape demonstrate that management options to maximise forest cover (‘i.e. the Certified’ and ‘Conservation Enhanced’ scenarios) can sustain substantial orangutan populations in areas under agricultural production for a period of at least 40 years. Oil palm plantations are productive for 25–30 years before replanting (Ng et al. 2013), and therefore, these landscapes are likely able to support stable populations of orangutans over the lifespan of a plantation. This is in line with several field observations that orangutans have survived in similar human-modified landscapes for long periods (> 20 years: Ancrenaz et al. 2021). To date, orangutan research has focused almost exclusively on intact landscapes, and no long-term empirical data are available from human-modified landscapes to investigate population dynamics (Voigt et al. 2018). Our models provide longitudinal insights into the increasingly recognised potential of human-modified landscapes to support orangutan populations. These results require validation from field observation and highlight the need for increased research focus in these highly modified landscapes.

A striking difference between the outcomes of the landscape scenarios was the proportion of dispersing individuals that either settled back in their original natal patch or transferred to a non-natal patch. Our simulations demonstrated that in the ‘Land Sparing’ scenario, where there is limited structural connectivity (e.g. the absence of riparian buffers and other remnant forest patches), individuals are largely confined to protected areas, with the majority of successfully dispersing females settling within their natal patches, and hence unable to disperse elsewhere. In the ‘Land Sparing’ scenario there was also a limited number of connections among patches (a mean of 1.9 connections over 100 iterations) and dispersing individuals moved over the shortest distances before settling (median 4.1 km). This short-distance dispersal is supported by observations of females from areas of continuous forest, that are generally philopatric and will move and settle close to their mothers’ home range (Goossens et al. 2006; Ashbury et al. 2020). The inclusion of Malaysia’s current minimum ecological requirements in the ‘Uncertified’ scenario (retention of 20 m riparian buffers on permanent rivers and forest on slopes above 25°) increased the connectedness of the landscape, with an average of 12.7 connections and 30% of successfully dispersing individuals settling in non-natal patches. However, this rise in animal movement across the landscape did not greatly increase the total population size (a population increase of 0.2% compared to the ‘Land Sparing’ scenario) and produced the largest median Dispersal Distance at 7.4 km, compared to just over 6 km in the ‘Certified’ and ‘Conservation Enhanced’ scenarios. The retention of additional natural forest through conservation set-asides in the ‘Certified’ and ‘Conservation Enhanced’ scenarios resulted in higher proportions of individuals settling in non-natal patches (58% and 63% of all successful dispersers respectively) and a substantial number of connections amongst patches (an average of 40.3 connections in ‘Certified’ and 52.4 in ‘Conservation Enhanced’). The retention of natural forest also led to a considerable increase in the final population size (by 20% in the ‘Certified’ scenario, and 24% in the ‘Conservation Enhanced’ compared to the ‘Land Sparing’ scenario).

In fragmented landscapes, facilitating orangutan movement between isolated populations will be a key conservation strategy to ensure the long-term genetic health of populations (Templeton et al. 1990; Bruford et al. 2010) and to facilitate range shifts in response to climate change (Struebig et al. 2015; McGuire et al. 2016). Our simulations suggest orangutans are more likely to move short distances across agricultural matrixes when connected by smaller fragments or riparian remnants, rather than undertaking long excursions into plantations (see Fig. SI.4). Similar behaviour has been observed with wild female orangutans in oil palm landscapes, where individuals have been reported to move between areas of natural forest to cross plantations (Ancrenaz et al. 2015). Thus, small patches can function as stepping-stones and are likely to be of high importance in facilitating orangutan movement across human-modified landscapes, whilst not necessarily supporting high levels of biodiversity on their own (Deere et al. 2019). A phenomenon we captured in our simulations were female orangutans occupying and reproducing in remnant forest patches in oil palm. This is consistent with field observations, where practitioners have identified female orangutans residing in remnant forest in oil palm dominated landscapes for a sufficient amount of time for an offspring to reach maturity and leave to establish its own home range (Ancrenaz et al. 2021).

These individuals born in remnant forest patches will have limited opportunities to establish their own home range and this will be contingent on the amount of remaining natural forest. However, the highly philopatric nature of female orangutans means movements over long distances for this sex are likely to be extremely rare. Our simulations suggest that orangutans in heavily modified landscapes may have to change their behaviour, moving over larger distances to find suitable areas to settle, particularly when there is minimal remaining natural forest cover. In the ‘Uncertified’ scenario, three of the remnant forest patches were occupied by dispersing individuals but failed to produce successful emigrants, representing potential population sinks. Thus, we highlight the need for not only providing linear forest corridors (such as riparian buffers) and/or stepping stones to increase the permeability of the landscape, but also to increase the potential for females to establish their own home ranges, reducing the need for long distance dispersal and maximising the population size.

Ensuring sufficient habitat remains or is restored in agricultural landscapes will be vital to promote functional connectivity, but this is only half the battle. Connectivity may be overestimated if negative human–wildlife interactions (such as crop foraging or hunting) are not considered (Day et al. 2020; Ghoddousi et al. 2020; Bleyhl et al. 2021). Our sensitivity analysis suggests changes in mortality risk during dispersal are likely to have a large effect on the ability of the agricultural landscape to support orangutan and allow movement between patches. As part of the certification process itself, the High Carbon Stock methodology includes a risk assessment to prioritise patches for conservation based on proximity to anthropogenic features (such as roads and settlements), aiming to limit the impact of human–wildlife interactions. Despite this, human accessibility to forest patches is high in oil palm estates and hunting can be common (Azhar et al. 2013; Deere et al. 2019). Here, we assumed all HCS forest areas to be suitable for orangutans, however, other factor such as hunting may limit the ability of these areas to support orangutans. In Kalimantan, hunting of orangutans for meat and persecutory killings following negative orangutan–human interactions is widespread (Meijaard et al. 2011). Although hunting of orangutans in Sabah appears to be relatively low (Ancrenaz et al. 2007), further research is needed to fully establish the extent to which human–wildlife interactions will inhibit the ability of oil palm certification to facilitate movement for orangutans.

A conservation strategy routinely employed in anthropogenic landscapes is the translocation of orangutans from forest fragments to areas of continuous forest or rescue centres (Sherman et al. 2020a). However, most orangutans translocated from agricultural landscapes appear to be in good health (Sherman et al. 2020b), indicating that these individuals have been surviving in these landscapes. In addition, our models demonstrate the potential importance of these individuals in maintaining movement and connectivity across the landscape (Ancrenaz et al. 2021). Although orangutan translocation can be justified when the life of an individual is threatened, such as during fire events or direct conflict with people, the blanket deployment of translocation may in fact be detrimental to the species (Sherman et al. 2020b). Indeed, we showed that maintaining orangutan numbers within remnant forest in plantations will increase the overall population size and increase movement across the landscape, potentially providing vital transfer of genetic information between isolated individuals or groups of individuals.

There is mounting evidence that orangutans can survive and reproduce within plantation-dominated landscapes (Ancrenaz et al. 2015; Spehar and Rayadin 2017; Seaman et al. 2019) and our models show these animals are likely to play an important role in maintaining movement between otherwise isolated populations. However, there remains little known regarding orangutan movement behaviour in human-modified landscapes. Therefore, we used our current understanding of orangutan dispersal and expert judgment to parametrise the model and create an approximation to the same pattern as observed in the field. Our use of an individual based model allowed the incorporation of individual variability in both movement and demography, providing a more realistic process than other modelling approaches. There is unavoidably a level of subjectivity to this approach, and to address this we performed an extensive sensitivity analysis. Encouragingly, the model projections, especially in how the alternative plausible scenarios ranked, seem robust to permutations in almost all emigration, movement, and settlement parameter values, giving us confidence in the conclusions we draw here. The orangutan’s slow life history and low reproductive rate makes the species highly vulnerable to even small rises in mortality rates above natural levels (~ 1%: Leighton et al. 1995). This sensitivity highlights the critical need to ensure mortality in human-modified landscapes is kept to a minimum though reducing hunting and conflict killings, raise the level of acceptance of people sharing the same habitat for peaceful co-existence via targeted awareness campaigns and capacity building, as well as increasing conservation focus towards these areas.

Conclusion

With increasing environmental degradation, protected areas alone will be insufficient to secure conservation goals for much of the earth’s biodiversity (Dinerstein et al. 2020). Maintaining wildlife populations and ensuring connectivity in agricultural landscapes is therefore now essential, particularly to sustain viable populations of large-bodied and wide-roaming terrestrial mammals (Carroll et al. 2015). The extent to which this can be achieved will depend on the overall design of the landscape and how land-use practices meet the ecological needs of species, as well as the acceptance of people to coexist in proximity with wildlife. We need a paradigm shift about how we conserve wide-ranging species such as orangutans—to embrace landscape-level management in human-modified habitats, as well as staunch protection in intact forest areas (Kremen and Merelender 2018). There is currently a paucity of research from these landscapes and additional research will be vital to better inform land-use policy and focus conservation efforts. As further agricultural expansion is unavoidable, our modelling suggests that maximising natural forest cover in farmland landscapes through conservation set-asides, will lead to improved long-term conservation outcomes for critically endangered species such as orangutan.