Long-term growth data were sourced from the Department of Agriculture and Fisheries, Forestry Science database, combining permanent sample plots in State Forests (described by Ngugi et al., 2015 and tested for representativeness by Beetson et al., 1992) with long-running thinning trials in Queensland State Forests (Table 1). This mix ensured we captured forests with and without harvesting effects (including non-commercial thinning) measured over a minimum of 5 years, and up to 79 years. The earliest measurements began in 1931 and the latest measurements were made in 2011. Individual trees were tagged with a unique identifier. Plot sizes varied from 0.0405 to 0.5 ha (mean ± standard deviation 0.179 ± 0.145 ha) and were generally rectangular in shape. Given the age of many plots, their dimensions used imperial units, now rarely used, e.g., 0.1 acre (0.0405 ha) is limited by a 1 chain × 1 chain plot with 1 chain, equal to 66 feet or 20.1168 m. No plots in the dataset were completely free of tree removal since European settlement, as far as we can determine. We calculated the contribution of tree species to stand basal area (Eq. 1) and grouped the plots by dominant species (that contribute most to stand basal area in Table 1).

We calculated basal area BA (m² ha^–1) and stem density N (ha^–1) as the sum of individual trees with DBH of at least 10 cm (Eqs. 1, 2).

Where DBH is the diameter at 1.3 m height above ground (cm) measured over-bark, A is the plot area (ha) and “count” is the number of trees per plot. DBH was measured at 4.5 feet above ground until 1974, which is equivalent to 1.37 m. As there are no suitable taper functions available for all tree species to correct this, we decided to not modify our data. Our paired analysis (see section “Effect of Management on Stand Structure”) ensures that this effect has no effect on our analysis.

We also calculated quadratic mean diameter (DG, cm) as a measure of mean tree size.

For every period (years) with subsequent repeated measurements, we calculated the basal area increment of existing trees (BAI), recruitment (BA₁ = zero), and losses due to mortality and harvesting (BA₂ = zero) based on basal area (m² ha^–1 year^–1) (Eq. 4).

Sapling trees were counted as recruits when they reached 10 cm DBH. Measurement protocols changed over time (Beetson et al., 1992) and differed among the different types of plots (i.e., permanent sample plots vs. thinning trials). In permanent sample plots (i.e., those described by Ngugi et al., 2015) the minimum tree size for recruitment used from 1975 to 1991 was 7.5 m in height, which is approximately equivalent to a DBH of 10–15 cm. After 1991, on all plots all trees ≥ 10 cm DBH were measured. Such differences in measurement protocols were accounted for in our analysis by splitting the results by first and last measurements (see section “Effect of Management on Stand Structure”). Further, stems < 10 cm DBH were measured in some plots, but were excluded from analysis to ensure we retained comparability among plots. Trees ≥ 10 cm DBH are generally considered to be a permanent part of the stand and able to survive effects of low intensity fire, browsing and drought (Hurst et al., 2011; Lawes et al., 2011).

Mixed-species native forests in southern Queensland cover an area of around 9.5 million hectares (Figure 1). About one quarter (2.6 million ha) are state forests and the remainder are private forests, national parks and reserves (Ngugi et al., 2015; Lewis et al., 2020). For each plot, we identified annual rainfall and annual average temperature between 1970 and 2000 (Fick and Hijmans, 2017). The great majority of plots were located south of the Tropic of Capricorn and all plots received annual rainfall greater 500 mm (Figure 1). The regions has a humid sub-tropical climate and the plots have an average daily temperature of 20°C (min-max 16–22°C) and annual rainfall ∼ 960 mm (min-max 572–1620 mm year^–1; Fick and Hijmans, 2017).

A large number of plots were unthinned and had not been harvested during the measurement period (Table 1). Using three models for predicting self-thinning (i.e., similar to the approach used by Trouvé et al., 2017), we sought to identify maximum stocking using only unmanaged plot data. First, we assessed predicted self-thinning on the basis of basal area (BA = N × π/40,000 × DG²) applying a self-thinning line with slope of –2.000 (ST_BA). Secondly, we assessed predicted self-thinning on the basis of stand density index (SDI = N × (DG/25)^1.605; Reineke, 1933), assuming a slope of –1.605 (ST_Reineke). For both models we calculated the 95th quantile (P95) of BA and SDI for the four most frequent species, E. pilularis, E. crebra, Callitris glaucophylla, and Corymbia citriodora ssp. variegata (see Table 1 and Figure 1). Finally, we used quantile regression to fit self-thinning lines with variable slope (ST_variable) to the point cloud, again using 95th quantile. We compared these lines to those with fixed slopes derived from either density (ST_Reineke) or basal area (ST_BA). We thus evaluated self-thinning without prior assumption of fixed slopes for self-thinning lines (Trouvé et al., 2017). The 95th quantile of BA is an easy-to-measure proxy for maximum potential stocking.

