Introduction

Human societies are intrinsically embedded in their natural environment. Throughout history, we have been relying on nature to collect food, hunt animals, procure wood and mine ores. Our need for energy and resources inescapably induces a certain impact on the environment. Today, human impact is so widespread that global climatological and environmental changes are threatening the persistence of the earth system in a Holocene-like state (Steffen et al. 2015). In the past as well, human societies often had to walk the tightrope between sustainability and overexploitation. For societies to be resilient over longer time scales, sustainable resource exploitation strategies were needed. With resource exploitation strategies, we denote all of society’s aims and practises to target and collect resources for subsistence or other purposes. For these to be sustainable, a balance between society and nature must be found in which the consumption of resources does not exceed availability and/or regeneration rate (Muys 2013).

It is often axiomatically assumed that societies primarily affect their immediate surroundings in search for the required foodstuffs and material resources. This assumption has been conceptualised under the moniker of the ‘catchment area’, denoting a given area suitable for exploitation to sustain a community (Vita-Finzi & Higgs 1970). The catchment area of a community is defined by projecting a fixed spatial extent around a settlement location, determined by walking time limits observed in contemporary rural economies or ethnographic research. The underlying assumption is that, all else being equal, the greater the distance from the site, the less likely an area is to be exploited, and the less rewarding its exploitation due to the necessary energy expenditure to move to and from the site.

The extent of the exploitable territory is only one of the parameters that defines the sustainability of a given community. The nature of the resources available within that extent is equally important. The overall quality of the landscape in providing the necessary foodstuffs and resources to sustain a given number of organisms can be quantified as its ‘carrying capacity’. Hailing originally from ecology, this term has been defined as ‘the maximum population of a given organism that a particular environment can sustain’ (Allaby 2010: 67). Despite (or perhaps because) its quantitative orientation, practical difficulties in measuring carrying capacity have led the concept to be labelled not only deficient in theory, but also unrealistic in implementation and impossible to measure (Hayden 1975: 11).

In recent years, major progress has been made in the study of past human-environment interactions, resulting in an ever-growing intertwining of social and ecological data (Schoon & van der Leeuw 2015). Yet, the study of past resource exploitation strategies has, so far, typically relied on the separate analysis of (paleo-)environmental datasets and the subsequent overlaying of archaeological data such as settlement patterns and material culture as proxies for human activity. In such a setup, different scientific fields approach a common research question by conducting separate analyses of the appropriate datasets and comparing findings post hoc (Haldon et al. 2018). While these analyses are informative and useful, they inherently approach human-environment interactions through the impacts resulting from cumulative effects over time.

To address this, we advocate for agent-based modelling (ABM) as a suitable bottom-up approach to disentangle strategies of decision-making on an individual and collective level. ABM offers a platform for interdisciplinary cooperation, capable of simulating human-environment interactions and their effects as continuous processes (Schulze et al. 2017). It provides a formal computational framework that allows a more dynamic approach to complex issues such as the delineation of catchment areas and calculations of carrying capacity. ABM has been successfully applied in archaeological studies on human-environment interactions and resource exploitation strategies for topics such as subsistence strategies (Gunaratne & Garibay 2020), risk-decreasing strategies driving social organisation (Shultz & Costopoulos 2019), and resource distribution and exploitation strategies (Janssen & Hill 2016).

In this paper, we present the first iteration of an agent-based model of resource exploitation and subsistence strategies during the Iron Age to Hellenistic period (1000–25 BCE) in the area of Sagalassos in southwest Anatolia (modern Turkey). For most of this period, hilltop sites on strategic locations constituted the focal points of the settlement pattern in inland Anatolia (Vanhaverbeke et al. 2010). Such hilltop settlements were agricultural communities dependent on their immediate surroundings for resource procurement and subsistence activities (Daems & Poblome 2016). It has been suggested that these settlements were significant drivers of local environmental change during this period (Daems et al. 2021). Their exact functioning, given location-specific opportunities, constraints and legacies from past land use, however, remains understudied. Our simulation-based approach allows us to explore the sustainability of these settlements in light of the available archaeological and environmental data. We will specifically investigate how differential availability of, access to, and demand for resources in a realistic landscape gives rise to catchment areas for known settlements. We will also explore how the sustainable operating space for these settlements is shaped through the interaction between the landscape’s productive capabilities and their proximity to one another. Through this work, we aim to advance the usage of computational modelling and simulations in Classical archaeology in general, and Anatolian archaeology in particular.

Methodology and Data

The main methodological framework used in this paper hails from computational modelling, more specifically agent-based modelling (ABM). ABM is a powerful tool for hypothesis testing and is well-suited to simulating processes of human agency, reactive environments, and the feedbacks between them (Brughmans et al. 2019). ABM also enforces conceptual clarity by requiring all assumptions in theory and data to be made explicit. It therefore offers an ideal medium to facilitate interdisciplinary collaboration (Schwarz et al. 2020). In this paper, we discuss the construction and first results of SAGAscape, a new model of subsistence and resource gathering in the territory of Sagalassos written in the agent-based modelling software NetLogo v6.1.1 (Wilensky 1999). In this section we discuss the different inputs, settings and assumptions implicit in the model. Afterwards, we give an overview of a single simulation step.

1. Initiating simulation environment

SAGAscape simulates the procurement of three different resources: food, wood and clay. These were chosen as ‘token’ resources representing the three main material requirements for any community: subsistence, fuel and artisanal production. Settlements and their associated households, which constitute the agents in the model, harvest these resources in the surrounding landscape according to their availability and accessibility. Resource abundance and accessibility depend, on the one hand, on settings inherent to the simulation environment: the initial quantities present, the landscape’s capacity for regeneration after exploitation, and, specifically for wood, the occurrence of forest fires. On the other hand, local resource availability will change due to the exploitation decisions of other settlements nearby (see the next subsection on Model Functioning). Following the gathering of resources by settlements, the landscape partially regenerates its food and wood resource pools within the same timestep.

