Passive acoustics and sound recognition provide new insights on status and resilience of an iconic endangered marsupial (koala Phascolarctos cinereus) to timber harvesting

Retention forestry aims to mitigate impacts of native forestry on biodiversity, but data are limited on its effectiveness for threatened species. We used acoustics to investigate the resilience of a folivorous marsupial, the koala Phascolarctos cinereus, to timber harvesting where a key mitigation practice is landscape exclusion of harvesting. We deployed acoustic recorders at 171 sites to record male bellows (~14,640 hours) for use in occupancy modelling and for comparisons of bellow rate (bellows night-1). Surveys targeted modelled medium-high quality habitat, with sites stratified by time since logging and logging intensity, including old growth as a reference. After scanning recordings with software to identify koala bellows, we found a high probability of detection (~0.45 per night), but this varied with minimum temperature and recorder type. Naïve occupancy was ~ 64% across a broad range of forests, which was at least five times more than expected based on previous surveys using alternative methods. After accounting for imperfect detection, probability of occupancy was influenced by elevation (-ve), cover of important browse trees (+ve), landscape NDVI (+ve) and extent of recent wildfire (-ve, but minor effect). Elevation was the most influential variable, though the relationship was non-linear and low occupancy was most common at tableland elevations (> 1000 m). Neither occupancy nor bellow rate were influenced by timber harvesting intensity, time since harvesting or local landscape extent of harvesting or old growth. Extrapolation of occupancy across modelled habitat indicates that the hinterland forests of north-east NSW support a widespread, though likely low density koala population that is considerably larger than previously estimated. Retention forestry has a significant role to play in mitigating harvesting impacts on biodiversity, including for forest specialists, but localised studies are needed to optimise prescriptions for koalas.

Introduction with low elevation [26]. Moreover, koala scats are correlated positively with the number of selective harvesting events, indicating that koala populations are resilient to historical, lowintensity, harvesting [31,32]. But scats are also associated with structurally complex, unevenaged forests with some mature and old-growth elements, a large basal area and mixed species associations dominated by preferred browse species [32]. These results led [32] to suggest that high intensity timber harvesting that creates extensive gaps, especially those lacking mitigation measures, are incompatible with koala conservation. In the last 20 years, specific koala prescriptions have been implemented to protect koalas in forests harvested for timber, though they have relied upon detection of koala scat accumulations (high-use areas) during pre-harvest surveys [33].
Our study assessed koala occupancy and bellow rate across the extensive forested area of north-east NSW, Australia. We specifically targeted areas with different timber harvesting intensities and times since harvesting, as well as old growth forest. Harvesting ranged from recent modern practices where varied habitat prescriptions aim to mitigate environmental impact [7], to forests regenerating after historical harvesting that pre-dated standard environmental prescriptions. We also surveyed koala high-use areas that previous pre-harvest surveys identified based on scat accumulations, which then triggered a small 'high-use' exclusion area specifically for koalas. We surveyed koalas passively with acoustic recorders as these are effective for recording the male mating bellow during the spring breeding season [19,29]. Occupancy and bellow rate were compared across harvesting treatments after first controlling for imperfect detection [34]. A second aim of our study was to provide an updated assessment of the status of koalas in the hinterland north-east forests of NSW based on occupancy estimates. Occupancy is an ideal population parameter for studies encompassing regional scales, especially as koalas typically occur at low densities [32,34].
We predicted that if timber harvesting impacted koala populations, occupancy and calling rate would be greater in old growth forests and selectively harvested forests than heavily harvested forests, but that this impact would be moderated in areas with recently implemented prescriptions designed to protect koalas. We also predicted that occupancy and calling rate would increase with time since harvest.

Study area
The study area spanned native forests of the hinterland, ranges and tablelands of north-east NSW bounded by the Hunter River in the south, the Queensland border in the north and Armidale in the west (~8.5 million ha) (Fig 1). The coastal strip was generally excluded as koalas in this area are exposed to multiple threats associated with urbanisation. Both State forests and National Parks were sampled, while private land (including private forestry) was excluded. The study region supports a rich diversity and mosaic of forest types, from tall moist eucalypt forest to dry sclerophyll forests and woodlands, including a range of important koala browse tree species [17,35,36]. A scientific license for our research was provided by the NSW national Parks and Wildlife Service (SL 102086). No animal ethics approval was required as only passive acoustic recording was undertaken.