Individual tree information included measures of status (alive, standing dead or removed/ring-barked). For standing dead trees, the reason for tree death was assessed. We assigned all plots, as being either managed (n = 130) or unmanaged (n = 175) and used only the first and last measurement for each plot. We excluded plots with measurement period less than 10 years (n = 10), to allow for time needed to respond to management. Managed plots were subject to non-commercial thinning of small trees or ring-barking of undesired trees up until 1990, as well as selective harvesting on a 20–30 year cycle. Managed plots benefitted from fire protection. We cannot preclude that harvesting of nearby trees outside the plots changed growth conditions within the plots. Most plots (76%) had a buffer to accommodate such effects. No trees had been removed or killed in unmanaged plots (as far as we can determine). For every plot we calculated stem density and assessed the skewness of DBH distributions. We assumed the same year-to-decade variation due to climate change, CO₂ fertilization or nitrogen deposition (Kirschbaum, 2004; Ruiz-Benito et al., 2014; Adams et al., 2021) for all plots. We checked for differences in rainfall (Supplementary Figure 2) and sampled area (Supplementary Figure 3) between managed and unmanaged plots. Average rainfall on unmanaged E. pilularis plots was 180 mm greater than that on managed counterparts (Supplementary Figure 2). There were no other significant differences. Temporal variations in climate or fire frequency were beyond the scope of this study.

We used QGIS, version 3.2.2, for preparing Figure 1 and R statistical software, version 4.0.3 (R Development Core Team, 2021) for processing and analysis. We used analysis of variance to determine whether management has had an effect on stand structure, increment and recruitment and whether managed plots have different plot area or precipitation than unmanaged plots (see previous sections). We assessed both parametric statistical tests (Student’s t-tests with logarithmic transformation) and non-parametric tests (Kruskal-Wallis tests). Both accommodated skewness in the data and produced similar results; results from Kruskal-Wallis tests were reported unless otherwise noted.

Sample plots used in this study were distributed over a large proportion of southern Queensland (an area of about 9.5 million hectares, slightly greater than the area of Austria; Figure 1) and span a large gradient in rainfall (Figure 2). Eucalyptus pilularis is a dominant tree near the coast, where rainfall exceeds 1,200 mm per year. Rates of growth of E. pilularis (average basal area increment BAI 0.84 m² ha^–1 year^–1) were much faster than those of other species (see also Supplementary Figure 1). Corymbia citriodora ssp. variegata dominates at higher rainfall in coastal areas, but also on sites further inland with annual rainfall < 600 mm. Rates of growth of C. citriodora ssp. variegata decreased (mean BAI 0.23 m² ha^–1 year^–1) with precipitation (Supplementary Figure 1). Eucalyptus crebra (mean BAI 0.09 m² ha^–1 year^–1) dominates large areas in the 600–700 mm rainfall zone, but grows slowly. Callitris glaucophylla (mean BAI 0.23 m² ha^–1 year^–1) is increasingly dominant as annual rainfall declines from 700 mm. Callitris glaucophylla maintains rates of growth similar to those of C. citriodora ssp. variegata (at higher rainfall). Species distributions are clearly influenced by climate, as well as topographic and soil conditions (not analyzed here).

Fifty-five plots recorded negative BAI (Figure 2); mostly associated with drought years between 1979 and 1983 (the Eastern Australian drought or “great Australian drought,” associated with a strong El Niño event, Gibbs, 1984; Allan et al., 1991).