In accordance, a GIS-environment imbedding known settlements from the study period in a realistic landscape was constructed. All GIS data was prepared in QGIS v. 3.0.3-Girona (QGIS Development Team 2018). Time-dependent settlement locations and population estimates are derived from survey and excavation data compiled by the Sagalassos Project (see Table 1). Figure 1 shows these settlements in an ASTER GDEM V1 elevation raster, resampled to a resolution of 1 ha (100 × 100 m), which formed the basis for all subsequent raster calculations.

Table 1

Overview of archaeological sites included in the model, showing ID number in Figure 1, site name assigned from archaeological surveys or excavations, starting period, the type of settlement and the settlement altitude.

ID	SITE	STARTING PERIOD	TYPE	ALTITUDE (m.a.s.l.)
1	Sagalassos	Achaemenid	town	1472
2	Düzen Tepe	Achaemenid	town	1400
3	Körüstan	Iron Age	town	1270
4	Sandalion	Iron Age	town	897
5	Aykırıkça	Iron Age	town	1185
6	Hısar	Iron Age	town	1054
7	Belören	Achaemenid	town	1283
8	Çatal Pınar	Achaemenid	town	1184
9	Kayiş Kale	Iron Age	town	1399
10	Kepez Kalesi	Iron Age	town	1156
11	Kökez Kale	Iron Age	town	1840
12	Bereket	Hellenistic	town	1501
13	Hacılar	Iron Age	town	1037
14	Kozlüca	Iron Age	town	971
15	Düver Yarımada	Iron Age	town	940
16	Seydiköy	Iron Age	town	822
17	KoçaPınar	Iron Age	village	919
18	Düver Çay 1	Iron Age	hamlet	922
19	Düver Çay 3	Iron Age	hamlet	910
20	Düver SE 1	Iron Age	hamlet	949
21	South of Saga	Achaemenid	hamlet	1376
22	East of Çatal Oluk	Hellenistic	hamlet	1269
23	2004 site	Hellenistic	hamlet	1298
24	Susaklı	Hellenistic	hamlet	1302
25	Akyamaç	Hellenistic	hamlet	1154
26	F085	Achaemenid	hamlet	1159
27	Kale Mevkii	Hellenistic	hamlet	1513
28	Kapıkaya	Hellenistic	town	1430
29	Çıngıraklı	Hellenistic	village	1240
30	Ören	Hellenistic	hamlet	961
31	Taşkapı Kale	Hellenistic	village	1589
32	Düver East	Iron Age	hamlet	920

Figure 1 

Map of the study area with settlement locations. Numbers refer to entries in Table 1.

The archaeological dataset used in the simulations contains 32 sites dated to the Iron Age (1000–550 BCE), Achaemenid (550–334 BCE) and Hellenistic (334–25 BCE) periods. Accordingly, Iron Age sites are initiated at the start of a model run, with Achaemenid and Hellenistic sites entering the simulation, respectively, after 450 and 650 timesteps of one year. For most sites, the starting date is defined by specialist studies of the attested pottery material. As the goal of the model is to assess the general level of human impact on the environment, we only included sizable sites identified as towns, villages or hamlets, leaving out farms and other types of small settlements whose footprint would likely have been rather limited in space and scope. The identification of site types is defined by the extent of the artefact scatters encountered during the archaeological surveys. The surface area of a hamlet is estimated between 0.8–4 ha, for villages 4–7 ha and larger than 7 ha for towns.¹

Population estimates were allowed to vary between different model runs by sampling every settlement’s population size from random normal distributions with mean 50, 500 and 1000 for hamlets, villages and towns, respectively, and standard deviations set at 20% of the mean. These numbers are based on estimates for Düzen Tepe (Cleymans, Daems & Broothaerts in preparation) and by transposing estimates for settlement sizes in the subsequent Roman Imperial period (Boogers, Muys & Poblome accepted).

To construct resource availability from the natural environment in our simulations, we used realistic GIS data to assign properties of soil fertility, maximum wood standing stock and clay availability. For all patches situated in a lake, on a stream or at the site of a settlement, these were permanently set to zero.

The potential of the landscape for food production was defined by initial soil fertility in terms of wheat harvest per ha and per year based on the PhD research of Van Loo (2018). In this work, a model of modern-day winter wheat productivity was constructed for a part of our study area. However, the model failed to explain a large part of measured crop yield variance, did not cover the entire spatial extent of our study area, required substantial data input, and was run on modern crop data. We therefore spatially extrapolated the modelled values over the study area using altitude as the sole predictor. Afterwards, these modelled values were rescaled to more realistic values for Antiquity (see Supplementary Materials for details).

The resulting map expresses initial soil fertility, which forms the basis for calculating current patch fertility (see Supplementary Materials for the implementation). Current fertility is set to zero when a patch is used for agriculture, after which it gradually regenerates to the initial value over the course of a user-determined number of timesteps, following a logistic growth function. If a patch is harvested without proper time for regeneration, its current productivity will lie between zero and the initial soil fertility value. The model runs shown in the Results section are based on settings of two and three years of regeneration time.