Survey design
Surveys targeted modelled medium-high quality koala habitat (MaxEnt model >0. 45), which covers~1.66 million ha of forest in north east NSW [19]. Thus forests of low predicted habitat quality and non-habitat (e.g. rainforest, heath) were excluded. We stratified our survey by time since harvesting and harvest intensity and sampled eight treatments (Table 1). Between 20-25 acoustic recorders (SongMeter SM2 or SM4, Wildlife Acoustics, Maynard, USA) were set each sampling trip, with one recorder located near the centre of an allocated treatment, though within a local landscape comprising a mosaic of forest ages and timber harvest exclusion zones. We sampled two levels of harvest intensity: light-moderate selective (<80 m 3 timber removed per compartment (~250 ha)) and heavy harvesting (>80 m 3 timber removed). Harvest intensity was designated based on recent Geographic Information System (GIS, ArcMap 10.4.1, ESRI) layers (GISO.EventPoly_State; FCNSW unpubl. data) or, for older forests (prior to 2001), using a Management History layer containing volumes of timber removed from compartments (FCNSW unpubl. data). Heavy harvesting, also referred to as regeneration harvesting, targeted dry hardwood forests, especially E. pilularis where site disturbance aids regeneration [37]. Within each harvest treatment, we sampled three different times since harvesting; recent (2-10 years since harvesting), medium (11-25 years since harvesting) and old (> 25 years since harvesting). Older harvested forests had few environmental protection measures at the time of harvest as environmentally sustainable practices had yet to develop [7]. In comparison, recent heavy harvesting incorporated exclusions for environmental protection, including for riparian zones, rainforest, old growth forest and specific protective measures for threatened species [7]. However, this treatment lacked any specific prescriptions for koalas because they were not previously detected by scat surveys. We also sampled koala high-use areas, which were areas where scat accumulations were detected by forestry surveys, and which were subsequently excluded from harvesting (in addition to other landscape exclusions). High-use areas are designated during pre-harvest surveys by either sighting a koala, identifying trees with > 20 scats or 30% of searched trees containing a scat [32]. Such areas were typically small patches (mean = 4 ha) with surrounding harvests occurring mostly < 10 years prior to acoustic sampling. Our final treatment was old growth forest, which was sampled for comparison with harvested sites. Old growth is mapped as a GIS layer (High Conservation Value Old Growth Forest) and although some of these sites had been lightly harvested historically, all sites represented long undisturbed forests. Within these treatments we sampled a range of topographic locations (gullies to ridge tops), forest types and elevations (10 to 1327 m above sea level (ASL)). Site locations were initially selected within GIS to sample the full range of forest treatments in each sampling session, and modified in the field if GIS mapping proved inaccurate. The mean minimum distance between sites was 5.6 km, which is considerably greater than the diameter of male koala home ranges (e.g. 60 ha) [21].

Acoustic sampling
At each site, we deployed one Song Meter (SM2+ in 2015, SM4 in 2016/17-Wildlife Acoustics, Maynard USA) to record koala bellows. Song Meters were programmed to record from sunset until sunrise, the peak calling period of koalas [29], with a sampling rate of 22 kHz, and resolution of 16 bits per sample. We deployed a total of 171 Song Meters for at least seven consecutive nights (7-14 nights) over three breeding seasons (September-December) in each year between 2015 and 2017. Table 1. Harvest treatments surveyed with acoustic recorders for koalas and number of replicate sites. Heavy harvests were defined as compartments treated by 'heavy single tree selection' or where more than 80 m 3 per ha of timber were removed in an operation. Koala high-use areas were previously identified based on scat accumulations and resulted in patches (mean = 4 ha) being excluded from harvesting. All treatments were based on GIS mapped layers and field assessment. See S1 The distance at which koala calls can be detected is likely to vary with environmental conditions and topographic position. Studies using call playback at 75 dB in forests have found that call amplitude attenuates to background noise levels in recordings (e.g. SM2 recorders) within 100-150 m [38,39], although we expect under ideal conditions and with later model acoustic recorders/microphones (e.g. SM4) this could extend to~300 m based on detection distances for other species, such as owls [40]. [40] also noted that recording devices detect sounds over shorter distances than human ears and detection algorithms are less efficient for faint calls. Given variation in acoustic recorder models, background noise levels (e.g. during rain) and topographic locations of recorders, we assessed weather, topography and recorder unit effects on probability of detection and accounted for important effects in occupancy analyses (see below).

