This paper presents a new reconstruction of prehistoric population of Australia for the last 50 ka, using the most comprehensive radiocarbon database currently available for the continent. The application of new techniques to manipulate radiocarbon data (including correction for taphonomic bias), gives greater reliability to the reconstructed population curve. This shows low populations through the Late Pleistocene, before a slow stepwise increase in population beginning during the Holocene transition (approx. 12 ka) and continuing in pulses (approx. 8.3–6.6, 4.4–3.7 and 1.6–0.4 ka) through the Holocene. These data give no support for an early saturation of the continent, although the estimated population following initial landfall was probably greater than previously allowed (comparable with the Early Holocene). The greatest increase in population occurred in the Late Holocene, but in contrast to existing intensification models, changes in demography and diversification of economic activities began much earlier. Some demographic changes appear to be in response to major climatic events, most notably during the last glacial maximum, where the curve suggests that population fell by about 60 per cent between 21 and 18 ka. An application of statistical demographic methods to Australian ethnographic and genetic data suggests that a founding group of 1000–2000 at 50 ka would result in a population high of approximately 1.2 million at approximately 0.5 ka. Data suggests an 8 per cent decline to approximately 770 000–1.1 million at the time of European contact, giving a figure consistent with ethnographic estimates and with historical observations of the impact of smallpox, and other diseases introduced by Macassans and Europeans during and after AD 1788.
Quantifying population levels in prehistoric Australia has been at the heart of long-standing debates in Australian archaeology, with most controversy surrounding how rapidly hunter–gatherer populations increased following initial colonization of the continent. The key studies were undertaken between the 1950 and 1990s [1,2], and highlight two main issues: (i) the likely shape of the population growth curve, especially following Late Pleistocene colonization of the continent, and during the Late Holocene when archaeological field evidence suggests rapid and sustained population growth; and (ii) absolute population estimates at either end of this population curve, at initial colonization and at European contact.
Birdsell  was the first to propose a model for saturation of the prehistoric continent by a hunter–gatherer population that rapidly achieved ethnographic population densities. Using modern analogues, he argued that a small founding population could have colonized the continent in only 2000 years—populations ‘budding off’ into new areas as carrying capacity was reached. This model was criticized by Smith  a demographer who considered the rate of population increase (4% per year) too high, and observations that suggested hunter–gatherers rarely exploited their environment to its maximum potential. Using archaeological evidence from north Queensland, Beaton  proposed an alternative model, suggesting populations were consistently low in the Pleistocene, before exponentially expanding in the Mid- to Late Holocene.
Lourandos [5,6] put forward an influential and widely quoted alternative, arguing that the increasing number, density and widespread distribution of archaeological sites observed in the Late Holocene was a result of a socio-economic intensification leading to rapid population growth at this time. In his model, increased competition between groups led to pressure for larger more elaborate ceremonies and intergroup gatherings and alliance networks, and triggered a shift to intensified food production to sustain these networks. Although Lourandos argued for very specific ‘social’ triggers, later work has blurred the distinction between social and economic intensification. Most archaeologists, however, accept that there is evidence for rising populations and a range of related economic and social changes in the Late Holocene (see  for review).
These early models treated demography as unidirectional growth through time. More recent work suggests that the actual population history could be more complicated [8–10]. Davidson , for example, demonstrated that even a small continuous growth rate (0.001%) over the last 40 ka would produce values well in excess of several billion people (2.5 × 1017), a value significantly above the carrying capacity of the continent. In reality, population change can be expected to be both nonlinear and reactive to a range of climatic (e.g. hazards like flooding, cyclones, etc.), and ecological conditions. O'Connor et al.  suggested the last glacial maximum (LGM) would have caused a significant reduction in prehistoric population. The magnitude, rate and direction of these changes are still poorly understood, despite several decades of research providing no definitive answer [7,11–14]. This work, however, has also raised questions about the role of taphonomic loss and/or homogenization of diverse temporal and spatial records in structuring the field evidence [8,15–17].
Less attention has been given to absolute population estimates. Founding populations are generally considered to have been small (less than 150), reflecting a family group or band of people [1,3,18]. However, recent mitochondrial DNA research suggests that a large population (greater than 1000) would have been required to produce the genetic diversity observed [10,19]; and O'Connell & Allen  argue that initial colonization of the continent would have required deliberate organized sea travel, involving hundreds of people.
Population estimates at the other end of the time scale, at European contact (AD 1788), are widely divergent. These are primarily based on extrapolation between a few detailed regional censuses, to arrive at overall estimates of people per km2 in various environments. An early study by Radcliffe-Brown  arrived at the widely quoted figure of 250–300 000 for the continent. These numbers were largely corroborated by Smith , who suggested a maximum population of approximately 314 000. However, Butlin  argued that smallpox introduced by the early Europeans in 1789 would have significantly reduced Aboriginal populations by up to 80 per cent prior to formal ethnographic documentation, and that a maximum number of 1.2 million was more realistic; Campbell  reached similar conclusions, but proposed the introduction of the disease from Macassan contact in northern Australia prior to European contact. Mulvaney & White  reviewed the evidence and suggested a continental population at AD 1788 in the order of 750 000, but this number is largely a compromise between Radcliffe-Brown and Butlin's estimates.
It is now possible to revisit these issues using a comprehensive continent-wide dataset of archaeological radiocarbon data—allied with better techniques for correcting for taphonomic loss  and estimating prehistoric population levels . In this paper, I want to determine: (i) the shape of the long-term population curve, and whether ethnographic population densities were rapidly established following colonization or whether these were delayed until the Late Holocene, and (ii) derive more rigorous estimates of prehistoric population following first landfall and at European contact.
2. Key assumptions
Several assumptions underpin this analysis. I have assumed that the colonization of Australia occurred between 46 and 50 ka [10,20]; while earlier dates have been proposed , these values represent the most conservative age accepted in Australian archaeological literature (approx. 46 ka) and the temporal limits of the dataset (less than 50 ka).
To apply Peros et al.'s  techniques, the population of the continent must have been isolated between two points: Pleistocene landfall and European contact. I therefore, assume that colonization took place in one or several closely timed waves of migration, and no further colonization occurred after this. Some studies hypothesize that multiple waves of migration may have occurred through the Pleistocene  and/or Holocene , but a large body of genetic research suggests no new migration since 42–50 ka [19,29,30] and supports an isolated system. I also assume there has been no significant continent-wide hiatus in the archaeological data over this 50 ka period—either through climatic impacts such as the LGM or through sea-level change. Comparison of the complete dataset with a range of subsets including rockshelter site data (which Surovell et al.  argues is relatively unaffected by taphonomic loss) and archaeological data near Pleistocene coastlines (see the electronic supplementary material) shows a close correlation, indicating that systematic gaps within the dataset are not present.
This paper only deals with the Australia landmass. Prehistoric Sahul, of course, also included the island of New Guinea and the extensive, now flooded, continental shelves between Australia and New Guinea. A working assumption of this paper is that these will have separate demographic histories.
Time-series radiocarbon data is widely used as a proxy for human activity or prehistoric population [13,17,25,31–34]. Analysis and interpretation of this form of data is complex and has several limitations (; see the electronic supplementary material). In this paper, I assume that these radiocarbon data reflect levels of human activity, specifically demographic change (cf. [11,15]).
I have adopted the taphonomic correction of radiocarbon data developed by Surovell et al.  and modified by Williams . These procedures are based on time-decay correction of the dataset using an equation based on a global volcanic dataset, the relevance of which is untested in an Australian context. A specific correction curve has yet to be developed for Australia. However, the good correlation between uncorrected and corrected data (see the electronic supplementary material, figure S12) suggests that the results presented here are robust and not an artefact of correction factors. The electronic supplementary material provides further analysis and discussion on this issue.
3. The dataset
This paper uses the most comprehensive radiocarbon dataset for the continent, assembled by the author. This has been published sequentially (AustArch 1, AustArch 2 and IDASQ [35–37]) and the remaining unpublished data is presented in the electronic supplementary material. The complete dataset comprises 5041 radiocarbon dates from 1750 archaeological sites (figure 1). Of the 5041 dates, 4575 are used in this analysis; 466 are excluded as they lack full information or are regarded as anomalous by their authors.
4. Material and methods
(a) Developing the basic times-series curve
The methods adopted here are based on earlier work by Williams  and Peros et al.  (see the electronic supplementary material, for further discussion). Williams explored the application, use and limitations of radiocarbon time-series data, including the relationship of detrital charcoal to human activity; minimum sample size; the role of the calibration process and protocols for correcting for taphonomic loss. Peros et al. developed techniques to convert radiocarbon data into demographic estimates.
As an initial step, radiocarbon data was calibrated using Oxcal (v. 4.1) . Terrestrial dates were calibrated using Intcal09 and marine dates using Marine09 , with ΔR values after Ulm . Oxcal was used to both obtain a median value for each radiocarbon date (95.4% confidence) and to create sum probabilities. To remove some of the calibration anomalies and allow subsequent analysis, each calibrated date was then ‘data binned’ into 200-year intervals based on its median value.
(b) Taphonomic correction and testing
Identification of archaeological site loss through time and the application of taphonomic correction to time-series data is becoming increasingly common in archaeological literature [16,24,25]. The underlying rationale for correcting the dataset is the expectation that time-decay progressively removes older evidence —distorting these time-series curves. One drawback, however, is that correction factors can also remove real structure in these data.
In this study, correction for taphonomic loss followed procedures in Williams . This involves correction of the actual number of dates per 200-year bin using a decay curve created from a volcanic radiocarbon dataset. The correction equation is 4.1where nc = taphonomically corrected number of radiocarbon dates for the dataset of interest and na is the actual number of radiocarbon dates present at a specific time (t) in the dataset of interest. After Surovell et al. , I ran a series of tests which indicated that while it was necessary to correct open site data, the data from rockshelters was best used without such corrections, which tend to artificially flatten any trends (see the electronic supplementary material). Therefore, the final time-series curve combined ‘corrected’ open site data and ‘uncorrected’ data from rockshelters.
A series of tests were run to identify whether the time-series curve was distorted by heightened erosion during the LGM or by inundation of coastal zones. Tests included comparison of the overall time-series curve with subsets of data from rockshelters, midden sites and those situated on Pleistocene coastlines—each considered to highlight stratigraphic disconformities if present and explore the influence of prehistoric coastal exploitation. Cross correlation of time-series curves showed these factors have little significant effect on the shape of the curve or the time-series trends it described.
(c) Estimating prehistoric populations
Numerical estimates of prehistoric population were derived by taking the ethnographic population of 750 000–1.2 million as an endpoint and fitting the time-series curve to this.
To determine prehistoric population levels, a smoothing spline was first put through the data (degrees of freedom; d.f. = 25) to remove noise and zero. Using the values interpolated by the spline, I determined the annual percentage growth rate (GRAnn), using the equation after Peros et al. : 4.2where in a given pair of consecutive (interpolated) 200-year bins, d2 is the number of radiocarbon dates in the younger bin and d1 is the number of dates in the older bin. Each GRAnn value was multiplied by 0.5 to convert to a percentage and to produce an annual rate from the 200-year bins.
The GRAnn value was then integrated into a compound interest equation to determine quantitative values of palaeo-populations: 4.3where P is final population, f is initially the founding population followed by the P value from each preceding 200-year data bin and t is number of years. Founding populations, colonization dates and the detail of the smoothing spline were fitted to observed ethnographic population levels.
The resulting population curve for prehistoric Australia is shown in figure 2. This shows that although there was widespread dispersal of the hunter–gatherer population following the first landfall, population stabilized at levels much lower than Holocene levels throughout the Late Pleistocene and only began to increase in the Early–Mid-Holocene.
The uncorrected data-binned radiocarbon data of the entire dataset shows a positive curvilinear shape over the last 50 ka (see figure 2 and electronic supplementary material, figure S14). Notable features of the data include low numbers between 50 and 20 ka; a significant increase between 20 and 18 ka, before a decline to some of its lowest values between 16 and 12 ka; and finally constant increase from 12 to 0.4 ka. From 12 ka, the curve begins initially as a steady increase, within which there are a series of pulses (at approx. 8.3–6.6, 4.4–3.7 and 1.6–0.4 ka). While not as clear, several declines or plateaus are also evident through the Holocene (at approx. 5.6–4.8, 2.6–2.2 and 0.4–0 ka).
The taphonomically corrected data (figure 2) effectively exaggerates or subdues the numbers of dates in any given time interval. For reasons explained in , the correction procedure can only be undertaken up to 40 ka, and is increasingly problematic in periods greater than 25 ka. The taphonomic correction plot shows divergence from the uncorrected dataset in the form of a saw-tooth pattern between 40 and 24 ka; a significant peak at 20–18 ka; and short declines at approximately 36 and 26 ka. However, several similar features are evident in both plots, including the steep decline to low values from approximately 18 to 13 ka, followed by the steady stepwise increase from approximately 12 to 0.4 ka. The pulses within the Holocene are presented more clearly in the taphonomically corrected data, showing defined peaks at approximately 8.3–6.7, 4.7–3.7 and 0.6 ka, interspersed with declines at approximately 5.5, 2.6–2.2 and 0.4–0 ka. Effectively, the two curves show the same trend, the continent was rapidly filled during the Pleistocene, but to a level well below levels during the Mid-to Late Holocene.
Further analysis can be undertaken by dividing the data by longitude and latitude (see the electronic supplementary material, figures S15 and S16). Overall the data indicate the greatest amount of human activity was within the eastern half (133–153°) of Australia. There appears to be little difference between latitude, with data evenly distributed across the continent. Both plots indicate low, or ephemeral, levels of human activity across Australia by 45 ka. Between 40 and 20 ka, increasing populations are found in lower and higher latitudes, and focused on the eastern seaboard. The large peak evident in figure 2 between 20 and 18 ka is localized to southeast Australia, and probably reflects the dating of the terminus ante quem of several sites in Tasmania and the switching on of the Willandra Lakes system (; see the electronic supplementary material, for discussion). The decline at approximately 18–13 ka is not particularly evident in any one region, although hiatuses and low numbers are observed in the northern latitudes and western longitudes. In most areas, steady increase in human activity through to the Late Holocene is evident by approximately 12 ka, with some areas (especially in the north and west) not initiating until approximately 8 ka. A highly variable period through the Mid-Holocene (including declines and/or hiatuses in most latitudes at approximately 5 ka) gives way to a stepwise increase from approximately 4 ka, peaking between approximately 1.6 and 0.4 ka. A decline in data from approximately 0.6–0.4 ka is evident in all regions, although this may be an artefact of sampling (e.g. few archaeologists use radiocarbon to date the top levels of their sites).
The GRAnn analysis of the uncorrected dataset indicates an average annual growth rate of 0.01 per cent over the last 50 ka, with a range of between 0.07 and −0.03 per cent (figure 3). This average is robust. Experiments changing the variability of the spline consistently produced an average of 0.01–0.02% over 50 ka, however, increasing the variability of the spline (by increasing the d.f. value) significantly affects the extreme average annual growth rates, with observed values up to 1.84 per cent during some periods.
Following an initial increase in the GRAnn, the Late Pleistocene (50–20 ka) appears highly variable, but generally reflects positive steady growth (figure 3). A notable decline is evident at approximately 36 ka. Following this initial growth, a significant double-dip drop in growth rates occurs in the Terminal Pleistocene at approximately 18.1 and 14.1 ka. From 12 ka, there is strong positive stepwise growth through the Holocene, culminating in the highest growth rates at approximately 2.3 ka. Three noticeable reversals are also evident in the Holocene centred at approximately 6.3 ka, 4.1 ka and the last 2000 years. Taphonomically correcting the GRAnn equation has little effect on the shape of the growth trends over the last 50 ka (figure 4), however, several of the peaks and troughs are intensified (most notably those at approximately 25, 18.1 and 14.1 ka); overall average population growth remains the same (0.01%), but the range of extremes values increases (0.09 to −0.05%).
This analysis implies a founding population size of greater than 1000 and less than 3000 (figure 5). Estimates of effective population size at various time slices are given in table 1 (see also the electronic supplementary material). Figure 5a applies equation (4.3) to a range of founding populations. The initial founding population has a significant effect on the eventual proto-historic values; a starting value of 50 at 50 ka will result in populations of approximately 19 000 at time of contact for example. Based on the observed proto-historic populations [2,3,22], initial populations of between 1000 and 3000 would be required to produce populations of between approximately 380 000 and 1.15 million at time of contact. If the colonization date is closer to 46 ka , initial populations of 2000–5000 would be required to produce the same results (approximately 407 000–1.01 million; electronic supplementary material). Owing to the low variability within the spline (d.f. = 25) used to create the estimates, the detail is too coarse to identify recent declines from disease or other factors. Increasing the variability of the spline (d.f. = 50) incorporates greater extremes from the data providing a more realistic picture of population change (figure 5b). It demonstrates a stepwise increase through the Holocene; a peak in populations at approximately 0.5 ka; and a decline in the most recent period. Using this spline, initial populations of 500–2000 would result in peak populations of between approximately 299 000 and 1.19 million at 0.5 ka before reaching proto-historic values of approximately 275 000–1.10 million.
Observed palaeo-populations could not be duplicated with the taphonomically corrected data, a finding similarly published by Peros et al. . Founding populations of between 30 and 50 000 would be required to reproduce proto-historic values (see the electronic supplementary material, figure S17), numbers well in excess of those considered realistic in archaeological and genetic literature. This failing of the technique probably relates to the younger colonization date necessitated by the temporal barrier of equation (4.2) (40 ka), and increasing uncertainties of the correction method in time intervals greater than 25 ka .
6. Discussion and conclusions
I present here a synthesis of all published (and extensive unpublished) radiocarbon data from documented archaeological sites within Australia. For the first time to my knowledge, this allows the exploration of the population history of the entire continent over the last 50 ka. Attempts have been made in the past to undertake similar analysis on local and regional scales, but to date such studies have been based on small spatial and/or temporal datasets [13,42,43]. The time-series analysis presented here provides a new perspective on prehistoric population changes across Australia over the last 50 ka and it suggests greater changes to prehistoric populations in the Early and Mid-Holocene than previously considered.
Taphonomic correction of the data suggests a significant number of sites by 40 ka (figure 2) and lends some support to Birdsell's rapid model of colonization, but it does not support saturation of the continent to ethnographically observed population densities during the Late Pleistocene (see figure 5 and electronic supplementary material, figure S17). Williams  demonstrated that the correction function had significant errors (up to 34%) in periods greater than 25 ka, and the failure of equation (4.3) to taphonomic data indicates that there may be unresolved issues with the technique. The high fluctuations in data to more consistent levels in values less than 32 ka, also suggest that there may be over-correction in the data until this point. Even with these limitations, the results for the Late Pleistocene (from 32 ka onwards) population levels were comparable with the Early Holocene. Using the uncorrected dataset, population estimates during this time are still low—less than 20 000 (approx. 1 person per 385 km2)—and suggest widespread settlement of the continent would probably have been unfeasible prior to the LGM.
The largest population increase began as a stepwise increase from approximately 12 ka and peaking in the Late Holocene (figures 2 and 3). Importantly, these trends can be observed in both the corrected and uncorrected data, indicating that the pattern of Late Holocene growth is not a result of taphonomic loss, a finding similarly demonstrated by Johnson & Brook . I further highlight that while the dataset is relatively patchy in some parts of the Late Pleistocene, these data are strongest over the critical period spanning the Terminal Pleistocene and Holocene, indicating these trends are likely to be reliable.
Although there is evidence of exponential growth in the Late Holocene along the lines proposed by Beaton and Lourandos—the data here suggests it forms part of a longer process of population increase beginning in the Terminal Pleistocene. The growth during this period appears to be primarily driven from southeastern Australia, and demonstrates increasing spatial use of the landscape and diversification of economic activities, with greater divergence from rockshelters to a range of other sites types (see the electronic supplementary material, figure S18). By the Early Holocene, there is evidence for increasing population—with activation of new sites ; dense occupation deposits [6,44] and rock engravings ; appearance of complex technology (e.g. fish traps) and food processing techniques (e.g. Macrozamia plants) . The initiation for this growth at 12 ka is yet to be proved, but coincides with increasing sea level  and strengthening of the northern monsoon as continental shelves were flooded [47–51].
As population growth coincides with Holocene stabilization of climate, it may be the case that high-amplitude environmental fluctuations kept population relatively low throughout the Late Pleistocene (as O'Connell & Allen  has suggested). Demographic change over the last 50 ka reveals close correlation to a range of climatic and environmental factors (see the electronic supplementary material, figure S19); it appears that climatic instability may have had significant influence over population change, with notable declines at the LGM (approx. 18 ka), Antarctic Climate Reversal (approx. 14 ka), Mid-Holocene Climatic Optimum (approx. 8–5 ka) and the onset of the El Niño Southern Oscillation (approx. 4–2 ka). The greatest impact to populations appears to have been the LGM, where declines of up to 61 per cent are observed between 21 and 18 ka. Before this, the GRAnn indicates steady increase since colonization, which was abruptly stopped at the onset of the LGM. Following the LGM, populations did not recover to pre-LGM numbers until approximately 9.6 ka.
The development and application of GRAnn to founding populations (equation (4.3)) produces coherent results. The average annual growth (0.01%) and the maximum and minimum extremes (up to 1.84%) over 50 000 years are comparable with those observed in modern ethnographic societies: Pennington  demonstrates global average growth rates of between 0.008 and 0.14 per cent over the last 90 ka, with observed extremes of up to 3 per cent. Modelling experiments similarly demonstrate Holocene growth of approximately 0.04 per cent , well within the ranges produced here. While definitive resolution on the numbers of people at colonization or at European contact cannot be provided, results here generally concur with previous studies described above. Specifically, founding populations of between 1000 and 2000 people at 50 ka (as suggested by DNA evidence) would result in a population at time of contact in the order of 770 000–1.1 million—values considered realistic by Mulvaney & White  and Butlin . The highest population estimates in the order of 1.2 million were observed at approximately 0.5 ka.
I thank Mike Smith and Sean Ulm for their extensive guidance on the paper and provision of unpublished data; Matthew Peros for guidance on the methods; and E. Williams, P. Thomson and D. Berrada for statistical support. I thank L. Robin, P. Hiscock, W. Steffen, C. S. M. Turney, P. Veth, P. Douglas, S. Mooney and N. Williams for suggestions and review. I thank T. Surovell, C. Johnson and an anonymous reviewer for formal review.
- Received February 22, 2013.
- Accepted March 26, 2013.
- © 2013 The Author(s) Published by the Royal Society. All rights reserved.