To assess wood harvesting, we used the GREFOS gap succession model, which is specifically adapted to the study area, to calibrate forest productivity in units of m³ of wood per ha and per year (Kint et al. 2014). In order to reduce computation time, we opted for an exponential growth function that simply cuts off when the standing stock reaches a certain maximum value, rather than employing logistic growth. This maximum stock value, along with annual growth rate, was determined on the basis of patch altitude, exposition and climatic zone. The full calibration procedure is explained in the Supplementary Materials.

Initial forest standing stock in NetLogo is set for each patch by returning the growth function value for a random forest stand age between 200 and 400 years, representative of mature stands. To account for deforestation in the period preceding the simulation, we used a fossil pollen-based land cover reconstruction (López 2020), which indicated that initially, about 40% of the area would have been under forest cover. Forests in the model were therefore cleared in circles around the different settlements, using population estimates as weights for distribution of the total deforested area. During model runs, forest harvesting is equivalent to clearing the entire patch, after which the forest gradually grows back unless the land is converted (see below).

Clay availability was determined on the basis of soil clay content and bulk density data accessed through SoilGrids (Poggio et al. 2021). Bulk density values for the upper two metres of soil were multiplied by their respective clay content values, which results in a mass estimate of clay per ha. Only data on the upper two metres of soil were available, but this depth can probably serve as a first approximation of the maximum excavation depth in a non-industrial context.

Within this environment, all settlements are required to exploit their surroundings for the three resources just explained. Availability of these resources varies throughout space and time. Settlements have to take into account not only availability, but also accessibility of resources. In order to incorporate the highly variable topology of the area, Tobler’s (1993) hiking function was used to quantify the time needed to cross every patch (see Equation 1). We assume that inhabitants of settlements would follow the same path back and forth between a resource and their home. We therefore ignore slope directionality in Equation 1. Travelling time across lakes was rendered impossible to ensure travel by land. Travel across streams however was allowed without penalization. The shortest-time pathway from each settlement to all patches within a 5 km radius was then calculated through use of the “nw” extension of NetLogo.

$t = \frac{100}{6000}\times e^{3.5 \times \left|S + 0.05\right|}$

M1 \documentclass[10pt]{article} \usepackage{wasysym} \usepackage[substack]{amsmath} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage[mathscr]{eucal} \usepackage{mathrsfs} \usepackage{pmc} \usepackage[Euler]{upgreek} \pagestyle{empty} \oddsidemargin -1.0in \begin{document} \[t = \frac{{100}}{{6000}} \times {e^{3.5 \times \left| {S + 0.05} \right|}}\] \end{document}

With:

• – t the travelling time across the patch [hrs]
• – S the directional slope of the patch [rad]

Equation 1: Adaptation of Tobler’s (1993) hiking function, quantifying the time needed to cross a patch with a side of 100 metres.

2. Model functioning

The constructed model world described above was used to run simulations of 1000 timesteps, each representative of one year. In this section, we describe the consecutive steps constituting a single timestep.

2.1 Resource procurement

All communities attempt to fulfil their resource needs in the following order: food, wood, clay. An analogous procedure is followed for all three: settlements rank patches within a five km radius according to the ratio “actual resource availability/patch accessibility”. They then harvest these patches in descending order of the ratio until their resource requirement is met or until the settlement runs out of a predetermined amount of working days. Figure 2 illustrates the decision-making process. Each settlement runs through this entire process once, with the order of settlements decided randomly every timestep. The amount of working days available to a settlement is calculated as a user-defined percentage of the total population, multiplied by 365 days. Our model settings considered an active population percentage of 10 and 25%. The amount of working time spent on each harvesting action is described per resource in the following subsections.

Figure 2 

Flowchart of main steps in the resource procurement algorithm.

2.1.1 Food procurement

For agriculture, we used estimates derived from Roman sources to inform working days and energetic needs. It must be noted that these estimates are not necessarily directly transposable to the Anatolian highlands during the Iron Age to Hellenistic periods, but they provide a useful ballpark figure to inform the model. The maximum number of working days was limited to 250 instead of 365 to account for the length of the agricultural working season according to Columella (cited in Goodchild 2007: 301). When a patch is harvested, its current fertility value is transferred as food to the settlement and the patch is entirely cleared, as explained above. Any initial wood stock present on the patch is also transferred to the settlement to account for the full clearance of the patch. Daily per capita food requirement was set to three different levels: 0.75 kg, sufficient to fulfil the energetic needs of an average active Roman citizen (Van Limbergen 2018: 1051), 1.50 and 2 kg. The latter two settings were chosen to probe the parameter space and account for surplus production. Two additional sources of food requirement are added to these amounts: Firstly, every settlement attempts to stockpile enough food to last up to one additional year. Secondly, a grain-per-grain yield of 6/1 was presumed (Goodchild 2007: 252), meaning that one sixth of all harvested wheat is lost to resowing and needs to be replaced with additional harvesting. It was further assumed that 42 working days were spent on each ha (Goodchild 2007: 301), to which the time needed to travel to and from the field during these 42 days was added, with 10 hours accounting for one working day.

2.1.2 Wood procurement

Basic wood requirement was set at 1 kg per capita per day, to which 180 grams were added per 1000 metres of elevation above sea level (Boogers, Muys & Poblome accepted). The model was also run using 4 and 10 kg as base requirements for the sake of parameter exploration. As with food, when a wooded patch is harvested, the whole current standing stock is transferred to the settlement. The amount of working days required per harvested patch is calculated by estimating the number of headloads required to carry off all wood and attributing time to their gathering and hauling. The weight of a single headload was sampled from a random normal distribution with mean 29.21 and standard deviation 14.14 kg, while respecting a minimum value of 4.54 kg (Kefa et al. 2018: 476). Conversion of wood mass to volume assumed a mass density of 695 kg/m³, the average wood density of the species present in the study area (sd 189) (see Supplementary Materials). Based on field measurements by Preston (2012: 30–33), gathering time of a single headload was set to 4.2 minutes per kg. Hauling time was calculated as the number of headloads multiplied by the shortest travelling time to the patch and back to the settlement.