Automated analysis of koala bellows
Recordings were scanned by acoustic software and a koala recogniser [41]. Recordings matched by the koala recogniser were checked for false positives by manually visualizing spectrograms of the audio and listening to recordings. Random checks were carried out for false negatives, revealing that very faint bellows were not recognised by the software. A single koala bellow comprised multiple event triggers. We defined a koala bellow as sequential event triggers that were <60 s apart. Elsewhere the mean length of a koala bellow has been reported as 36 sec and few last more than a minute [28,42]. The number of koala bellows was tallied per site per night.

Habitat and GIS covariates
A rapid assessment of habitat variables was undertaken within a 50 m radius around each recorder, which represents the central zone of detection of the recorder. Projected foliage cover of the canopy was measured using a smart phone application ('Habitapp' V1.1, Android application). Per cent cover was then apportioned based on a visual estimate to the different tree species comprising the canopy based on a visual estimate of their percentage contribution. Understorey cover was scored on the Braun Blanquet scale (<5%, 5-25%, 26-50%, 51-75%, > 75%) and the dominant contributor (rainforest, Acacia, vines, lantana, eucalypt regeneration) recorded. The presence of trees with hollows was scored as 1. common, 2. rare or 3. absent and stand structure was classed as 1. even or 2. uneven age. Harvest intensity and time since harvest were assessed from the presence of stumps and tree size to confirm the GIS classification of the site. Topographic position was scored on a scale of 1-12 (1 = summit, 12 = swamp). GIS was used to classify the surrounding landscape of each site within a 1 km buffer ( [22]; Table 2) based on a number of GIS layers. This included the extent of wildfire and area harvested in the last 10 years, as well as the extent of old growth and cleared land ( Table 2). We also derived site productivity by calculating normalized difference vegetation index (NDVI, [43]) values using MODIS MOD13Q1 granules acquired for September-December for each year between 2005 and 2015, at a 250 m scale and a 1 km site buffer.

Occupancy modelling
We used single-season occupancy modelling to account for imperfect detection of koala bellows and to estimate probability of site occupancy [34]. We considered single-season occupancy to be appropriate for modelling because, although our study extended over three years, the low metabolic lifestyle of koalas, and the consequent restrictions on activity and reproduction [44], suggested little likelihood for significant change in populations over this short time frame. Moreover, each site was sampled in only one year and we included year as a covariate for occupancy in our modelling (see below). To restrict the number of models, we employed a hierarchical approach [45] whereby we first modelled probability of detection (ρ) from 7-14 consecutive nights of sampling while using a global model for occupancy that included all individual covariates for site occupancy (C). Model selection was carried out using PRESENCE (V10.5, [46]) by comparing AIC scores of each model [47] to identify the most supported covariate for ρ, which was carried forward and used to model all subsequent parameters. Covariates used to model variation in ρ are listed in Table 2. Covariates modelled for C included site-based disturbance (harvesting treatment, wildfire, clearing), habitat (foliage cover of important browse trees) and broad environmental (elevation, NDVI) variables, as well as landscape assessments of disturbance based on GIS layers (Table 1). Supported candidate models (dAIC < 2 points) were model-averaged to provide estimates of all parameters. All continuous covariates were standardised prior to analysis.
A Canonical Analysis of Principal coordinates (CAP) was conducted using PRIMER 6 (V6.1.18) and PERMANOVA+ (V1.0.8) (PRIMER-E Ltd, Plymouth, UK) to explore the relationship between different levels of occupancy and cover of individual tree species at each site, as well as supported covariates from modelling. Conditional occupancy for each site (probability that a site is occupied, given its particular detection history) was classified as probably absent (< 5%), low likelihood of occupancy (5-25%) and occupied (100%). The data were normalised and a Euclidean distance dissimilarity matrix was used.