Stocking of unmanaged plots was clearly affected by self-thinning (Figure 3). Increases in tree size (based on quadratic mean diameter, Eq. 1) were accompanied by strong reductions in stem density. All variants of self-thinning models (ST_Reineke, ST_BA, ST_variable) showed similar results. The slope of the self-thinning line (–1.662) of E. pilularis was similar to Reineke’s SDI. Self-thinning was thus weakest (density declined more slowly with increasing diameter) in E. pilularis among our four target species (Figure 3A). Stronger self-thinning was observed in plots in drier areas (cf. Figure 2; slope for Callitris glaucophylla –2.021; slope for C. citriodora ssp. variegata –2.253). Moderate self-thinning in plots dominated by E. crebra (slope –1.795) is conditional of fewer observations for this species (Figure 3D). While many managed Callitris glaucophylla plots have already reached the self-thinning line (potential stocking), managed E. pilularis plots are still below potential stocking (Figures 3A,C).

Basal area (BA) provides a good index of maximum stocking, taking into account self-thinning (see Figure 3). As with BAI (Figure 2), P95 BA also varied by species. P95 BA was 56.3 m² ha^–1 for E. pilularis, 19.0 m² ha^–1 for C. citriodora, 16.8 m² ha^–1 for Callitris glaucophylla and 14.5 m² ha^–1 for E. crebra.

In Figures 4–6, we contrast first and last measurements of repeatedly measured sample plots (n = 305, with a measurement period of at least 10 years). There were 10 plots without remeasurements. Initially, stem density was greater (p < 0.001) in managed plots (n = 130) than unmanaged plots (n = 175), but there was no difference in basal area (Figures 4A,C). After an average interval of four decades (min. 10, max. 79 years), stem density had clearly declined in E. pilularis plots, but was constant or had increased in plots dominated by other species (Supplementary Figures 2–5). At the time of last measurements, basal area of unmanaged plots was greater than that of managed plots (p = 0.032; Figure 4B). Stem densities were not statistically different between managed and unmanaged plots (Figure 4D). E. pilularis and Callitris glaucophylla were major contributors to this pattern (Supplementary Figures 4, 6) as the two other species showed no differences in basal area or stem density (Supplementary Figures 5, 7).

DBH distributions varied markedly between the first and last measurements. There were significant differences between managed and unmanaged plots (Figure 5) albeit data were heavily skewed. The skewness at the last measurement was greater for the unmanaged plots (1.225 vs. 0.839 for managed plots), while initially greater in managed plots (0.657 unmanaged vs. 0.722 managed). Initial significant differences in stem densities (10–20 cm, 30–40 cm, > 60 cm DBH) between managed and unmanaged plots (p < 0.05, Figure 5B) were not sustained to the last measurement (Figure 5A). At the last measurement, only the 40–50 cm DBH class showed statistically greater stem density (p < 0.001) in unmanaged plots. A similar trend (although less significant) was evident for the 50–60 cm DBH class. These results suggest that early stand structure is more important than management for determining long-term density of large trees. The density of larger trees (> 60 cm DBH) increased between the first and last measurements for all plots (Figure 5). Hence, our first hypothesis, that the density of larger trees would be reduced in harvested stands, was not supported.

Figure 4 shows change in basal area of stands that increase due to growth and recruitment and decrease due to mortality or harvesting. BAIs (i.e., growth rate of repeatedly measured trees) were significantly smaller in unmanaged plots than managed plots at the last measurement (p < 0.001; Figure 6B). Negative BAIs were observed in some plots (Figures 6A,B). While there was evidence of greater recruitment in managed plots—in particular at the last measurement—the data were highly skewed and differences between managed and unmanaged plots were not highly significant with p = 0.003 (Figure 6D). BAIs at the last measurement were considerably less than those recorded initially, irrespective whether plots were managed or unmanaged. Our results supported our second hypothesis, that management increased recruitment. Our third hypothesis (i.e., management reduced increment) was not supported by the results of our study. We note the importance of initial stand conditions. As differences in BAI and recruitment were less pronounced at the initial measurement relative to the last measurement (based on levels of significance), we conclude that management had (unexpected) effects on stand structure and dynamics.

DBH distributions in Queensland forests were increasingly skewed over the study period. Recruitment increased the density of trees < 20 cm DBH. Management, most likely through thinning, result in more evenly distributed tree sizes from 10 to +60 cm DBH and reduced the skewing of density (skewness for managed plots = 0.839, skewness for unmanaged plots = 1.225; see Figure 7). Fifty years without active management resulted in changes in distribution of tree size (Figure 7B).