2.1.3 Clay procurement

Annual clay requirement was arbitrarily set to 1 and 10 kg per person. A threshold in clay availability was established to enforce a larger degree of selectivity for this resource. The reported model runs are based on thresholds of 0.25 and 0.50 tons per m³. When a patch is excavated for clay, 10 m³ of pure clay is gathered. The amount of working days required for the excavation, transport and firing of one such unit is based on a number of assumptions. First, quarrying and placing the clay in baskets takes 0.193 days per m³ (Delaine 1992: 182). Every 50-kg basket then needs to be carried back to the settlement, for which the shortest travelling time to the patch is charged twice, further taking a slowing factor of 6.5 into account (estimated on the basis of Delaine 1992: 182). Estimates for Roman Imperial tableware production in Sagalassos attribute an additional 4.6 working days per ton of clay (Janssen et al. 2017).

As with food sourcing, any wood present on a patch selected for clay excavation is added to the settlement’s pool, following the same procedure as described in section 2.1.2. The patch then remains permanently void of forest cover. Conversely, the firing of clay constitutes an added wood requirement: brick kilns efficiencies fall in the range of 0.8–8 MJ/kg (Vide supra). We therefore selected the mean value of 4.4 MJ/kg, equivalent to 0.3 kg of wood per kg of clay (see Supplementary Materials).

2.2 Resource consumption and regeneration

The next step in a model run constitutes the actual consumption of the resources. For this, the annual requirement is simply deducted from the stock for all three resources. Forest standing stock or agricultural fertility then regenerate according to the processes delineated in section 1. Some flexibility between the two land uses is introduced by allowing agricultural patches to reforest after an abandonment period of 6 years. In contrast, clay does not regenerate and patches used for clay procurement can no longer be used for agriculture or produce wood.

2.3 Disturbances

SAGAscape generates two types of events with potentially negative effects to the settlements: bad harvests and forest fires. A bad harvest year is instantiated when a sample drawn from the Poisson distribution, with the expected rate of occurrence set to a user-defined value, is at least one. Our model considers an arbitrary rate of occurrence of six years. A single sample is drawn every timestep and, if a bad harvest year occurs, food stocks of all settlements will be halved after the food procurement step during the next year.

The forest fire module was derived from the ILand model (Seidl, Rammer & Spies 2014) and separates the processes of ignition and spread (see Supplementary Materials for details). Every forested patch is asked to ignite if a number of conditions relating to minimum standing stock, fire return rate and average fire size are met. A number drawn from the standard uniform distribution is further used to express the severity of drought during every timestep, with higher values leading to a higher incidence of fire ignition. Fire is then allowed to spread to random neighbouring patches until a maximal fire size, dependent on a known average fire size, is reached. When a patch is affected by fire, its standing stock is set to zero.

Results

1. Landscape

Figure 3 shows the model output with regards to the impact on the landscape after 1000 steps. It is clear that all settlements have cleared the forest landscape in their immediate surroundings for agriculture and wood procurement. Forest structure is uneven-aged across the study area because of forest fires impacting destruction and regeneration across the landscape. Due to the combination of the fertility regeneration procedure and original soil fertility, actual fertility varies across the patches used for agriculture. The band pattern surrounding several settlements indicates the alternation between zones of harvested, regenerating and fully regenerated agricultural land. A large agricultural area in the valley below Sagalassos, Düzen Tepe and Körüstan is apparent. This can be explained by the close proximity of major settlements as well as smaller hamlets in the area. The intensity of survey activities in this area over the last 30 years should be considered an important explanatory factor for the observed occupation density. The geographical axis of Kökez Kale-Bereket displays a similarly extensive human impact. Again, this is not surprising given the existence of these two extensive settlements within this specific and relatively limited ecological and geographical niche as part of the broader study area.

Figure 3 

Typical model output after 1000 steps using the low settings. Orange specks are ongoing forest fires. Settlements are shown as house icons, with size proportional to population. Patches excavated for clay, only 72 in total, are not marked. Band patterns in cultivated land indicate different stages of fertility regeneration.

2. Land cover evolution

Figure 4 shows how the proportional composition of the landscape evolves over time in response to two different settings of food and wood demand (all runs with these settings are shown together). The large effect of the former and imperceptibility of the latter are apparent. A significant effect of fertility regeneration time can also be noted. The introduction of new settlements at the start of the Achaemenid (after 450 timesteps) and Hellenistic (650) periods is clearly noticeable through the upshots in the agricultural/forest ratio. Yet, the amount of forest area remains vastly larger than the amount of agricultural area, regardless of settings or timing. In the simulations without settlements (see figure S1 in Supplementary Materials), the amount of land patches under forest cover remains approximately constant at 99%, with fires keeping 1% of land cover cleared of forest. The difference in forest cover between the low and high food & wood scenarios does not exceed 5% of total land cover.

Figure 4 

Evolution of the agricultural/forest area ratio, split up between two levels of food and wood demand. Columns represent two different levels of daily food demand, 0.75 and 1.5 kg per person, while the two rows distinguish between 1 and 4 kg of daily wood consumption per capita. Notice the sudden jumps in the curve at timesteps 450 and 650.

3. Resource stocks