Results
Passive acoustic surveys were completed over three years at 171 sites resulting in~14,640 hours of recording across 1,464 nights. Site and landscape attributes for each of the eight forest disturbance treatments are shown in S1 Table. Across the 171 sites, 2513 validated koala bellows were recorded at 62% of sites (naïve occupancy). On average, 1.6 bellows night -1 were recorded across the 171 sites surveyed. Koala high-use areas supported nearly three times the bellow rate (3.1 bellows night -1 ) as other treatments, but an ANCOVA found that the difference among treatments was not significant (F 7,162 = 0.82; P = 0.6). Minimum nightly temperature was a significant covariate (-ve), indicating more frequent bellowing at lower temperatures (F 1,162 = 5.93, P = 0.02). Bellow rate also did not vary with topographic position (F 1,162 = 0.47, P = 0.62); i.e. among ridges, mid-slopes or gullies.

Detection probability
In all, 40 candidate models were fitted for detection probability (S3 Table). Four models were supported, with model 1 allowing detection probability to vary with year and minimum nightly temperature, whereas models 2-4 also allowed detection to vary with nightly rainfall, moon phase and topographic position (Table 3).
Minimum nightly temperature had a major negative influence on detection probability ( Fig  2). Detection probability declined from an estimated 0.57 per night at 3˚C minimum temperature to 0.32 at 23˚C. Minimum nightly temperature was also correlated negatively with month of survey, indicating lower detectability in December compared with September. Detection probability was also influenced by year of sampling, which was approximately 10% lower in 2015 (SM2 units) compared with 2016 and 2017 (SM4 units), respectively (Fig 2). The other covariates (influence of nightly rainfall (-ve), moon phase (-ve on full moons) and topographic position (+ve on ridges)) made only minor contributions to the models. For rainfall this is likely to be because little rain fell during the study.

Occupancy
The covariates in the most supported model for detection probability (i.e., year and minimum nightly temperature) were carried forward to account for imperfect detection when modelling occupancy. In all, 35 candidate models were fitted (S2 Table). Forest treatment (trt) did not feature in any of the supported occupancy models (Table 4). For example, after holding other supported covariates at their mean, estimated occupancy across the full range of treatments, including old growth, was 0.64±0.04 (±SE). Similarly, landscape assessments of the extent of disturbance within the 1 km buffer did not influence occupancy. In addition, there was no support for occupancy to have been influenced by land tenure or to have changed over the three survey years.
Modelling revealed support for six models with DEM 2 included in each model (Table 2). Other supported covariates included as additive or interactive effects were important feed tree cover, NDVI 2 and wildfire extent. The data were considered to be a reasonable fit to the  Of the six supported models for koala occupancy, DEM 2 showed the strongest association (negative) with occupancy, declining from an estimated occupancy of 0.75 at 10 m ASL to 0.13 at 1,327 m ASL (Fig 3). However, DEM 2 plus the foliage cover of important browse trees had almost twice the AIC weight as DEM 2 alone. All other associations (browse tree cover, NDVI 2 and wildfire) were positive, though the strength of associations varied (Fig 3). Occupancy increased slightly (~10%) with both cover of browse species and NDVI. Occupancy was variable at low values of NDVI (generally from open forest in the Richmond River valley-e.g. Carwong S.F.), but was consistently higher at sites with high NDVI. Occupancy also declined slightly (~5%) with greater extent of wildfire, though we sampled few sites with a large wildfire extent.
Two interactive effects were also supported in modelling, although neither were in the top model. The interaction of DEM 2 and browse tree cover was such that probability of occupancy for koalas was high (>0.7) at low elevation even when there was low cover of important feed trees. As an example, some sites (e.g. Whian Whian SCA) had almost 100% cover of blackbutt E. pilularis (low importance browse tree), but high numbers of bellows. In contrast, at high elevation no sites were sampled with high cover of important browse trees, and occupancy was  predicted to be close to zero when these were absent (Fig 4). Similarly, the interactive effect of NVDI 2 and DEM 2 indicated that at low elevations very low values of NVDI 2 were not recorded, and that low-moderate values were associated with moderate koala occupancy (0.5-0.6) compared to negligible koala occupancy under both low and medium NDVI 2 at high elevation (Fig 4).