Assessed across all measurements, BAIs in managed plots (0.294 m² ha^–1 year^–1) were on average less than those of unmanaged plots (0.333 m² ha^–1 year^–1, Table 2). This pattern matches the analysis of first measurements (Figure 6A), but is contrary to that based on last measurements (Figure 6B). Existing trees (BAI) contributed at least an order of magnitude more to BA than recruitment (i.e., > 90% of basal area gains were due to growth of existing trees). Managed stands lost significantly more BA than unmanaged stands (p < 0.001) with harvesting accounting for 74% of the total loss (Table 2). There was a positive net change in basal area (gains due to increment and recruitment exceeded losses due to mortality and harvesting, Table 2) irrespective of whether forests were managed or unmanaged. Increases in BA were almost twice as great for unmanaged compared to managed forests (p < 0.001), largely due to harvesting removals in the latter. In turn, the rate at which unmanaged forests approached potential stocking (see self-thinning lines in Figure 3), was twice as fast as that of managed counterparts.

Native forests dominated by Eucalyptus and Callitris are limited to Australia. BAIs recorded here are mostly small by comparison with other Australian forests and forests elsewhere. Eucalyptus pilularis grew most quickly with an average BAI of 0.84 m² ha^–1 year^–1. This species is often dominant where precipitation is at least 1,400 mm per year. Average BAI for 3–7 year old regrowth of Eucalyptus regnans in Victoria was 4.8 m² ha^–1 year^–1 (Attiwill, 1992). BAI for Eucalyptus obliqua in Victoria was recorded (Attiwill, 1979) as 1 m² ha^–1 year^–1, while that for E. fastigata in New South Wales was 1.14 m² ha^–1 year^–1 (Bi, 1994). Coastal mixed eucalypt forests in New South Wales (mean BAI 0.27 m² ha^–1 year^–1; Curtin, 1970) and coastal spotted gum (Corymbia maculata, mean BAI 0.11–0.54 m² ha^–1 year^–1; Bauhus et al., 2002) show similar rates of growth to the forests studied here. Within Queensland, BAIs of dry sclerophyll spotted gum-ironbark (Corymbia citriodora—Eucalyptus siderophloia) forest ranged from 0.22 to 0.34 m² ha^–1 year^–1 (Grimes and Pegg, 1979). Higher temperatures largely account for reduced rates of growth in Queensland compared to more southern states of Australia. Average annual temperatures for study plots are mostly greater than 16^°C (average 19.7^°C) and growth increments generally reduce sharply once average temperatures exceed 15^°C (Prior et al., 2011). Internationally, there are limited published examples of native forests that accumulate basal area at similar (slow) rates to those reported here. Native Prosopis forests in Texas and Argentina (Felker and Guevara, 2003) and the miombo woodlands of sub-Saharan Africa (Pelletier et al., 2018; Chidumayo, 2019) share some similarities.

BAIs also varied with precipitation (Figure 2, Eq. 5). Callitris glaucophylla—a native conifer—showed much faster rates of growth (0.23 m² ha^–1 year^–1) than E. crebra (0.09 m² ha^–1 year^–1) at intermediate rainfall (600–650 mm, Figure 2). Callitris glaucophylla is often dominant on deep sandy soils (Horne and Robinson, 1987; Whipp et al., 2012) and the closely related species Callitris columellaris has been shown to more effectively reduce transpiration during drought than sympatric eucalypts (E. microcarpa) (Attiwill and Clayton-Greene, 1984). Ability of drought-resistant C. glaucophylla (Prior et al., 2012) to “conserve” water and maintain foliage in drought periods (in contrast to eucalypts shedding their leaves) may underpin its relatively faster overall rates of growth.

Natural thinning (“self-thinning”) in the absence of disturbances is a widely acknowledged feature of forests as shown in numerous studies (Sterba and Monserud, 1993; Bi et al., 2000; Berger et al., 2004; Trouvé et al., 2017). Self-thinning limits maximum stocking. Such limits also feature in deliberations of critical species abundance and carbon uptake by vegetation (Wilson and Puettmann, 2007; Finkral and Evans, 2008; Verschuyl et al., 2011; Blakey et al., 2016).