Figure 5 shows the evolution of the three resource stocks over time in the low and intermediate settings. Specifically, the stock/requirement ratio is shown, with the horizontal line indicating the point of sufficiency. For food, there appears to be little difference between the two scenarios as a surplus is maintained along the entire timeline and for every settlement. Different settlements also show very similar behaviour. Switching to wood however, differences between settlements become more apparent. In the low scenario, the Iron Age towns Kökez Kale, Körüstan and Kayiş Kale are able to maintain a stock level well above the sufficiency level, with the latter two even becoming more productive after the introduction of Achaemenid and Hellenistic settlements. During the Achaemenid period, Sagalassos and Düzen Tepe are also noticeably above sufficiency level, with Sagalassos outperforming Düzen Tepe. At intermediate settings, the higher wood demand (four instead of one kg per capita) pushes stock levels to the sufficiency level. Sagalassos is able to sustain a higher level of stock and is the only Achaemenid settlement to do so. Regardless of the settings, clay levels are generally kept just at the sufficiency level for most settlements, with only hamlets (such as the neighbouring Düver Cay 1 and 3, Düver SE and Düver East) consistently outperforming this general trend. In the intermediate scenario, Sagalassos struggles to gather enough clay.

Figure 5 

Evolution in resource stocks, expressed relative to demand, with each colour representing a settlement. Left: low scenario, right: intermediate scenario. IA = Iron Age, ACH = Achaemenid, HELL = Hellenistic, KöKa = Kökez Kale, K = Körustan, KaKa = Kayiş Kale, S = Sagalassos, DT = Düzen Tepe, DH = Düver hamlets. The 1-line represents sufficiency. Only in the case of clay is this level not consistently reached by certain settlements, for instance by Sagalassos in the intermediate scenario. Notice the high levels of wood collection in Kökez Kale, Körustan and Kayiş Kale in the low scenario and the high amounts of clay collected in both cases for the Düver hamlets.

Figure 6 again shows the evolution of the different stock levels over the years, this time however for the two high settings. The left side of Figure 6 shows the effects of high food and wood demand while the right illustrates the consequences of high clay demand in the face of low labour availability and high selectiveness for clay sources. At both settings, food stocks are least affected as these are sourced first, when the highest amount of working days is still available. High F&W settings still allow for sufficient food to be collected, while high C settings cause the Hellenistic settlement of Bereket to become incapable of maintaining its food stock level. Given that food and wood are sourced before clay and that the demand level for these is the same as at intermediate settings, this shortcoming must be due to the lowered active part of the population. Effects of high clay demand are more pronounced for Hellenistic settlements in general. Interestingly, Düzen Tepe outperforms Sagalassos for food in both scenarios. Many settlements struggle to gather enough wood at both high settings, with more pronounced effects in the high C scenario and for the Hellenistic period. Kayiş Kale, Körüstan and Kökez Kale, which were outperforming their peers in the low scenario, are unable to sustain their wood supply in the high C scenario and are even surpassed by all contemporary settlements. Both Sagalassos and Düzen Tepe are also incapable of sourcing enough wood at either high setting. At the high F&W settings, both settlements’ supplies sink to approximately the same level, though the rate at which this occurs is higher for Düzen Tepe. Faced with the high C scenario, Düzen Tepe is however able to surpass Sagalassos. Clay sourcing is highly disrupted by the high settings, with more obvious effects at high C settings. At high F&W settings, quite a few settlements are still able to gather enough clay, with Körüstan and Kökez Kale being notable exceptions in the Iron Age and Sagalassos and Düzen Tepe in the Achaemenid period. In the high C scenario, most settlements fail to gather any clay at all, with the group of Düver hamlets being marked exceptions in the Iron Age, as well as the isolated hamlet of Ören in the Hellenistic period.

Figure 6 

Evolution in resource stocks, expressed relative to demand, with each colour representing a settlement. Left: high F&W scenario, right: high C scenario. IA = Iron Age, ACH = Achaemenid, HELL = Hellenistic, DT = Düzen Tepe, S = Sagalassos, B = Bereket, KöKa = Kökez Kale, K = Körustan, KaKa = Kayiş Kale, DH = Düver hamlets, Ö = Ören. The 1-line represents sufficiency. All settlements except Bereket are able to gather sufficient food in both scenarios. Few settlements can however maintain adequate stocks of wood in either scenario. The Düver hamlets and Ören consistently outperform other settlements for clay.

4. Catchment area

Figures 7 and 8 show the annual average distance travelled to the location of harvest for, respectively, food and wood (see Supplementary Materials for the calculation methodology and an alternative graph with confidence intervals for wood). For the sake of comparison, we only show the results for Düzen Tepe and Sagalassos. Overall, food is always harvested closer to home compared to wood. The agricultural catchment area ranges from a 23 minute- to a 75 minute-walk, while the exploitation zone for wood mostly lies between 30 minutes and 3 hours, depending on the settings. Distance to the settlement is often correlated with the amount of food and wood gathered, though this is not an absolute relationship. Regardless of the settings, Sagalassos requires a slightly larger catchment area for food than Düzen Tepe, despite their likeness in elevation (1472 and 1400, respectively) and size. There are no common trends over time for wood, with Düzen Tepe generally decreasing its catchment area over time and Sagalassos either showing increasing, decreasing or constant behaviour, depending on the scenario. Comparisons between wood catchment areas are not very straightforward, with the outcome dependent on the settings and the point in time. In the high F&W scenario, catchment area is consistently higher for Sagalassos while the opposite is true in the low scenario; the differences however remain small. The first few timesteps are always associated with large amounts of wood gathering due to agricultural fields being cleared.

Figure 7 