Koala browse species
More than 42 tree species were identified at the survey sites (S2 Table). The most widely distributed tree species was tallowwood Eucalyptus microcorys, which is a primary browse species for koalas in the study region, being recorded at 120 sites, both occupied (75%) and likely absent sites (52%) (S2 Table). The CAP indicated that occupied sites were associated with higher cover of a wide range of tree species, such as E. microcorys, E. saligna and grey gum species, as well as important timber species (e.g., E. pilularis, ironbark and spotted gums) (Fig 5). Occupied sites were commonly at lower elevations and with greater NDVI 2 . Occupied sites at higher elevations were associated with the cover of 'other' browse and supplementary browse tree species (see also S2 Table). Sites where koalas were likely absent (<5% probability of occupancy) were associated with higher elevations (DEM 2 ), low cover of important browse tree species and high cover of E. nobilis, E. obliqua and E. campanulata and other non-browse species (Fig 5). Sites with low likelihood of koala occupancy (5-25%) were associated with a greater extent of recent wildfire (within a 1 km buffer) and also more cover of spotted gum species (Fig 5; S2 Table).

Discussion
Our study used innovative technology to automate identification of passively recorded koala calls and has provided new insights into the resilience of koalas to timber harvesting at a local landscape scale, as well as their status in hinterland forests of north-east NSW. Much higher detection rates (naïve occupancy) were recorded than previously known for the study region [25,26]. After correcting for imperfect detection, occupancy did not vary with timber harvesting treatments, nor did bellow rate. Although we targeted modelled moderate-high quality koala habitat for our surveys, this represents an extensive area of north-east NSW, covering 1.66 million ha, approximately 25% of which is public timber production estate [19]. Koalas were broadly distributed throughout this area across many forest communities, comprising different browse tree species, and on the full range of topographic positions. Three covariates, and their interactions, were supported as influencing koala occupancy, including elevation, the cover of important browse trees and site productivity (NDVI). Extent of recent wildfire also had a minor influence on occupancy. Outside of modelled moderate to high quality habitat, a range of factors further reduce habitat quality for koalas such as high frequency wildfire, lack of browse species and low productivity soils on steep terrain [19].

Bellow detection and interpretation of occupancy
Two primary covariates influenced koala detection probability. We found minimum nightly temperature had the greatest (negative) influence on detection probability and this was  correlated with month of survey during the breeding season, with lower detectability in December than September. Similarly in Queensland, most male bellows are recorded during peak mating season in spring and there is a negative relationship with maximum nightly temperatures [29]. This pattern could reflect an energetic constraint of bellowing at warmer temperatures or bellow rate may simply reflect a correlation with the spring breeding period of koalas, where koalas call more frequently at the beginning of the season when temperatures are cooler and less at the end when temperatures are warmer [29]. We also found detection probability was influenced by year of sampling, being 10% lower in 2015 when SM2 units were deployed compared with 2016 and 2017 when SM4 units were used. SM2s were recently confirmed to have lower detectability and detection distances among a range of sound recorders tested due to less sensitive microphones [40]. Other covariates had a minor influence on bellow detection. We adjusted occupancy for the effects of variable detection probability (temperature and year/ acoustic recorder) and this is a key step needed when estimating occupancy. Occupancy is an appropriate surrogate for population assessment when density is low as is expected for koalas [32,35]. For example, given SM2/SM4s sample a 150-300 m radius of forest for koalas [38,39], and male koala density in better quality forests of north-east NSW is 0.03 per ha [32], then just a single male would be expected in the acoustic sample area where a bellow was recorded. In addition to occupancy, bellow rate is correlated with koala density [39]. Bellowing also serves both to avoid male-male interactions by functioning as a signal of body size [48,49] and to attract females [29], and thus an increase in bellow rate likely reflects an increase in breeding activity in a population [29,50].