While Trouvé et al. (2017) studied self-thinning in temperate Australian forests, ours is among the first studies for sub-tropical forests (see Vanclay and Sands, 2009). We also found clear species differences in slopes of self-thinning lines. For temperate forests, slopes range from –1.90 to –1.52 (Trouvé et al., 2017). The range recorded here for sub-tropical forests is greater (–2.336 to –1.662), likely due to a greater range in rainfall (Figure 2). E. pilularis growing in high-rainfall coastal regions shows similar self-thinning to Victorian eucalypt forests (slope of –1.662; c.f. Trouvé et al., 2017) and to temperate coniferous and broadleaved forests in Europe and North America (Reineke, 1933; Monserud et al., 2004; Pretzsch and Mette, 2008). While reductions in increment and basal area with rainfall are expected, the increasing slope of self-thinning lines are less well appreciated (Figures 2, 3). The steepest slope we recorded was –2.253 for C. citriodora ssp. variegata growing in open forests and woodlands. Absent from our data are unmanaged plots with large diameters for C. citriodora ssp. variegata, that if present may have changed the slope of the self-thinning line. Basal area (slope = –2.000) is arguably a better metric for delineating species “boundary curves” than Reineke’ SDI (slope = –1.605) for C. citriodora ssp. variegata. Reineke used a hand-drawn line, while MacKinney and Chaiken (1935) using the same data and statistical methods, determined a slope of –1.707. Further analysis and more data for forests with rainfall between 1,200 and 1,400 mm (see Figure 2) are needed to confirm that self-thinning slopes in drier Australian forests show consistent differences with forests elsewhere. We note that basal area is a potent and simple to measure index for self-thinning, using Bitterlich sampling (“sampling by diameter” or “angle count sampling,” Eastaugh and Hasenauer, 2013).

Our study showed that unmanaged stands have higher overall stem densities than managed plots, particularly where DBH was < 20 cm and > 40 cm (Figures 4, 5, 7). Density of trees with DBH 20–40 cm increased more on managed plots. These differences are likely due to thinning prior to ∼1980 in State Forests (Norman et al., 2004), removing mostly trees < 20 cm DBH. The large increase in stem density of the 10–20 cm DBH class in unmanaged plots has to be viewed considering the combined effects of abandoning silvicultural thinning as common management practice and changes in policy from 1999. Prior to 1999 trees >60 cm DBH were selectively harvested (Norman et al., 2004). The 1999 South-East Queensland Forest Agreement extended harvesting to trees > 40 cm DBH. This would have presumably promoted recruitment of trees < 20 cm DBH (Henry and Florence, 1966). Data reported here show recruitment was greater in managed plots relative to unmanaged plots (Figures 6C,D) and this is supported by studies of other forest ecosystems (Espelta et al., 1995; Albrecht and McCarthy, 2006; Ouédraogo et al., 2011; Seiwa et al., 2012; Hu et al., 2020). We add a note of caution in that recruitment data are skewed.

Increasing density of trees < 20 cm DBH will likely result in dense stands, slow growth and increased mortality of individual trees (Vanclay, 1988; Bauhus et al., 2002; Lewis, 2020). Selective thinning of small to middle sized trees, can reduce competition and hence promote growth of larger trees that may grow to become important habitat trees (e.g., Jones et al., 2015; Lewis et al., 2020; Bhandari et al., 2021). The presence of large and tall habitat trees and diverse canopy cover are particularly important for native forest management (Franklin et al., 2000; Norman et al., 2004; Bollmann et al., 2005; North et al., 2017) with implications for timber and habitat value (Wormington et al., 2002; Norman et al., 2004). In Queensland, arboreal mammals have been the focus for determining the minimum number of habitat trees in relevant “Codes of Practice” (Department of National Parks Recreation Sports and Racing, 2014; Queensland Parks and Wildlife Service, 2020). Suggested lower limits for arboreal mammals in Queensland are four hollow-bearing trees ha^–1 (Wormington et al., 2002). The presence of hollow-bearing trees was not assessed on the plots used here. Wormington and Lamb, 1999) reported that tree hollows become frequent in trees with a DBH greater 60 cm, depending on species. In this study, density of trees > 60 cm DBH was 7.8 ha^–1 for unmanaged plots and 3.6 ha^–1 for managed plots, on average (Figure 7). Both managed and unmanaged plots had similar initial densities of 1.0 ha^–1. Thus, although we rejected our hypothesis that harvesting reduced density of potential habitat trees, we note that there was a trend of increased density of these trees in unmanaged plots.

Unmanaged plots exhibited a large net increase in basal area (+0.321 m² ha^–1 year^–1) over five decades. This suggests that these plots were previously managed or disturbed and were certainly not in equilibrium (i.e., growth and recruitment balanced by mortality) during the study period. This underpins the importance of self-thinning lines (at least as rough guides) to quantify deviations from the maximum stocking.