Average walking times to the agricultural plots for Düzen Tepe and Sagalassos over time. Model settings are indicated on the right. The red lines are fitted linear regression curves for the average walking time, using timesteps as the sole predictor. Walking time increases over time in all scenarios, and walking times are slightly higher for Sagalassos than for Düzen Tepe. Sagalassos generally harvests greater amounts of food, except in the high F&W scenario. The highest walking times overall occur in the high F&W scenario.

Figure 8 

Average walking times to the exploited forest stands for Düzen Tepe and Sagalassos over time. Model settings are indicated on the right. The red lines are fitted linear regression curves for the average walking time, using timesteps as the sole predictor. Walking time is similar across scenarios for Düzen Tepe, and slightly decreases over time. For Sagalassos, the overall value and trend over time depend on the scenario, with the highest values occurring at high F&W settings. Sagalassos generally harvests greater amounts of wood than Düzen Tepe.

Discussion

The current version of SAGAscape emerged as an approach to studying catchment area within a specific spatiotemporal context. As with all simulation approaches, we made clear choices to focus on the core processes of past resource exploitation strategies by simplifying real-world dynamics. While the model thematically overlaps with other ABMs of resource exploitation and subsistence in the past, we consciously opted to start from a simplified basis, disregarding household composition, demographic processes or migration. Instead we rely on known archaeological sites and population estimates. Only a simple scalar was used to arbitrarily set the working proportion of the population for each model run. This stands in contrast with several other ABMs where randomly generated settlements with alterable, migrating populations engage in agriculture and/or material sourcing. For examples, see MayaSim (Heckbert 2013, Heckbert et al. 2016), Indus Village (Angourakis et al. 2020, 2022) and the family of models produced by the Village Ecodynamics Project (VEP) (Kohler et al. 2012, Crabtree et al. 2017).

Furthermore, little consideration is given to climate, except through bad harvest events and increased occurrence of forest fire. Conversely, VEP models rely on high-resolution tree ring data for precipitation and temperature reconstructions and ensuing land productivity (Crabtree et al. 2017). Indus Village is based on a number of interconnected submodels determining weather, terrain characteristics, soil water availability and crop productivity (Angourakis et al. 2022). The MedLand modelling laboratory (MML) further integrates land use and erosion/deposition processes by coupling models across programming languages (Barton et al. 2016). Finally, direct interactions between settlements are absent from the current version while trade and/or competition over resources are central to most of the aforementioned models. Settlements in SAGAscape do compete over available land, with attribution happening on a simple first come, first served basis without a tenure mechanism between timesteps.

Our approach is thus significantly less elaborate than these pre-existing ABMs. This is a direct result of the simplicity of our research question and the limited data at our disposal. Most of the existing models are data-intensive and require substantial computational effort. By limiting the number of exploitable resources to three, working on the basis of known settlements and by keeping exploitation and regeneration procedures relatively simple, we were able to curb computation time and memory load. Moreover, the effects of particular combinations of settings can still be traced, even though non-linear mechanisms are clearly present in our results. This tracing of cause and effect is especially interesting in the comparison of different settlements, such as Sagalassos and Düzen Tepe. While further developments to SAGAscape are planned (see section “Planned model development”), the current version allows to describe a certain base scenario of self-directed behaviour that can serve as a future point of reference. We highlight the most noteworthy features in the rest of this section.

The simulated landscape generally appears to produce results consistent with empirical observations. It consists of anthropogenic zones characterised by various degrees of agricultural and wood exploitation, embedded in a forest matrix. The extent of woodland exploitation seems not to depend greatly on wood demand, at least not at the levels used in our model runs. This is in agreement with the modern observation that conversion for agriculture is globally responsible for significantly more tree cover loss than forestry (Curtis et al. 2018).

Some caveats need to be stated regarding the resource abundances constituting the model environment. The methodologies employed for resource estimation are internally consistent, but require extensive processing, meaning that we are sometimes forced to employ heuristics in the absence of adequate data. For example, the procedure utilised for attributing forest growth parameters already entails a large degree of simplification. Furthermore, our implementation of the fire module only considers standing stock for potential fuel while the original took additional pools of woody biomass (e.g. dead branches) into account. Compounded with the absence of GREFOS results for the specific climate of our study period, the accuracy of our forest productivity data is reduced. The same could be said about the extrapolated and adapted agricultural productivity data, despite the underlying real-world measurements. Even so, the average agricultural area per capita over all model runs and timesteps amounts to 0.83 ha (sd 0.30), which is roughly in agreement with the 0.7 ha value proposed for Greek Sicily by De Angelis (2000). In summary, the environment generated from these base data layers will be off by a constant degree of magnitude in all cases, meaning that, at the very least, qualitative results are still likely to hold relevance as a useful basis for discussion.

Our approach to simulating human agency in the landscape also has to deal with quantitative lacunae. Specifically for the usage of clay in artisanal production, the time cost of the different steps in processing remains only a rough approximation. Better knowledge of the chaîne opératoire of artisanal production, including reconstructed figures on kiln size, clay throwing times, firing duration, personnel requirements etc. is required for improvement in this regard. Clay was also an important resource for vernacular architecture, yet clay requirements for construction activities have not been incorporated in the model. Fixing the active part of the population on a single number for all three resource-gathering activities introduces a further problem: members of a household probably contributed in different capacities to different economic activities. For instance, wood gathering is often a task for women and children, with higher proportions of men joining in during the agricultural off-season (Cooke, Köhlin & Hyde 2008). The ancient author Hesiod in fact specifically advises the month of September to be given over to wood-cutting (Osborne 1987: 38). Again, these simplifications are common to all simulation runs and, while they can definitely be improved, should not disproportionately affect any single settlement within a scenario.

