INTRODUCTION
The chronology of ancient cultures, such as Egypt or Mesopotamia, has been widely studied (Spence Reference Spence2000; Huber Reference Huber2011; Manning et al. Reference Manning, Griggs, Lorentzen, Barjamovic, Bronk Ramsey, Kromer and Wild2016; Manning et al. Reference Manning, Barjamovic and Lorentzen2017). However, currently, information about American archaeological cultures, i.e., Maya, Aztec, and Inca, is scarce (Culleton et al. Reference Culleton, Prufer and Kennett2012; Ogbourn Reference Ogburn2012; Huster et al. Reference Huster and Smith2015). Despite the widespread agreement that Inca culture is among the more important cultures of South America, there is little consensus about the chronology of the Incan period.
In the context of the historical chronology used to date, the Inca civilization dates from the period from the first imperial conquests of Pachacuti Inca Yupanqui to the moment of Atahualpa’s capture by Pizarro, a time span of no more than 100 years (Rowe Reference Rowe1945). According to this narrative, the imperial phase of the Inca State (called Tahuantinsuyu, Figure 1) started around 1438 CE, after Pachacuti Inca Yupanqui came to power. Thus, the period of the dynamic development of the Inca State, as postulated by the historical chronology, was brief. However, recent dating studies based on the radiocarbon method have revealed several discrepancies in the historical chronology, suggesting that the expansion of the Inca State occurred over a longer time frame than has heretofore been understood (Ogbourn Reference Ogburn2012; Marsh et al. Reference Marsh, Kidd, Ogburn and Duran2017; Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020; Burger et al. Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021). However, radiocarbon-based chronologies have presented discrepancies resulting from prior, different calibration curves. Thus, in recent decades, it has been challenging to find a consensus between the historical Incan chronology and the different, radiocarbon-based chronologies.
Figure 1 The Inca Empire modified from Ogburn (Reference Ogburn2012) and Marsh et al. (Reference Marsh, Kidd, Ogburn and Duran2017).
The historical chronology, according to Rowe (Reference Rowe1945), was based mainly on the Spanish chronicle written by Miguel Cabello Balboa (1586) after the fall of the Incan Empire. According to this chronology (Rowe Reference Rowe1945, Reference Rowe1990, Reference Rowe2003), the Cordillera Vilcabamba region was conquered by the Pachacuti Inca Yupanqui approximately 10 years after Pachacuti took power. As part of the occupation of the new territories, the Pachacuti began the large-scale urbanization of the Machu Picchu region by erecting the llaqta of Machu Picchu and dozens of new sites around it. Based on the historical sources, authors such as Chavez Ballon (Reference Chavez Ballón1971) and Salazar (Reference Salazar, Burger and Salazar2004) suggest that the construction of those sites may have begun in 1450 CE.
Recently published studies by Ziółkowski et al. (Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020) and Burger et al. (Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021) show that the sites could have been constructed about 30 years earlier than assumed in the historical chronology presented by Rowe (Reference Rowe1945). In both studies, radiocarbon-dated artifacts (bones, teeth, or charcoal) have been used to obtain the calendar date corresponding to the beginning of the occupation of the Machu Picchu (Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020; Burger et al. Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021), Choqesuysuy and Chachabamba (Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020). Thus, already at the level of the very beginning of the Incan imperial expansion in the Machu Picchu area, there are discrepancies between the historical dates and recent findings.
Similar suggestions have been made in the context of the Incan expansion northward, namely into the area of present-day Ecuador. The study presented by Ogbourn (Reference Ogburn2012) also pointed to some discrepancies between historical chronology and chronology based on Bayesian analysis of radiocarbon dates for the Incan occupation of Chamical in Ecuador. Ogbourn (Reference Ogburn2012) likewise suggests that high-precision radiocarbon dates could be used to develop local chronologies based on Bayesian modelling for historical sites associated with the expansion of the Inca Empire. Considering that the Tahuantinsuyu occupied an area from today’s southern Colombia to central Chile, which is under the influence of Northern Hemisphere (NH) and Southern Hemisphere (SH) air parcels (Ancapichún et al. Reference Ancapichún, De Pol-Holz, Christie, Santos, Collado-Fabbri, Garreaud, Lambert, Orfanoz-Cheuquelaf, Rojas, Southon, Turnbull and Creasman2021; Figure S1), some of the difficulties identified by Ogbourn (Reference Ogburn2012) are the differences in the calendar age obtained, depending on the calibration curve used. Relatively similar conclusions were reached by Marsh et al (Reference Marsh, Kidd, Ogburn and Duran2017) in an attempt to apply Bayesian modelling to the Argentinian region of Tawantinsuyu. Differences between the NH calibration curve (IntCal20; Reimer et al. Reference Reimer, Austin, Bard, Bayliss, Blackwell, Bronk Ramsey and Talamo2020) and the SH calibration curve (SHCal20; Hogg et al. Reference Hogg, Heaton, Hua, Palmer, Turney and Wacker2020) amount to about 20 years in the study period (Ogbourn Reference Ogburn2012). As there is a lack of information on the proportions of Northern and Southern hemisphere air masses in the study area, it was difficult to decide how to combine both calibration curves.
The radiocarbon concentration differences between the NH and the SH force the creation of separate calibration curves for both hemispheres. In the youngest period, those curves were built based on radiocarbon concentrations in tree rings that were previously dated dendrochronologically (Willis et al. Reference Willis, Tauber and Münnich1960). In the newest curves, IntCal20 and SHCal20 (Reimer et al. Reference Reimer, Austin, Bard, Bayliss, Blackwell, Bronk Ramsey and Talamo2020; Hogg et al. Reference Hogg, Heaton, Hua, Palmer, Turney and Wacker2020), many periods were substituted with 1-year-resolution measurements based on annual tree rings. For reliable radiocarbon dating, the selection of the appropriate calibration curve is dictated by the location of the sampling site as that location relates to the convergence of NH and SH air parcels. In some cases, it is necessary to use a mix of both calibration curves (Marsh et al. Reference Marsh, Bruno, Fritz, Baker, Capriles and Hastorf2018), but this raises this issue concerning the appropriate proportion of both calibration curves and their associated standard deviation. Additionally, is the calibration curve 50:50 ± 20, as reported by Ziółkowski et al. (Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020), or another ratio, e.g., 40:60 ± 5? Therefore, to estimate as precise and exact an Incan time frame as possible and to find a consensus between different radiocarbon-based Inca chronologies, it is crucial to study in-depth the proportions of NH and SH air parcels that arrived at the different Incan settlements discussed here.
The main goal of the present investigation is to develop a new methodological approach to reconstruct a radiocarbon-based Incan chronology, which is based on the determination of the proportions of NH and SH air parcels arriving at three relevant Incan settlements: Machu Picchu (13ºS, 72ºW), Chamical (3ºS, 79ºW), and Tiquischullpa (16ºS, 68ºW). Through the Hybrid Single-Particle Lagrangian Integrated Trajectory model v. 4 (HYSPLIT; Stein et al. Reference Stein, Draxler, Rolph, Stunder, Cohen and Ngan2015), we identified the geographical provenance (GP) of air parcels that arrived at these settlements and subsequently calculated the air parcel contributions from each hemisphere, producing recalibrated radiocarbon dates for Machu Picchu, Chamical, and Tiquischullpa.
MATERIAL AND METHODS
Samples
Radiocarbon dates were recalibrated with the adequate mix of IntCal20 and SHCal20 according to results obtained from the HYSLIT model (Section 2.3) for five relevant sites: 1) Machu Picchu (13.16ºS, 72.53ºW, 2800 m asl) (Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020, Burger et al. Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021) in Peru; 2) Chachabamba and Choquesuysuy (Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020) in Peru; 3) Chamical in Ecuador (3.48ºS, 79.26ºW, 1380 m asl) (Ogburn Reference Ogburn2012), and 4) Tiquischullpa in Bolivia (16.87ºS, 68.63ºW, 4060 m asl) (Pärssinen and Siiriäinen Reference Pärssinen and Siiriäinen1997). The lists of samples along with sample types are presented in Table 1. We used the chronology models presented by Ziółkowski et al. (Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020), Ogburn (Reference Ogburn2012), and Burger et al. (Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021) to show how the results differ depending on the contribution of each calibration curve.
Table 1 Radiocarbon dates of samples from Machu Picchu (Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020; Burger et al. Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021), Chachabamba and Choquesuysuy (Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020), and Chamical (Ogburn Reference Ogburn2012) used in our study.
Tropical Low-Pressure Belt (TLPB)
Previous studies have used the Inter-Tropical Convergence Zone (ITCZ) to represent the area where winds from the NH and the SH converge (Hua et al. Reference Hua and Barbetti2004, Reference Hua, Barbetti and Rakowski2013). However, the ITCZ concept represents a meteorological phenomenon occurring over the oceans (see Vuille et al. Reference Vuille, Burns, Taylor, Cruz, Bird and Abbott2012; Marsh et al. Reference Marsh, Bruno, Fritz, Baker, Capriles and Hastorf2018). Over continental masses, such as tropical South America, however, the ITCZ definition is not as straightforward. A Tropical Low-Pressure Belt (TLPB) covers tropical latitudes around the globe, accompanied by the convergence of winds from the NH and the SH (Ancapichún et al. Reference Ancapichún, De Pol-Holz, Christie, Santos, Collado-Fabbri, Garreaud, Lambert, Orfanoz-Cheuquelaf, Rojas, Southon, Turnbull and Creasman2021). Thus, recent studies have used the TLPB concept as a more reliable indicator of the atmospheric circulation-boundary between hemispheres (Ancapichún et al. Reference Ancapichún, De Pol-Holz, Christie, Santos, Collado-Fabbri, Garreaud, Lambert, Orfanoz-Cheuquelaf, Rojas, Southon, Turnbull and Creasman2021; Hua et al. Reference Hua, Turnbull, Santos, Rakowski, Ancapichun, de Pol-Holz, Hammer, Lehman, Levin, Miller, Palmer and Turney2022).
To determine the latitudinal limits of the TLPB influence, we used the reanalysis dataset (Kalnay et al. Reference Kalnay, Kanamitsu, Kistler, Collins, Deaven and Gandin1996) for above-sea-level pressure data from the National Center for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR). In addition, to obtain the same spatial resolution between our air parcel fields (see below) and the above sea-level pressure, the raw data from the reanalysis NCEP/NCAR was remodeled into a matrix of 1º latitude x 1º longitude. Considering austral spring, summer, and early autumn as the growing season of plants (October, November, December, January, February, and March; ONDJFM), which fixed atmospheric 14C concentration in their tissues, we studied the spring-summer months (ONDJFM) from the years 1949 to 2019. As in Ancapichún et al. (Reference Ancapichún, De Pol-Holz, Christie, Santos, Collado-Fabbri, Garreaud, Lambert, Orfanoz-Cheuquelaf, Rojas, Southon, Turnbull and Creasman2021), we considered 10ºN and 35ºS as the respective limits of the TLPB. In this way, after the evaluation of the above-sea-level pressure pattern, we established the zonal lowest-pressure values between the latitudes 10ºN–35ºS for each spring–summer month within each year studied. To minimize the outlier-effect produced by the low spatial-resolution of the reanalysis NCEP/NCAR on the latitudinal position of the TLPB, we smoothed our data using a moving average of 10 points. The average location of the TLPB was calculated using the spring-summer months (ONDJFM) from 1949 to 2019.
Geographical Provenance (GP) of Air Parcels
We used HYSPLIT (Draxler and Hess Reference Draxler and Hess1998; Draxler and Stunder Reference Draxler and Stunder1988; Stein et al. Reference Stein, Draxler, Rolph, Stunder, Cohen and Ngan2015) to assess the GP of the air parcels that potentially transported the 14CO2 that was subsequently fixed and assimilated by the local biosphere at Machu Picchu, Chamical, and Tiquischullpa (Figure 2). The HYSPLIT model v. 4 is a complete system designed to calculate simple air-parcel trajectories, complex transport, dispersion, chemical transformation, and depositional simulations. The NCEP/NCAR reanalysis data were used as the input data of the HYSPLIT model (Kalnay et al. Reference Kalnay, Kanamitsu, Kistler, Collins, Deaven and Gandin1996). On the other hand, the backward trajectories calculated using HYSPLIT are generated as a set of longitudinal, latitudinal, and altitudinal data, based on a window of an arbitrary number of hours into the past, for a given hour of the day (Draxler and Taylor Reference Draxler and Taylor1982). Considering that ONDJFM is the main growing season of plants, which fixed atmospheric 14CO2 concentration in their tissues, we modeled the spring-summer months from the years 1949 to 2019.
Figure 2 HYSPLIT backward trajectory analysis results for Machu Picchu (left panels; green dots), Chamical (middle panels; blue dots), and Tiquischullpa (right panels; yellow dots). Air parcel density arriving to Incan settlements. The colorbar scale indicates the total number of air parcels located on each pixel grid (1º x 1º) in: hour -120, hour -240, hour -360, hour -480, and hour -600 (top to bottom panels) during the study period (ONDJFM, 1949 to 2019). The solid white lines show the mean temporal shape of the TLPB.
Ancapichún et al. (Reference Ancapichún, De Pol-Holz, Christie, Santos, Collado-Fabbri, Garreaud, Lambert, Orfanoz-Cheuquelaf, Rojas, Southon, Turnbull and Creasman2021) demonstrated that modeling air parcel backward trajectory for 1 hour per day represents the GP of air parcels that arrived in tropical South America at this location and for this day (a 24-hour time period). Thus, to evaluate the contribution of air parcels of each hemisphere, this contribution was calculated using the modeling of one air parcel backward trajectory per day at midday (Figure 3).
Figure 3 Air parcel contribution (%) of each geographical provenance arriving in Machu Picchu (a), Chamical (b), and Tiquischullpa (c), during the studied period: Southern Hemisphere (white area) and Northern hemisphere (grey area) (Hour -600, ONDJFM from 1949 to 2019).
Each backward trajectory (12,922 for each Inca settlement) describes the spatial displacement of air parcels during a determinate time-window. To evaluate the GP variability relative to this time-window, we chose five specific hours into the past: –120, –240, –360, –480, and –600. To determine the number of air parcels that derived from the NH and SH, we counted the number of air parcels located northward and southward of the TLPB in a specific hour (–120, –240, –360, –480, and –600) by each year. Thus, we obtained 10 annual series (2 by each selected hour: –120, –240, –360, –480, and –600) of the number of air parcels arriving at a determined Inca settlement coming from the NH (5 annual series) and the SH (5 annual series). The atmospheric circulation-boundary between hemispheres is defined by the averaged position of the TLPB, calculated from each month and year studied.
To calculate the percentage of air parcels derived from the NH and SH, respectively, we averaged the annual series of air parcel GP and multiplied the result by 1.82 (182/100). For each studied location, 182 air parcels were modeled by year (one backward trajectory per day). Thus, we obtained the average percentage and standard deviation of the GP for each time of day: –120, –240, –360, –480, and –600 (Figure 3, Table 2). This procedure was done for each Incan settlement studied.
Table 2 NH and SH contribution of air parcel that arrived at Machu Picchu, Chamical, and Tiquischullpa.
Radiocarbon Calibration Curve Mixing
For the calibration curves, mixing probability distribution functions describing results obtained for percentages of NH and SH were checked. The normality of the results was first checked by the Shapiro-Wilk method. If the results were not normally distributed, their uniformity was checked by the Kolmogorov-Smirnov test. Neither natural nor uniform distribution described one of the results; therefore, the histogram was plotted, and an exponential distribution fitted to the histogram. This fit was then verified by a Kolmogorov-Smirnov test. All tests were performed with Scipy 1.7.0 and Python 3.10.0; Scipy stats.shapiro() and stats.kstest() were employed for testing.
Calibration and modeling of radiocarbon dating results were made with OxCal 4.4.4. r:1 (Bronk Ramsey and Lee Reference Bronk Ramsey and Lee2013). The IntCal20 radiocarbon calibration curve (Reimer et al. Reference Reimer, Austin, Bard, Bayliss, Blackwell, Bronk Ramsey and Talamo2020) was used for NH air, and the SHCal20 radiocarbon calibration (Hogg et al. Reference Hogg, Heaton, Hua, Palmer, Turney and Wacker2020) curve was used for SH air. Radiocarbon calibration curve mixing was used to obtain probability distribution functions describing the results of air parcel contribution percentages during the period 1949–2019. All OxCal scripts and Python code used are available in the Supplementary Material.
RESULTS AND DISCUSSION
Northern and Southern Hemisphere’s Air Parcel
The air parcels arriving at Machu Picchu, Chamical, and Tiquischullpa came from both north and south of the TLPB (Figure 2). For the Altiplano, a similar pattern was reported by Ancapichún et al. (Reference Ancapichún, De Pol-Holz, Christie, Santos, Collado-Fabbri, Garreaud, Lambert, Orfanoz-Cheuquelaf, Rojas, Southon, Turnbull and Creasman2021). These authors found that from the bomb period (1950–1985), atmospheric radiocarbon levels above tropical South America were affected by each of the following: the geographical provenance of air parcels, the carbon reservoirs acting on these air parcels, the altitude of the 14C records, and the bomb-14C pulses. At the global scale, during the bomb period, human-made events and carbon reservoirs had large effects on atmospheric radiocarbon levels at the following levels: fossil fuels (-7‰ year-1), oceans (-80‰ year-1), terrestrial biosphere (-44‰ year-1), nuclear weapon testing (+320‰ year-1), and cosmic production of 14C (+8‰ year-1) (Heimann and Maier-Reimer Reference Heimann and Maier-Reimer1996; Masarik and Beer Reference Masarik and Beer1999; Key et al. Reference Key, Kozyr, Sabine, Lee, Wanninkhof and Bulliste2004; Marland et al. Reference Marland, Boden and Andres2006; Naegler and Levin Reference Naegler and Levin2009; Levin et al. Reference Levin, Naegler, Kromer, Francey, Gomez-Pelaez and Steele2010). In contrast, during the preindustrial period (including the Tahuantinsuyu period), carbon reservoirs had a small effect on atmospheric radiocarbon levels because of the absence of human-made events: Oceans (mostly the Southern Ocean; -9‰ year-1), terrestrial biosphere (-0.4‰ year-1), and cosmic production of 14C (+9‰ year-1) (Heimann and Maies-Reimer Reference Heimann and Maier-Reimer1996; Masarik and Beer Reference Masarik and Beer1999; Key et al. Reference Key, Kozyr, Sabine, Lee, Wanninkhof and Bulliste2004; Naegler and Levin Reference Naegler and Levin2009). The cosmic production of 14C mainly occurs in the lower stratosphere, with the subsequent entry of 14C entry to the troposphere from mid and high latitudes during early spring; this 14C source occurs as a nearly symmetrical mechanism between both hemispheres (Masarik and Beer Reference Masarik and Beer1999; Levin et al. Reference Levin, Hammer, Kromer, Preunkert, Weller and Worthy2021). In fact, several studies have found that the main driver of the interhemispheric atmospheric 14C gradient from the preindustrial period (including the Tahuantinsuyu period) was the ocean-atmosphere 14CO2 exchange occurring in the Southern Ocean (Braziunas et al. Reference Braziunas, Fung and Stuiver1995; Rodgers et al. Reference Rodgers, Mikaloff-Fletcher, Bianchi, Beautien, Galbraith, Gnanadesikan, Hogg, Iudicone, Lintner, Naegler, Reimer, Sarmiento and Slatter2011). The 14C initially formed in the atmosphere was taken up by the ocean and transported through its interior via thermohaline circulation for centuries to millennia (Weaver et al. Reference Weaver, Marotzke, Cummins and Sarachik1993; Braziunas et al. Reference Braziunas, Fung and Stuiver1995). As a result of the natural radioactive decay of the 14C, the Southern Ocean deep water masses have a low 14C:12C isotopic ratio relative to the atmosphere (Key et al. Reference Key, Kozyr, Sabine, Lee, Wanninkhof and Bulliste2004; Graven et al. Reference Graven, Guilderson and Keeling2012). Thereby, the upwelling of 14C-depleted deep water to the Southern Ocean surface, along with the subsequent ocean-atmosphere CO2 exchange, dilutes the Δ14CO2 of the SH’s atmosphere. This mechanism causes the SH atmosphere to have lower atmospheric 14C levels than the NH atmosphere (Braziunas et al. Reference Braziunas, Fung and Stuiver1995; Rodgers et al. Reference Rodgers, Mikaloff-Fletcher, Bianchi, Beautien, Galbraith, Gnanadesikan, Hogg, Iudicone, Lintner, Naegler, Reimer, Sarmiento and Slatter2011). Based on the above, we assume that the GP (NH or SH) of air parcels arriving at Machu Picchu, Chamical, and Tiquischullpa was the main factor influencing the atmospheric 14C concentrations above these locations during the preindustrial period.
At –120 and –240 hours, we found that air parcels coming from north of the TLPB were concentrated mostly in the Amazon basin (Figure 2). Then, at –360, –480, and –600 hours, our analysis evidenced that air parcels came from tropical latitudes of the NH. At the same time, at –120 and –240 hours, we found that air parcels that had come from south of the TLPB were concentrated in western South America and the Eastern South Pacific Ocean (Figure 2). At –360, –480, and –600 hours, air parcels that came from south of the TLPB were dispersed in subtropical and Southern latitudes of the SH. Ancapichún et al. (Reference Ancapichún, De Pol-Holz, Christie, Santos, Collado-Fabbri, Garreaud, Lambert, Orfanoz-Cheuquelaf, Rojas, Southon, Turnbull and Creasman2021) found similar patterns for the Altiplano location. These results reveal the displacement, over time, of air parcels and prove that air parcels arriving at Incan settlements had come from both the Northern and Southern hemispheres.
To estimate the air parcel contribution from NH and SH for Machu Picchu, Chamical, and Tiquischullpa, and to evaluate the time variability of these contributions, we calculated the average and the standard deviation of the annual series of the number of air parcels that came from the NH and the SH for –120, –240, –360, –480, and –600 hours (Table 2). For all study cases, at the –120 hours, air parcels mostly came from SH. However, over the day, the backward trajectory of air parcels revealed that a large number of them came from the NH.
At both –480 and –600 hours, air parcel provenance seemed to stabilize, with a similar contribution between these. However, systematic subtractions between late and early hours of the average of air parcel geographical provenance percentages showed a minimum value when subtracted from the average of the geographical provenance percentage of the –600 hour and the –480 hour. This result suggests that at the –600 hour, the air parcel contribution of each hemisphere was definitive and has likely not changed over time. Moreover, if we used the –480 hour, the uncertainties (standard deviation) could be underestimated, causing an overestimation of the precision of the radiocarbon dating. Therefore, the average and the standard deviation of the annual series of the number of air parcels that came from NH and SH, calculated for the –600 hour, allow to estimate of the Northern and Southern Hemisphere’s average daily air parcel contribution. Finally, to test the consistency of geographical provenance percentages of air parcels estimated from the –600 hour, we used different growing season periods (for example, November to February) and found no large differences (Table S1). Thus, to recalibrate the radiocarbon dates from Machu Picchu, Chamical, and Tiquischullpa, we used NH and SH contributions of air parcel obtained from the –600 hour (Figure 3).
Our time series of air parcel contributions (for the –600 hour; Figure 3) revealed high interannual variability and decadal oscillations. These patterns are represented by the average and the standard deviation in Table 2. Despite these patterns, recent climate change has caused extreme events around the globe (IPCC 2019), and yet, we found no marked trends in our data that could be attributable to recent climatic changes. This suggests that the GP contribution of air parcels calculated for the –600 hour could represent the GP contribution of air parcels that occurred during the Tahuantinsuyu period.
Thompson et al. (Reference Thompson, Mosley-Thompson, Davis, Zagorodnov, Howat, Mikhalenko and Lin2013) demonstrated that oxygen isotopic ratios (δ18O) recorded in ice cores from low latitudes can provide a wealth of unique information about the past climate in the tropics. Considering the above, to compare the recent climate state with the Tahuantinsuyu climate state, we used annually resolved ice core records from the Quelccaya ice cap that extend back ∼1800 years and provide a high-resolution δ18O record of the climate variability from tropical South America (Figure S1). We found no statistical differences between δ18O values from recent decades and the Tahuantinsuyu period, suggesting similar climate states (Figure S2). Therefore, we assume that the GP contribution of air parcels calculated for the –600 hour represents the GP contribution of air parcels that occurred during the Tahuantinsuyu period.
Recalibration of the Radiocarbon Dates from Inca Settlements
For Machu Picchu, air parcels come from the NH and SH at a proportion of approximately 50% ± 24% (Table 2). The air parcel contributions from 1949–2019 modeled data are described by a uniform probability distribution (Figure 4a). Kolmogorov-Smirnov uniformity test statistics were equal to 0.117, and the p-value was 0.264 (n = 71). For the Chamical site, air parcels come from the NH and SH at a proportion of approximately 29%/71% ± 18% (Table 2). The 1949–2019 modeled data of air parcel contributions are described by an exponential probability distribution (Figure 4b). The Kolmogorov-Smirnov test of exponentiality for NH were equal to 0.131, and the p-value was 0.162 (n = 71). For the Tiquischullpa site, air parcels come from the NH and SH at a proportion of approximately 41%/59% ± 13% (Table 2). The 1949–2019 modeled data of air parcel contributions are described by a normal probability distribution (Figure 4c). Shapiro-Wilk normality test statistics were equal to 0.993, and the p-value was 0.961 (n = 71).
Figure 4 Histograms of NH air parcel contribution percentage (%) for Machu Picchu (a), Chamical (b), and Tiquischullpa (c).
The results of radiocarbon calibration and modeling for the Machu Picchu site are presented in Tables S2 and S3. When comparing the calibration results for Machu Picchu for all calibration approaches, one can see that the age results are shifted from a few years to a couple of decades. In general, mixed curves give ages older than IntCal20 and younger ages than SHCal20 calibration results. Mixed curves also provide broader probability ranges when comparing them to IntCal20-derived results, but sometimes, these ranges are narrower than in SHCal20-derived results. For example, Machu Picchu’s first date has been modeled to the period 1409–1484 CE (68%) using a mixed curve, whereas IntCal20 gave 1405-1470 CE as the actual range, and the SHCal20 curve modeled the time period as dating from 1403–1494 CE. The results of radiocarbon calibration for the Chamical site are presented in Table S4 and those for the Tiquischullpa site in Table S5. The pattern of the shifted dating results for two other sites is similar to that for the Machu Picchu site.
In light of the archaeological studies carried out in recent years in the Machu Picchu area (Bastante et al. Reference Bastante Abuhabda, Sieczkowska and Deza2020), the use of mixIntCal20 SHCal20 helps to explain the situation in the region. On the one hand, the functioning of Machu Picchu at this level of development is the most probable determination, especially if we consider the geospatial distribution of the other two sites of Chachabamba and Choqesuysuy (Table S2). It is likely that the establishment of Choqesuysuy took place even before the construction of Machu Picchu had begun, which indeed points to a period near the end of the 14th century. Nevertheless, it already seems unlikely that Choqesuysuy was inhabited after the Spanish Conquest period. There is also a lack of archaeological data to support such a hypothesis. As for Chachabamba, the use of mixed curves does indeed help to clarify this foundation, which most likely began only after the initial expansion of Machu Picchu, and Chachabamba was itself ritually enclosed by the Incas (most likely) before the Spanish reached the region (Bastante et al. Reference Bastante Abuhabda, Sieczkowska and Deza2020). These assumptions support the hypothesis already suggested by Ziolkowski (Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020) and Burger (Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021), according to which the historical chronology according to Rowe is outdated in the context of the urbanization of the Machu Picchu region.
The use of mixed curves in the proportions (50:50) proposed by Burger et al. (Reference Burger, Salazar, Nesbitt, Washburn and Fehren-Schmitz2021) leads to the same conclusion on the construction and occupation of the llaqta of Machu Picchu as suggested earlier in the literature (Ziółkowski et al. Reference Ziółkowski, Bastante Abuhadba, Hogg, Sieczkowska, Rakowski, Pawlyta and Manning2020). Importantly, based on the analyses described above, the currently proposed mix of the two calibration curves still finds that the mixing ratio of the curves for the Machu Picchu area is 50:50. In conclusion, given that both the earlier publications and the current text have the same proportions, the overall range of dates concerning the occupation of the region has not changed significantly. And by doing so, it is still incompatible with Rowe’s historical chronology.
For the northern areas of Tahuantinsuyu, the use of a mixed curve approach in the modelled proportions of SHCal20 and IntCal20 also seems to reflect the chronological situation in the region, at least for the Chamical site. Both the start date and the end date of the Inca presence slightly shift in the timeline towards the older ages, which is particularly interesting given the latest theories about the earlier appearance of this culture in this area. Nevertheless, it would be particularly important for the Chamical site to obtain new samples that would not have been abandoned by the Incas but would instead have been occupied at least until the European invasion.
As for the sites in the area of present-day Bolivia (Tiquischullpa), the comparison of the mixed curve with IntCal20 does not seem to provide significant differences. At the level of determinacy, possible via such dating methods, a difference of a few years is certainly difficult to grasp. One can see a slight shift on the time axis after applying the mixed curve with respect to the individual two calibration curves, where in fact the results obtained lie, in most cases, in the middle of the difference between IntCal20 and SHCal20. The further use of a mixed curve could possibly refine the individual results, especially if new samples could be obtained to perform new dating.
CONCLUSIONS
Air parcels from both NH and SH were observed to be of influence for all three Incan sites. The influence of NH air parcels was significant far to the Southern Hemisphere (between 28% and 50% for the sites studied). This is consistent with the TLPB variability and the proximity of the sites studied with the TLPB. Our investigation, which is based on an empirical method, yields three proportions for determining the mix of NH and SH 14C curves in radiocarbon dating studies above tropical South America. Therefore, this study presents a new and suitable methodological approach to date radiocarbon samples located in tropical latitudes, which must be supported by high-resolution paleoclimatic tracers such as δ18O. Finally, our study confirms that great attention should be paid when applying radiocarbon calibration to radiocarbon measurements of samples originated in regions under the influence of the TLPB.
ACKNOWLEDGMENTS
The research project to establish the chronology of the development of the Incan domination in the Vilcabamba region is funded by the National Science Centre of Poland, within the framework of the OPUS grant entitled Chronology of the Inca Expansion in Cordillera de Vilcabamba (Peru), number 2020/37/B/HS3/01622. This work was also supported by the ANID-BMBF project AVOID grant no. 180005 and FONDECYT grant no. 1201810. We thank two anonymous reviewers and Associate Editor Pieter Grootes for their valuable comments that improved the manuscript.
Supplementary material
To view supplementary material for this article, please visit https://doi.org/10.1017/RDC.2022.87