Koalas and forest disturbance
Our study sampled a broad range of timber harvest intensities and times since harvesting, at both site (~300 m radius) and a larger landscape scale (1 km buffer), together with old growth forests for comparison. Neither occupancy nor bellow rate was found to be influenced by any of these treatments. These results are consistent with previous studies that have suggested koalas tolerate selective harvesting [26,30,32]. We also found occupancy and bellow rate were not lower in recent, heavily harvested forests after a significant component of the canopy had been removed. At the time of our surveys these sites were dominated by dense regeneration of sapling eucalypts in the understorey (mean = 5 years post-harvest). Although intensive harvesting in mixed species E. pilularis forests can favour the dominance of this species in the regeneration [36], koala occupancy remained high decades later in the old, heavy harvested stands, including those dominated by E. pilularis. Intensive harvesting of wet sclerophyll forests dominated by E. microcorys and E. saligna does not influence tree species diversity in the regenerating forest [51,52].
Resilience of koalas to recent, heavy harvesting is most likely explained by the landscape mosaic of forest types and disturbance history in north-east NSW; especially the level of harvest exclusion in the landscape. Over the last 20 years exclusions averaged~40% of the State forest area in the region [7]. In our study, about 50% of the 1 km area surrounding our recent, heavy harvest sites received this treatment in the last 10 years. The remainder comprised temporary off-set zones, but also permanent riparian buffers, old growth and rainforest exclusion areas and habitat protection for owls. In addition, large trees (40-80 cm dbh) provide important shelter and browse for koalas [31,32]. Within the harvest area, scattered habitat trees, recruit/seed trees and feed trees for other species assist in providing a scattered uneven age structure, even where harvesting is heavy [7]. We recorded an average of 27% overstorey cover in recent, heavy harvest sites compared to 48% cover in old growth and old, heavy harvest sites (S2 Table).
Koala high-use exclusion areas represent one method of retaining patches (mean size = 4 ha) of large browse trees where an accumulation of koala scats had been identified prior to harvest, but we found occupancy was no greater in these areas than other treatments, and although bellow rate was greater, the difference was not significant. It is possible that retaining browse trees as clumps across a harvested zone might be more beneficial to koalas than focusing harvesting exclusion in a single area. Limited radio-tracking near Eden NSW has shown koala home ranges can comprise a mosaic of regrowth and unlogged habitat [53]. Foliage nutrient concentrations, including nitrogen, decrease with tree age in E. grandis and E. pilularis [54], and also tend to be higher in younger eucalypt foliage [55] or in foliage exposed to light [56]. In contrast, plant secondary metabolites are higher in foliage on small compared to large E. microcorys trees [36]. While this highlights the complex trade-offs between nutrients and toxins in browse resources, young trees, for example in plantations containing preferred browse species, can support high densities of koalas [20,57]. More detailed fine scale movement data of koalas in a post-harvest landscape is needed to assess the effectiveness of different retention approaches, as well as assessing the extent to which koalas use young regenerating trees and exclusion areas post-harvesting.
Wildfire may have a greater immediate impact on koalas than timber harvesting both through direct mortality and indirectly by burning the canopy [58]. We found that the local extent of a wildfire in the last 10 years (usually a single fire) had only a small negative effect on occupancy. This is consistent with [21] who found that resource depletion from a major wildfire is short term for koalas because their mobility allows rapid recolonisation of the burnt forest, and they can maintain home ranges within forest regenerating from fire. Wildfire frequency rather than a single major fire may be a greater threat to koalas, noting that the areas we sampled had infrequent wildfires in the last 10 years. Indeed, wildfire frequency was a major contributor to modelling of koala habitat suitability in our study region [19].