Our results allow the validation of one of the major assumptions underlying the start of every model run: the maximum size of every settlement’s catchment area. In the model, the radius of this area is limited to 5 km, corresponding to the distance walkable in approximately one hour on flat terrain. According to the theoretical model of Bintliff (1999), all main subsistence activities would take place within this territorial range, with agriculture taking place in the centre. Due to the local topology however, walking times within this radius reach up to five hours in our model. As can be deduced from Figures 7 and 8, simulations for Sagalassos and Düzen Tepe indeed produce significantly larger catchments, especially for wood gathering. At both high settings, these towns even fail to source enough wood within the provisioned ranges, even though this is mainly due to labour shortage in the high C scenario. This is however not the case for all settlements. The solitary town of Kepez Kalesi is at the centre of an agricultural range that stretches up to only 35 minutes of walking distance. This community is able to gather sufficient amounts of wood even in high F&W settings, with walking distances reaching up to 2 hours. It seems that the 5 km theoretical boundary can be taken as a solid first approximation, albeit that settlement concentrations can push the boundaries of catchment areas. Labour shortages similarly limit the sourceable area, albeit through a different pathway.

To validate the model outcomes, we must note that cross-checking our results with palynological data is not straightforward because of the existence of land use types other than forestry and agriculture in the available pollen-based reconstructions. We can however make a rough generalisation with regards to agriculture. The pollen data currently at our disposal suggest that the agricultural area during the second half of the Hellenistic period was 1.5 to 2.2 times higher than at the end of the Iron Age. By contrast, the intermediate Achaemenid period is characterised by a slight decline of agricultural area. Over the course of our model runs, the amount of agricultural area indeed roughly doubles, independent of the settings. Yet, as can be seen in Figure 4, this happens in stepwise fashion. This is a direct result of the introduction of new settlements, with concomitant simulated population increases. It has been suggested that the formation of new Iron Age and Achaemenid settlements may have been (partially) driven by fission dynamics, where the population of the parent settlement is partitioned and part of the population is relocated (Daems 2019). Whereas community fission and the associated human impact on new areas in the natural environment may increase the overall amount of agricultural area, our modelling approach does not take into account a potential reduction of agricultural exploitation in the parent community. Further investigations of settlement patterns, population estimates and fission/fusion dynamics in relation to land use changes in the study area will be needed to shed light on this issue.

One particular subset of results we can briefly highlight is the comparison between Düzen Tepe and Sagalassos. Both sites were located at the edges of the Ağlasun valley, at a distance of 1.8 km as the crow flies, and originated contemporaneously during the late Achaemenid period. Given their close proximity, it has long been unclear to what extent their catchment areas may have overlapped, resulting in competition over shared resources (Daems & Poblome 2016). The question becomes even more pertinent when we consider that Düzen Tepe was abandoned during the 2nd century BCE, more or less at the same time when Sagalassos transformed into a political and urban community or polis (Daems 2019).

All the more surprising then, that our model suggests that under high demand scenarios, Düzen Tepe generally outperforms Sagalassos and that the latter actually required a larger catchment area to sustain a similar population size. For all intents and purposes, it seems that Düzen Tepe had been dealt better cards for survival. Recent calculations suggest that both communities would have been able to sustain the endosomatic energy needs of their populations (estimated around 1000 people) with the land available to them (Cleymans, Daems & Broothaerts in preparation). It must be noted, however, that these calculations only pertain to the physical subsistence of the population and do not include the need for fuel, resources for construction and artisanal production or any other good obtained from the natural environment. Given our simulation results, it may appear counterintuitive then that it is actually Sagalassos that was able to increase its energetic footprint and become the major urban centre in the area, rather than Düzen Tepe.

Yet, this finding actually corresponds with earlier hypotheses that part of the underlying drivers behind the urbanisation of Sagalassos must be found outside of local endogenous processes. More specifically, it might be related to a stronger participation by Sagalassos in wider regional and inter-regional networks of exchange (Monsieur, Daems & Poblome 2017) and a more proactive response to the political and economic policies of the Seleucid kingdom in the 2^nd century BCE, stimulating urbanism and market formation across southwestern Anatolia (Daems & Poblome 2016; Daems 2019; Daems & Talloen 2022). While it is too early to say that this model corroborates these hypotheses based on archaeological and epigraphical data, it does provide an excellent example of the great potential of computational modelling for providing alternative perspectives and novel avenues in exploring archaeological research questions.

Planned model development Modelling is inherently an iterative affair. Starting from the fundamentals, we gradually refine and add functionality to explore more intricate aspects of the system under study, all the while taking care to avoid sinking in the swamp of undue detail in pursuit of evermore realism. Still, a number of refinements and additions can be made to improve the performance of our model.

It is clear that many of our estimates could be deemed “rough” and have room for improvement. On the side of productivity estimates, input data generated from both GREFOS and AquaCrop should provide better coverage of the study area and be temporally explicit. For GREFOS, this can be accomplished by repeatedly running the model for more classes of geographic variables, combined with the use of high-resolution palaeoclimatic data (e.g. Karger et al. 2021). Kint and co-authors (2014) also remarked that forest productivity in GREFOS is probably overestimated at higher altitudes; the first step would therefore need to consist of a slight recalibration. For fertility data, an added temporal dimension should also be taken into account. Application of the AquaCrop model to temporally explicit sedimentation data from the Gravgaz catchment shows that fertility in the study area was not constant over time, with the amplitude of change dependent on topographical location (Van Loo & Verstraeten 2021).