Key drivers of koala occupancy in NSW's north-east forests
Rather than disturbance, the main drivers of koala occupancy were elevation, browse tree cover and landscape NDVI, as well as their interactions. Koala occupancy declined non-linearly with elevation and remained >0.5 at 700 m ASL, such as on the Dorrigo plateau. Lowest occupancy was found above 1,000 m ASL on the New England tablelands, where modelling of koala habitat also predicted high suitability areas to be scarce [19]. An association with low elevations has long been known (e.g., [17,26,59]); however, high occupancy at mid-elevation and even some high elevations (e.g., Nowendoc) appears to be less widely appreciated (but see [60,61]). The New England Tablelands (and the north coast NSW) are predicted to provide climate refugia for koalas under climate change scenarios [62], but the interaction of elevation with important browse tree cover suggests many high elevation areas are unsuitable for koalas because less preferred browse species dominate (e.g. E. nobilis, E. obliqua, E. campanulata) and/or greater concentrations of plant secondary metabolites occur at higher, colder elevations [36]. While koalas were once considered to be "in great numbers" on the tablelands in the early 1900s [63], it is likely that many of the more suitable forests on fertile soils were preferentially cleared for agriculture and pine plantations.
Unsurprisingly, we also found occupancy increased with cover of important browse trees at a site. Koalas are obligate folivores that specialise on a small suite of eucalypts, though the use of these species may vary because of differences in site productivity or because the availability of desirable tree species varies between sites [17,18,64]. Browse tree selection at a local scale is complex and substantial differences are found in the amounts and types of chemical compounds in leaves of neighbouring Eucalyptus trees, even between trees of the same species [65]. We found an interaction between browse tree cover and elevation, which means occupancy was high at low elevation even where cover of important feed trees was low. At high elevation, occupancy was predicted to be almost zero when important feed trees were absent. This could mean that tree species not commonly recognised as important browse species are used more than currently expected at low elevations. Alternatively, low elevation sites with few important browse trees may be an artefact of our site-based assessment of browse tree cover in an area well known for its fine-scale mosaic of tree species and forest communities compared to that found at higher elevations [66]. Finally, koala occupancy also increased with NDVI, reflecting the importance of site productivity to browsers [43]. NDVI is also correlated with the presence and abundance of another eucalypt folivore, the greater glider Petauroides volans [67].

Estimating the koala population in the north-east forests
Expert opinion currently suggests that the north coast and tablelands of north-east NSW supports a koala population of 8,367 (4,048-14,618) and 2,771 (468-5,838), respectively, though with high uncertainty (72-92% uncertainty) [15]. A minimum estimate based on our data for public-owned forests alone, suggests the population is substantially greater at more than 14,250 koalas (assuming a 1:1 sex ratio). We base this conservatively on just one male koala being present at occupied sites, a larger than expected sample area of 500 m radius surrounding acoustic recorders (78.5 ha), our established negative relationship with elevation and a range area of 0.83 million ha of predicted moderate to high quality koala habitat in public forests [19]. This estimate does not consider additional very low density koala populations in the extensive area of predicted low quality habitat. More importantly, it does not include privately-owned forests, which we did not survey, though they represent 50% of the predicted koala habitat for north-east NSW [19]. Systematic surveys are needed to estimate occupancy in private forests. We expect occupancy would be lower in many private forests than our estimate for public forests as they are often more fragmented by areas of permanent clearing. In particular, koalas in small forests near urbanisation are especially prone to direct mortality from cars and dogs [16,22,57,68]. However, the value of coastal forests near urban centres, some of which are privately owned, should not be overlooked as they may support higher koala densities than hinterland forests when occurring on more fertile soils [32].

Conclusion
In conclusion, koalas persist at high rates of occupancy and have a similar bellow rate across different timber harvest intensities and time since harvest in north-east NSW. Retention forestry likely plays a significant role in mitigating harvesting impacts on koalas, and other biodiversity [4,5], but further research on movements would assist in revealing how exclusion areas and regeneration are used, and in optimising prescriptions for the species. In addition, passive recorders have revealed a widespread and large, but likely low density koala population in north-east NSW. This is consistent with the view that such populations have been poorly detected in less accessible parts of NSW and that koalas may be more resilient than suspected in those areas [69]. Our results also highlight that many National Parks and State forests are currently important custodians of koala populations in north-east NSW. The species continues to decline in many areas [16], but there are few data describing population trends. We suggest passive recorders have great potential for monitoring this and other low density populations of cryptic, but vocal animals.
Supporting information S1 Table. Site attributes for eight different forest disturbance treatments. Treatments were classified by harvest intensity and time since harvest. (DOCX) S2 Table. Tree species recorded in a 50 m radius around each acoustic recorder (n = 171 sites) in north-east NSW and their percentage occurrence at sites classified by conditional occupancy values of koalas. (DOCX) S3 Table. The candidate models fitted for detection probability and for probability of koala occupancy using acoustic recorders to detect koala bellows. (DOCX)