An animal husbandry component is also conspicuously absent in the current implementation. Yet, cattle, ovicaprines and swine must have constituted a significant part of societal metabolism in our periods of interest (De Cupere et al. 2017). Besides the immediate contribution to the diet, the rearing of animals impacts on agricultural fertility by means of manure (see Sallares 1991: 381–383 for a discussion on Classical Attica), on forest productivity through the effects of browsing (Hughes 2011) or foraging (Frémondeau et al. 2017) and on wood requirements due to the availability of dung as an alternative fuel. The addition of grazing land to our simulated landscape would render our simulations more complex, but also more in line with known palynological and archaeofaunal datasets.

In the absence of more exact estimates of human-related variables, such as resource requirements or initial maximum catchment area, it would be interesting to run our model in a full factorial design, where all combinations of variable settings are evaluated a certain number of times. Gradual differences between these settings, instead of the relatively large differences in e.g. food, wood and clay demand understood above, would allow for the delineation of thresholds and non-linear effects. For instance, our results indicate that somewhere between 10 and 25% of the population, a threshold exists at which wood collection becomes impossible for some settlements. Alternatively, starting off from a larger maximum catchment area could serve to make more desirable, if distant, resources attainable. Especially for clay, certain areas seem to have been preferable in the past (Poblome 2004: 498). Finally, a comparison of our current results with results for virtual settlements, randomly located and populated, would allow us to further distinguish the contributions of environment and human interaction to the establishment of catchment area.

The current version of SAGAscape provides the baseline for the assessment of socio-ecological sustainability in subsistence and resource exploitation of independent communities. While highly informative, we must go beyond endogenous factors alone to explore the underlying drivers of community formation, urbanism and social complexity trajectories that can explain the observed human-environment interactions in our study area. Development of SAGAscape will therefore continue in the near future, with progress focused along a few major lines of inter-community interaction. Firstly, trade between settlements constitutes an important mechanism that is absent from the current iteration. At this point, a community lacking one or more of the three resources has no way to remedy this shortage beyond its own exploitation activities. Trading surpluses in one resource for another should broaden our understanding of local sustainability. Secondly, to explore the full interaction space between communities, it is essential to account for aspects of territoriality, allowing communities to lay claim on their own catchment areas and excluding others from accessing their resources. This is particularly important to untangle interactions between nearby sites such as Düzen Tepe and Sagalassos. It will be interesting to explore how such political delineations will shape resource exploitation strategies and the spatial extent of catchment areas, as well as the resultant impact on each community’s sustainability. One issue still to be resolved is how to introduce territorial claims as arising from dynamic interactions between communities rather than a static top-down imposition. Thirdly, we currently have only a limited understanding of how Sagalassos attained its primary position as the main political and economic centre in the area. To assess how likely it is for a single community to attain a primary position, we aim to introduce the possibility for communities to establish settlement hierarchies through an assessment of site importance based on a composite metric of site type, size and available resource stocks. By exploring the probability of Sagalassos being the primary centre, we can assess the likelihood of this development being the result of endogenous or exogenous factors.

Conclusions

Despite decades of archaeological and interdisciplinary studies by the Sagalassos Project, our understanding of settlement patterns, community dynamics and human-environment interactions in the area of Sagalassos in Pisidia from Iron Age to Hellenistic times remains limited due to the limited amount of excavated sites. We set out this paper with the aim of exploring the role of Iron Age hilltop sites as major drivers of environmental changes in the area, setting in motion a chain of events leading from deforestation and soil erosion to a diversification of settlement patterns across a wide range of topographic niches. In our simulations, we aimed to explore the impact of agricultural cultivation, wood procurement and clay exploitation on the natural environment and the sustainability of the immediate catchment areas of these communities.

The SAGAscape model presented here supports the hypothesis of human impact as a major driver of environmental change. Our simulations indicate a widespread deforestation of the landscape in function of agriculture and wood procurement, especially in restricted areas such as river valleys (e.g. Bereket) or more densely populated areas (e.g. the Ağlasun valley). Moreover, the results show that sufficient food, wood and clay stocks remained available for all communities, especially for the low and intermediate settings, throughout the entire runtime of the model. This overall pattern of sustainability matches our observations of continued inhabitation of hilltop sites until Hellenistic times, where the reduction of soil thickness on the hillslopes due to erosion did not cause a wholescale system collapse, but rather induced a diversification of settlement patterns towards the lower valley areas.

It is a tantalising thought whether the disruption of wood and clay stocks observed in certain communities under high resource demand settings might be likened to the impact of the increased energetic footprint from Sagalassos following its transformation to an urban community in middle Hellenistic times. At this time, Sagalassos needed to markedly extend its core catchment area, now also including the nearby Çanaklı valley, perhaps resulting in a reorientation of local resource flows that led to the abandonment of nearby Düzen Tepe. Further refinement and extensions of the model will be aimed at shedding light on this particular conundrum in the future.

It is clear that the ABM approach advocated here has been very effective in providing a new perspective onto past and ongoing research efforts to disentangle individual and collective decision-making strategies, their resultant impact on the environment and the sustainability of long-term human-environment interactions. ABM as a formal hypothesis-testing tool has slowly been finding its way into archaeology. Still, it has a long way to go before attaining its rightful place as a standard tool in Classical and Anatolian archaeology, where it remains considerably under-utilised. Yet, we strongly believe that the potential of ABM for archaeological research will one day fully disseminate through these fields and we hope that this paper will be able to contribute a small step towards this end-result.

Data Accessibility Statement

The release version of the SAGAscape model, its input and output data and associated processing scripts are freely available in the KU Leuven Research Data Repository. DOI: https://doi.org/10.48804/H3PSC6

Additional File

The additional file for this article can be found as follows: