Raimbault Juste
email: juste.raimbault@polytechnique.edu

Indirect Evidence of Network Effects in a System of Cities

Keywords: Urban Systems, Urban Growth, Spatial Interactions, Network Effects, Empirical AIC

We describe a simple spatial model of urban growth for systems of cities at the macroscopic scale, which combines direct interaction between cities and an indirect effect of physical network flows as population growth drivers. The model is parametrized on population data for the French system of cities between 1831 and 1999, which strong non-stationarity in correlation patterns suggest to apply the model on local time windows. The corresponding calibration of the model using genetic algorithms provide the evolution of interaction processes and network effects in time. Furthermore, the fit improvement when adding network module appears effective when controlling for additional parameters, what confirms the ability of the model to unveil network effects in the system of cities.

Introduction

Cities are paradoxically unsustainable and source of negative externalities, but also the best chance to reach sustainability and resilience to climate change [START_REF] Glaeser | Triumph of the city: How our greatest invention makes us richer, smarter, greener, healthier, and happier[END_REF]. The dynamics of urban systems at a macroscopic scale, and more precisely drivers of urban growth, are crucial processes that need to be understood to meet these potentialities. A better knowledge of how cities differentiate, interact and grow is thus a relevant topic both from a theoretical perspective and for policy applications. Indeed, [START_REF] Pumain | Innovation cycles and urban dynamics[END_REF] suggest that cities are incubators of social change, their fate being closely linked to the one of societies.

Various disciplines have studied models of urban growth with different objectives and taking diverse aspects into account. For example, economics are still reluctant to include spatial interactions in models [START_REF] Krugman | Space: the final frontier[END_REF] but are extremely detailed on market processes, even for models in economic geography. On the other hand, geography focuses more on territorial specificities and interactions in space but will have more difficulties to produce general conclusions. The example of this two disciplines and their misunderstandings [START_REF] Marchionni | Geographical economics versus economic geography: towards a clarification of the dispute[END_REF] shows how it is difficult to make bridges, how it needs exceptional efforts to translate from one to the other (as P. Hall did for Von Thunen work [START_REF] Taylor | A polymath in city studies[END_REF]), and therefore how it is far from evident to grasp the complexity of urban systems in an integrated way.

The simplest model to explain urban growth, the Gibrat model, assumes random growth rates for cities. It has been shown by [START_REF] Gabaix | Zipf's law for cities: an explanation[END_REF] that this model asymptotically produce the expected rank-size law (Zipf's law) for system of cities, which is considered as one of the most regular stylized facts, at least in its generalized scaling law formulation [START_REF] Nitsch | Zipf zipped[END_REF]. Explaining urban scaling laws is closely related to the understanding of urban growth, as [START_REF] Bettencourt | Why are large cities faster? universal scaling and self-similarity in urban organization and dynamics[END_REF] suggest that these reflect underlying universal processes and that all cities are scaled version of each other. This approach however does not reflect the complex relation between economic agents for which [START_REF] Storper | Rethinking human capital, creativity and urban growth[END_REF] advocate. Using a bottom-up reconstruction of urban areas with dynamical microscopic population data, [START_REF] Rozenfeld | Laws of population growth[END_REF] show indeed that positive deviations to the rank-size law systematically exist, and that these must be an effect of spatial interactions between urban areas.

Complexity approaches are good candidates to integrate these interactions into models. [START_REF] Andersson | A complex network approach to urban growth[END_REF] introduce for example a model of urban economy as a growing complex network of relations. The Evolutive Urban Theory, introduced by [START_REF] Pumain | Pour une théorie évolutive des villes[END_REF], focuses on cities as co-evolving entities and produces explanations for growth at the scale of the system of cities. [START_REF] Pumain | An evolutionary theory for interpreting urban scaling laws[END_REF] show that scaling laws could be a consequence of functional differentiations and of the diffusion of innovation between cities. The positioning regarding universality of laws is more moderate than theories of scaling [START_REF] West | Scale: The Universal Laws of Growth, Innovation, Sustainability, and the Pace of Lifein Organisms, Cities, Economies, and Companies[END_REF], as Pumain (2012b) highlights that ergodicity can difficultly be assumed in the frame of complex territorial systems. One crucial feature of this paradigm is the importance of interactions between agents, generally the cities, to produce the emergent patterns at the scale of the system. [START_REF] Pumain | Theoretical principles in interurban simulation models: a comparison[END_REF] investigate the advantages of agent-based models in comparison to more classical equation systems, and this methodological aspect is in accordance with the theoretical positioning, as it allows to take into account the heterogeneity of possible interactions, the geographical particularities, and to naturally include emergence between levels and render multi-scale patterns.

In this paper we aim at exploring further the assumption, central to Pumain's Evolutive Urban Theory, that spatial interactions between cities are significant drivers of their growth. More precisely, we consider both abstract interactions and flow interactions mediated through the physical networks, mainly the transportation networks. We extend existing models accordingly. Our contribution is twofold: (i) we show that very basic interaction models based on population only can be fitted to empirical data and that adjusted parameter values can be directly interpreted; and (ii) we introduce a novel methodology to quantify overfitting in models of simulation, as an extension of Information Criteria for statistical models, which application to our calibrated models confirms that fit improvement is not only due to additional parameters, but that the extended model effectively capture more information on system processes. This unveils network effects in an indirect way. We first review modeling approaches to urban growth based on spatial interactions.

Urban growth and spatial interactions

First of all, we must precise that we consider only models at the macroscopic scale, ruling out the numerous and rich approaches at the mesoscopic scale, that include for exemple cellular automata models, models of urban morphogenesis or land-use change models. Some include explicitly spatial interactions: [START_REF] Rybski | Distance-weighted city growth[END_REF] show that a gravity-based aggregation of population in a lattice produces realistic metropolitan patterns. Such kind of models operate however at an other level of detail than the one we consider. We also naturally rule out economic models that do not explicitly include spatial interactions.

Several models of urban growth at the macro scale have insisted on the role of space and spatial interactions. [START_REF] Frasco | Spatially distributed social complex networks[END_REF] introduce a spatialized model of social network growth that provides a good agreement with observed data both for city size hierarchy and the relation between population and density. [START_REF] Bretagnolle | Long-term dynamics of european towns and cities: towards a spatial model of urban growth[END_REF] propose a spatial extension of the Gibrat model. The gravity-based interaction model that [START_REF] Sanders | Système de villes et synergétique[END_REF] uses to apply concepts of synergetics to cities is also close to this idea of interdependent urban growth, contained physically in the phenomenon of migration between cities.

A more refined extension with economic cycles and innovation waves is developed by [START_REF] Favaro | Gibrat revisited: An urban growth model incorporating spatial interaction and innovation cycles[END_REF], yielding a version of the core of Simpop models (Pumain, 2012a) in terms of dynamical systems. This family of models have started with a toy-model based on economic interactions between cities as agents, that yield hierarchical patterns at the scale of the system [START_REF] Sanders | Simpop: a multiagent system for the study of urbanism[END_REF]. The Simpop2 model, still based on spatial interaction for commercial exchanges, including successive innovation waves, unveils structural differences between the European and the US urban systems [START_REF] Bretagnolle | Comparer deux types de systèmes de villes par la modélisation multi-agents. Qu'appelle t-on aujourd'hui les[END_REF]. The SimpopLocal model (Pumain and Reuillon, 2017a) simulates the emergence of initial settlement patterns. The Marius model [START_REF] Cottineau | L'évolution des villes dans l'espace post-soviétique[END_REF] couples population and economic growth with interactions between cities, allowing to accurately reproduce real trajectories on the former Soviet Union after calibration with processes multi-modeling.

Urban Growth and Transportation Networks

Under similar assumptions of previously reviewed models, the inclusion of transportation networks has been rarely pursued, contrary to the mesoscopic scale at which relations between networks and territories have been widely studied by LUTI models for example [START_REF] Chang | Models of the relationship between transport and land-use: A review[END_REF]. Transportation network growth models [START_REF] Xie | Modeling the growth of transportation networks: A comprehensive review[END_REF], prolific in economics and physics, are designed to reproduce network topology and do not focus mainly on the influence of networks on cities. [START_REF] Bigotte | Integrated modeling of urban hierarchy and transportation network planning[END_REF] study an optimization model for network design combining the effects of urban hierarchy and of transportation network hierarchy. [START_REF] Baptiste | Interactions entre le système de transport et les systèmes de villes: perspective historique pour une modélisation dynamique spatialisée[END_REF] models the dynamical interplay between network links capacity and city growth on a subset of the French city system. The SimpopNet model [START_REF] Schmitt | Modélisation de la dynamique des systèmes de peuplement: de SimpopLocal à SimpopNet[END_REF] goes a step further in modeling the co-evolution between cities and transportation networks, as it allows new network links to be created in time. These examples shows the difficulty of coupling these two aspects of urban systems in models of growth, and we will for this reason take into account network effects in a simplified way as detailed further.

The rest of this paper is organized as follows: our model is introduced and formally described in the next section; we then describe results obtained through exploration and calibration of the model on data for French cities, in particular the unveiling of network effects significantly influencing growth processes thanks to a novel methodology introduced. We finally discuss implications of these results.

Model Description

Rationale Some confusion may arise when surveying stochastic and deterministic models of urban growth. To what extent is a proposed model "complex" and is the simulation of stochasticity necessary ? Concerning the Gibrat model and most of its extensions, independence assumptions and linearity produce a totally predictable behavior, which is thus not complex in the sense of exhibiting emergence, in the sense of weak emergence [START_REF] Bedau | Downward causation and the autonomy of weak emergence[END_REF]. In particular, the full distribution of random growth models can be analytically determined at any time [START_REF] Gabaix | Zipf's law for cities: an explanation[END_REF], and in the case of studying only the first moment, a simple recurrence relation avoids to proceed to any Monte-Carlo simulation. Under these assumptions, it is natural to work with a deterministic model, as it is done for example for the Marius model [START_REF] Cottineau | L'évolution des villes dans l'espace post-soviétique[END_REF]. We will work with that hypothesis, capturing complexity through non-linearity.

We work on simple territorial systems assumed as regional city systems, in which cities are basic entities. The time scale corresponds to the characteristic scale associated to this spatial scale, i.e. around one or two centuries. Spatial interactions are captured through gravity-type interactions, this simple formulation having the advantage of being simple and of capturing the first law of Tobler, namely that interaction strength fades with distance. Other approaches introduced recently perform similarly at this scale [START_REF] Masucci | Gravity versus radiation models: On the importance of scale and heterogeneity in commuting flows[END_REF].

Model description

We consider a deterministic extension of the Gibrat model, what is equivalent to consider only expectancies in time. Let P (t) = (P i (t)) 1≤i≤n be the population of cities in time. Under Gibrat independence assumptions, we have Cov[P i (t), P j (t)] = 0. A linear extended version would write P (t + 1) = R • P (t) where R is an independent random matrix of growth rates (identity in the original case). This yields directly thanks to the independence assumption

that E P (t + 1) = E[R] • E P (t).
We generalize this linear relation to a non-linear relation that allows to be more consistent with model interpretation and more flexible. Denoting µ(t) = E P (t) , we take

µ(t + 1) = ∆t • f ( µ(t)) (1)
Note that in that case, stochastic and deterministic versions are not equivalent anymore, precisely because of the non-linearity, but we stick to a simple deterministic version for the sake of simplicity. The specification of the interdependent growth rate is thus given by

f ( µ) = r 0 • Id • µ + G ( µ) • 1 + N ( µ) (2)
where 1 is the column vector full of ones, and

G = G ij = w G • Vij <Vij
> such that the interaction potential V ij follows a gravity-type expression given by, with d ij distance between i and j (euclidian or network distance),

V ij = µ i µ j ( k µ k ) 2 γ G • exp (-d ij /d G ) (3)
The network effect term N is given by N i = w N • Wi <Wi> where the network flow potential W i reads

W i = k<l    µ k µ l j µ j 2    γ N • exp (-d kl,i /d N ) (4)
where d kl,i is the distance of city i to the shortest path between k, l computed in the geographical space, which can be through a transportation network or in an impedance field of the euclidian network. All seven model parameters are detailed below.

The first component corresponds to the pure Gibrat model, that we obtain by setting the weights w G = w N = 0. The second component captures direct interdependencies between cities, under the form of a separable gravity potential such as the one used in [START_REF] Sanders | Système de villes et synergétique[END_REF]. The third term is aimed at capturing network effects by expressing a feedback of network flow between cities k, l on the city i. Intuitively, a demographic and economic flow physically transiting through a city or in its surroundings is expected to influence its development (through intermediate stops e.g.), this effect being of course dependent on the transportation mode since a high speed line with few stops will skip most of the traversed territories. Note that we don't use exactly gravity flows in the network term, since there is no decay of interactions generating flows with distance, but a decay of the effect of the flow as a distance to the network: it is equivalent to assuming a long-range use of the network on average in time, and is this way complementary to the first gravity term.

Model Parameter Space

We give in Table 1 the description of model parameters, detailing the associated processes and parameter ranges. Both direct interaction and second order network flows effect have the same structure, namely separability between the effect of distance and the influence of population, an exponential decay parameter and a hierarchy parameter expressing the inequality of contribution depending on relative sizes of cities: the highest the exponent, the more contribution of smaller cities will be negligible regarding larger cities.

We propose to interpret the distance decay parameter the following way. Let fix an arbitrary fraction α and typical spatial ranges for a local urban system d L and for a long range urban system d R , consider a city i and two neighbors j, j with same population µ j = µ j , at distances d L and d R of i respectively. If we want to answer the question to what distance difference is equivalent an attenuation of α of the interaction potential with i, we obtain

d L -d R = -d G • ln α. Therefore, d G is exactly the proportionality coefficient answering this intuitive request.
Finally, we will consider only positive weights w G and w N , to follow empirical observations as detailed below. Numerical values for the weights will be given normalized by number of cities implied in the process, i.e. w G = w G /n and w N = w N /(n(n -1)/2).

Data

Our model is assumed as hybrid as it relies on semi-parametrization on real data. It could be possible to study it as a full toy-model, initial configuration and physical environment being constructed as synthetic data. We however aim at unveiling stylized facts on real data rather than on model behavior in itself, and setup therefore the model from the data we now describe.

Population data

We work with the Pumain-INED historical database for French cities [START_REF] Pumain | Le fichier de l'ined. Espace, populations[END_REF], which give populations of Aires Urbaines (INSEE definition) at time intervals of 5 years, from 1831 to 1999 (31 observations in time). The latest version of the database integrates Urban Areas, allowing to follow them on long time-period, according to the long time ontology for cities given by [START_REF] Bretagnolle | Villes et réseaux de transport : des interactions dans la longue durée[END_REF], that constructs a functional definition of cities as entities with boundaries evolving in time. We work on the 50 bigger cities in 1999. We furthermore isolate periods of similar length excluding wars, obtaining 9 periods of 20 years1 on which adjustment of the model will be done, taking into account non-stationarity in time.

Physical flows

The model takes into account stylized physical flows, in order to be rather independent of network effective shape. These are therefore assumed to follow the geographical shortest path according to terrain slope. This assumption avoids to obtain geographical absurdities such as cities with a difficult access having an overestimated growth rate. Using a 1km resolution Digital Elevation Model, we construct an impedance field, following Collischonn and Pilar (2000), which is given by

Z = 1 + α α 0 n0
where Z is the impedance of links the 1km grid network in which each cell is connected to its eight neighbors. α is the terrain slope computed with elevation difference between the two cells. We take fixed parameter values α 0 = 3 (corresponding to approximatively the real world value of a 5% slope) and n 0 = 3 which yield more realistic paths than smaller values2 .

Indicators of model performance

We work on an explanatory rather than an exploratory model. Therefore, indicators to evaluate model outputs are not directly linked to intrinsic properties of trajectories or obtained final states, but rather to a distance to the phenomenon we want to explain, i.e. the data. Given real population p i (t) (historical realizations of P i (t)) and simulated expected populations µ i (t) obtained with µ(t 0 ) = p(t 0 ) on a period of length T , we can evaluate two complementary aspects of model performance:

• Overall model performance, given by logarithm of the mean-square error in space and time

ε G = ln 1 T t 1 n i (p i (t) -µ i (t))
2

• Average local model performance, given by the mean-square error on logarithms

ε L = 1 T t 1 n i (ln p i (t) -ln µ i (t)) 2
Both are actually complementary, as using only ε G as it is generally done will focus only on larger cities and give poor results on medium-sized and small cities (for France only Paris will have reasonable fit as it strongly dominates other urban areas and cities). ε L allows therefore to take into account model performance in all cities simulated by the model.

Results

Stylized facts

Basic stylized facts can be extracted from such a database, as it has already been widely explored in the literature (Guérin-Pace and [START_REF] Pumain | 150 ans de croissance urbaine[END_REF]. We retrieve better fits of log-normal distributions of growth rates at all dates compared to normal fits, and also the fact that growth rates are mainly positive, on the cities we consider and when removing wars.

An interesting feature to look at in relation with our considerations on spatial interactions are correlations between time-series, and more particularly their variation as a function of distance. We consider 50 years overlapping time-windows to have enough temporal observations, finishing respectively in (1881,1906,1931,1962,1999), and estimate on each, for each couple of cities (i, j), the correlation between log-returns ρij = ρ [∆X i , ∆X j ] with a classical Pearson estimator, where

∆X i = X i (t) -X i (t -1)
and

X i (t) = ln P i (t) P i (t 0 )
This method, used mainly in econophysics [START_REF] Mantegna | Introduction to econophysics: correlations and complexity in finance[END_REF], reveals dynamical interactions without being biased by sizes. We show in Figure 1 the smoothed correlations curves as a function of distance, for each time period. First of all, the strong differences between each time period confirms the non-stationarity of growth rates over the whole time period, and justifies the use of local fit in time for the model. We can also interpret these patterns in terms of historical events for the system of city and the transportation network. The dynamic of the system begins with a flat correlation in 1881, around 0.2, that could be spurious due to simultaneous similar growth for all cities. It then stays flat but goes to zero, witnessing strong differentiations in growth patterns between 1881 and 1906. After 1931, the effect of the distance is clear with decreasing curves, starting between 0.4 and 0.5. We postulate that this evolution must be partly linked to the evolution of the transportation network: considering railway network for example [START_REF] Thévenin | Mapping the distortions in time and space: The french railway network 1830-1930[END_REF], the initial overall development may have fostered long range interactions flattening thus the correlation curves, whereas its maturation over time has conducted to the return of more classical interactions decreasing quickly with distance.

Model Exploration

Implementation Data preprocessing, result processing and models profiling are implemented in R. For performances reasons and an easier integration into the OpenMole software for model exploration [START_REF] Reuillon | Openmole, a workflow engine specifically tailored for the distributed exploration of simulation models[END_REF], a scala version was also developed. The question of a trade-off between implementation performance and interoperability is a typical issue in this kind of model, as a fully blind exploration and calibration can be misleading for further research directions or thematic interpretations. A NetLogo implementation, allowing interactive exploration and dynamical visualization, was also developed for this reason. Source code for models, cleaned raw data, simulation data, and results used in this paper are available on the open repository of the project at https://github.com/JusteRaimbault/InteractionGibrat.

We show in Figure 2 an example of model output. Cities color give city-level fit error and their size the population. Outliers can therefore easily be spotted (as Saint-Nazaire having the worst fit in the example shown) and possible regional effects identified. We illustrate in pink an example of geographical shortest path, from Rouen to Marseille, which reasonably corresponds to the actual current shortest path by highway. Top right plot shows trajectory in time for a given city, whereas the bottom right plot gives overall fit quality in time, by plotting simulated data against real data. The closest the curve is to the diagonal, the better the fit.

Behavior Patterns First model explorations, by simply sweeping fixed grids of the parameter space, already suggest the presence of network effects, in the sense that physical flow effectively have an influence on growth rates. We show in Figure 3 a configuration of such a case. At fixed gravity parameters and growth rate, we study variations of the parameters w N , d N and γ N and the corresponding response of ε G and ε L . At fixed values of γ N , we observe similar behaviors of the indicators when w N and d N change. The existence of a minimum for both as a function of d N , that becomes stronger when w N increases, shows that introducing the network feedback terms improves local and global fits compared to the gravity model alone, i.e. that the associated process have potential explanatory power for growth patterns.

Calibrating the Gravity Model

We now use the model to indirectly extract information on processes in time. Indeed, under assumption of nonstationarity, the temporal evolution of locally adjusted parameters shows the evolution of the corresponding aspect of the processes.

In a first experiment, we set w N = 0 and calibrate the model with four parameters on the 9 successive 20 years time windows. The optimization problem associated to model calibration does not present features allowing an easy solving (closed-form of a likelihood function, convexity or sparsity of the optimization problem), we must rely on alternative techniques to solve it. Brute force grid search is rapidly limited by the dimensionality curse. Classical methods [START_REF] Batty | The calibration of gravity, entropy, and related models of spatial interaction[END_REF] such as gradient descent fail because of the rather complicated shape of the optimisation landscape. Calibration using Genetic Algorithms (GA) are an efficient solution to find approximate solutions in a reasonable time. OpenMole embeds a collection of such meta-heuristics fo different purposes: [START_REF] Schmitt | Half a billion simulations: Evolutionary algorithms and distributed computing for calibrating the simpoplocal geographical model[END_REF] demonstrate the capabilities of these methods to calibrate models of simulation. In our case, it furthermore allow to do a bi-objective calibration on (ε G , ε L ).

We use the standard steady state GA provided by OpenMole (Pumain and Reuillon, 2017b), distributed on 25 islands, with a population of 200 individuals and 100 generations. We show in Figure 4 the calibration results on successive periods, by plotting final populations in the objective space. As expected, Pareto fronts that corresponds to compromises between the two opposite objectives are the rule. It means that the model cannot be accurate both globally and locally, and an intermediate solution has to be found. This may due to the fact that interaction range changes with city size (i.e. that terms in the potential are no longer separable), that we keep as a possible model development. The shape of the Pareto front are revealing the chaotic optimisation landscape, as for some periods such as 1921-1936 or 1962-1982 fronts are not regular and sparse. The change in shapes also translates different dynamical regimes across the periods: for 1881-1901, the quasi-vertical shape followed by an isolated front at high ε G values reveals a quasi-binary model behavior in the optimal regimes, in the sense that improving ε L under the limit is only possible through a qualitative jump at a high price for ε G . The values taken by d G for periods 1892-1911 and 1921-1936 show that larger cities have longer interaction ranges, as higher values of d G give better values of ε G .

We show in Figure 5 the values of fitted parameters in time, averaged over the Pareto front and for best single-objective solutions. First, the two peaks patterns for r 0 corresponds roughly to the patterns observed in average growth rates. The evolution of w G has a similar shape but lagged by 20 years: it can be interpreted as a repercussion of endogenous growth on interaction patterns in the following years, which is consistent with an interpretation of the interaction process in terms of migrations. The values of d G , with an increase until 1900 followed by a progressive decrease, is consistent with the behavior of empirical correlations commented above: the first 50 years windows in Figure 1 have a higher interaction range corresponding to flat correlation curves. Finally, the level of hierarchy γ G has regularly decreased, corresponding to an attenuation of the power of large cities that can be understood in terms of the progressive decentralization in France that has been fostered by the administration.

Unveiling Network Effects

We now turn to the calibration of the full model on successive periods, in order to interpret parameters linked to network flows and gain insight into network effects. The full calibration is done in a similar way with seven parameters being free. We plot in Figure 6 the fitted values in time for some of these parameters.

The behavior of growth rate and of the gravity weight relative to growth rate, that is similar to the gravity model only, confirms that network effects are indeed at the second order and that endogenous growth and direct interactions are main drivers. Network effects are however not negligible, as they improve the fit as shown before in model exploration, capturing therein second order processes. The evolution of d N , corresponding to the range on which network influences the territories it goes through, shows a minimum in 1921-1936 to stabilize again later, but at a value lower that past values. This could correspond to the "tunnel effect", when high-speed transportation does not stop in territories it goes through. Indeed, the evolution of railway has witnessed a high decrease for local lines at a date similar to the minimum, and later the emergence of specific high speed lines, explaining this lower final value. Hierarchy of flows have slightly decreased as for gravity, but are extremely high. This means that only flows between larger cities have a significant effect. This way, the model gives indirect information on the processes linked to network effects.

Evaluating Model Performance

We focus in this last experiment on quantifying the "performance" of the model, taking into account not only its predictive capabilities, but also its structure. More precisely, we want to tackle the issue of overfitting, which has been for long recognized in the field of machine learning for example [START_REF] Dietterich | Overfitting and undercomputing in machine learning[END_REF], but for which there is a lack of methods for models of simulation. We need to introduce a tool to confirm that the improvement in model fit is not only artificially due to additional parameters.

The Akaike Information Criterion (AIC) provides for statistical models for which a likelihood function is available the gain in information between two models [START_REF] Akaike | Information theory and an extension of the maximum likelihood principle[END_REF], correcting fit improvement for number of parameters. Similar methods include the Bayesian Information Criterion (BIC), which relies on slightly different assumptions. [START_REF] Biernacki | Assessing a mixture model for clustering with the integrated completed likelihood[END_REF] propose an integrated likelihood as a generalization of these criteria for unsupervised classification. [START_REF] Pohle | Selecting the Number of States in Hidden Markov Models -Pitfalls, Practical Challenges and Pragmatic Solutions[END_REF] show that in the case of selecting the number of states in Hidden Markov Models, real cases induces too much pitfalls for standard methods to work robustly, and suggest pragmatic selection based on their results and expert judgement. In our case, the problem is that it is not even possible to define the criteria. The method we propose is based on the intuitive idea of approaching models of simulation by statistical models and using the corresponding AIC under certain validity conditions. [START_REF] Bastani | Interpretability via Model Extraction[END_REF] use a similar trick of considering the models as black boxes and approaching them to gain insights, in their case to extract interpretable structure as decision trees.

Let (X, Y ) be the data and observations. We consider computational models as functions (X,

α k ) → M (k)
α k (X) mapping data values to a random variable. What is seen as data and parameters is somehow arbitrary but is separated in the formulation as corresponding dimensions will have different roles. We assume that the models have been fitted to data in the sense that an heuristic has been used to find an approximate optimal solution

α * k = argmin α k M (k)
α k (X) -Y , and we write

ε k = M (k)
α k (X) -Y 2 the corresponding mean-square error. For each optimized computational model, a statistical model S (k) with the same degree of freedom can be fitted on a set of realizations:

M (k) α * k (X) = S (k) (X), with an error s k = M (k) α * k (X) -S (k) (X) 2 .
We follow the intuition that if statistical models are good approximations of the simulation models compared to models discrepancy to reality, namely s k ε k , then the gain of information between the two should mostly capture the gain of information between simulation models3 . We define therefore an Empirical AIC measure between two simulation models by I M (1) , M (2) = ∆AIC S (1) , S (2) (5)

In practice we calibrate the gravity only model and the full model on the full time span, and choose two intermediate solutions giving M (1) at r 0 = 0.0133, d G = 4.02e12, w G = 1.28e -4, γ G = 3.82 with ε G = 31.2375, ε L = 302.89 and the full model M (2) at r 0 = 0.0128, d G = 8.43e14, w G = 1.230e -4, γ G = 3.81, w N = 0.60, d N = 7.47e14, γ N = 1.15 with ε G = 31.2366, ε L = 302.93. It is not clear how the empirical method is sensitive to the type of statistical model used. We use therefore several models for robustness, each time with the corresponding number of parameters (4 for the first and 7 for the second model): a polynomial model of the form a 0 + i>0 a i X i , a mixture of logarithm and polynomial as a 0 + a 1 ln X + i>1 a i X i and a generalized polynomial with real power coefficients that are optimized for model fit using a genetic algorithm a 0 + i>0 a i X αi . We fit the statistical models using successive years as different realizations.

Results are shown in Table 2. We give the value of s k /ε k and the ∆AIC. We also provide the ∆BIC, that correspond to the bayesian criterion differences on statistical models, to check the robustness regarding the information criterion used. We find a positive value for 5 criteria out of 6, confirming that information gain is indeed positive. The gain decreases when statistical model fit improves, and only the BIC for the optimized model fails to show an improvement. The assumption of negligible errors is always verified as the rate stays around 1%. This approach is of course preliminary and further work should be done for a more systematic testing and a more robust justification of the method. It suggests however that fit improvement in the model of simulation are effective, and that the model reveals therefore network effects.

Discussion

Multi-modeling

We obtain for most of periods and calibration points large values of γ N in the calibration of the full model, what suggests that a model with only major flows could be a possible alternative. An additional parameter determining the cutoff value for flows would have to be introduced, and to what extent this variation improves model performance should be tested in future developments. More generally, this raises the issue of more systematic multi-modeling approaches and of benchmarking models. [START_REF] Cottineau | A modular modelling framework for hypotheses testing in the simulation of urbanisation[END_REF] introduce a modular framework to compare concurrent hypothesis to explain urban growth. In the same spirit, an immediate extension of our work would be to test different functional specifications and benchmark more diverse network processes, such as the diffusion of innovation as done by [START_REF] Favaro | Gibrat revisited: An urban growth model incorporating spatial interaction and innovation cycles[END_REF] or the creation of social ties [START_REF] Frasco | Spatially distributed social complex networks[END_REF].

Theoretical implications

Our results support the hypothesis that physical transportation networks are necessary to explain the morphogenesis of territorial systems, in the sense that some aspects are fully contained within networks and cannot be approximated by abstract proxies. We showed indeed on a relatively simple case that the integration of physical networks into some models effectively increase their explanative power even when controlling for overfitting. This can be understood as a direction to expand Pumain's Evolutive Urban Theory [START_REF] Pumain | Pour une théorie évolutive des villes[END_REF]. This theory considers networks as carriers of interactions in systems of cities but does not put a particular emphasis on their physical aspect and the possible spatial patterns resulting from it such as bifurcations or network induced differentiations. The development of a sub-theory focusing on these aspects is an interesting direction suggested by our empirical and modeling results.

Urban System Specificity

The model has not yet been tested on other urban systems and other temporalities, and further work should investigate which conclusions we obtained here are specific to the French urban system on this period, and which are more general to any system of cities. Applying the model to other system of cities also recalls the difficulty of defining urban systems. In our case, a strong bias should arise by considering France only, as Lille must be highly influenced by Brussels for example. The extent and scale of such models is always a delicate question. We rely here on the administrative coherence and the consistence of the database, but sensitivity to system definition and extent should also be further tested.

Towards co-evolutive models

Our focus on network effects remains quite limited since (i) we do not consider an effective infrastructure but abstract flows only, and (ii) we do not take into account the possible network evolution, due to technical progresses [START_REF] Bretagnolle | Long-term dynamics of european towns and cities: towards a spatial model of urban growth[END_REF] and infrastructure growth in time. An interesting development would be first the application of our model with real network data, using effective distance matrices in time, computed e.g. with the train network used by [START_REF] Thévenin | Mapping the distortions in time and space: The french railway network 1830-1930[END_REF].

Then, allowing the network to dynamically evolve in time, as a function of flows, would yield a model of co-evolution between cities and transportation networks for a system of cities. This co-evolution has been highlighted empirically by [START_REF] Bretagnolle | Villes et réseaux de transport : des interactions dans la longue durée[END_REF]. This kind of model is not frequent, and [START_REF] Schmitt | Modélisation de la dynamique des systèmes de peuplement: de SimpopLocal à SimpopNet[END_REF] provides with the SimpopNet model one of the few examples. It is shown by [START_REF] Raimbault | Models coupling urban growth and transportation network growth: An algorithmic systematic review approach[END_REF] that disciplinary compartmentalization may be at the origin of the relative absence of such type of models in the literature. Indeed, it would imply to include heterogenous processes such as economic rules to drive network growth, that are quite far from the approach taken. It would however allow to investigate to what extent the refinement of network spatial structure and network dynamics can improve the explanation of urban system dynamics.

Conclusion

We have introduced a spatial model of growth for a system of cities at the macroscopic scale, including second order network effects among endogenous growth and direct interaction growth drivers. The model is parametrized on real data for the French city system between 1831 and 1999. The calibration of the model in time provides interpretations for the evolution of processes of interaction within the system of cities. We furthermore show that the model effectively unveils network effects by controlling for overfitting. This work paves the way for more complicated models with dynamical networks, that would capture the co-evolution between transportation network and territories.

Figure 1 .

 1 Figure 1. Time-series correlations as a function of distance. Solid line correspond to smoothed correlations, computed between each pairs of normalized log-returns of population time-series, on successive periods given by curve color.

Figure 2 .

 2 Figure 2. Example of output of the model. The graphical interface allows to explore interactively on which cities changes operate after a parameter change, what is necessary to interpret raw calibration results. The map gives adjustment errors by city (color) and their population (size). We illustrate in pink the geographical shortest path between Rouen and Marseille. The plot in top right panel gives in time the trajectory of a selected city, comparing simulated population with real population. The bottom right plot gives for each date the simulated data against real data: the closest the curve is to the diagonal, the better the fit.

Figure 3 .

 3 Figure 3. Evidence of network effects revealed by model exploration. Left plot gives ε G as a function of d N for varying r 0 /w N , at fixed gravity effect and γ N = 3. Right plot is similar for ε L

Figure 5 .

 5 Figure 5. Calibrated parameters for gravity model only. Each plot gives fitted values in time t for each parameter: in order from left to right and top to bottom, growth rate r 0 , gravity weight w G , gravity decay d G , gravity gamma γ G . Red and Green curves correspond to best points for ε G (respectively ε L ), whereas the blue curves give the average value over the Pareto front with standard deviation.

Figure 6 .

 6 Figure 6. Calibrated parameters for the full model. We plot values of r 0 , w G /r 0 , d N and γ N in time t, for single-objective optimal points (Red and Green curves) and averaged over the Pareto front (Blue).

Table 1 .

 1 Model parameters summary. We give the parameters names, notations, associated processes, possible interpretations, typical variation ranges and units.

	Parameter	Notation Process	Interpretation	Range	Unit
	Growth Rate	r0	Endogenous growth Growth rate	[0, 1]	1
	Gravity weight	wG	Direct interaction	Max average rate	[0, 1]	1
	Gravity gamma	γG	Direct interaction	Level of hierarchy	[0, +∞] 1
	Gravity decay	dG	Direct interaction	Interaction range	[0, +∞] km
	Feedback weight	wN	Flows effect	Max average rate	[0, 1]	1
	Feedback gamma γN	Flows effect	Level of hierarchy	[0, +∞] 1
	Feedback Decay	dN	Flows effect	Network effect range [0, +∞] km

Table 2 .

 2 Empirical AIC results.

	Statistical Model	M (1) Relative fit M (2) Relative fit ∆AIC ∆BIC
	Polynomial	0.01438	0.01415	19.59	3.65
	Log-polynomial	0.01565	0.01435	125.37 109.43
	Generalized Polynomial 0.01415	0.01399	11.70	-4.23

The time periods are more precisely:1831-1851, 1841-1861, 1851-1872, 1881-1901, 1891-1911, 1921-1936, 1946-1968, 1962-1982, 1975-1999.

More precisely, paths between main cities such as Paris-Lyon, Lyon-Marseille, Lyon-Bordeaux, were visually inspected for α 0 ∈ {2,

3, 4}and n 0 ∈ {2, 3}, and the parameters giving the most realistic paths were chosen.

There is to the best of our knowledge no way to establish this relation formally. According to[START_REF] Navakatikyan | A model for residence time in concurrent variable interval performance[END_REF], we haveAIC (k) = n log (s k /n) + 2K k if nis the number of data points and K k the number of parameters for model k. It follows thatI = n∆ log (s k /ε k ) + n∆ log (ε k /n) + 2∆K k .The two last terms can however not be related to the gain of information in simulation models as these have no likelihood function.

Acknowledgments

Results obtained in this paper were computed on the vo.complex-system.eu virtual organization of the European Grid Infrastructure ( http://www.egi.eu ). We thank the European Grid Infrastructure and its supporting National Grid Initiatives (France-Grilles in particular) for providing the technical support and infrastructure.


Morteza Djamali
email: morteza.djamali@imbe.fr

Alireza Askari Chaverdi

Silvia Balatti

Frédéric Guibal

Coralie Santelli

On the chronology and use of timber in the palaces and palace-like structures of the Sasanian Empire in "Persis" (SW Iran)

Keywords:

Timber in an archaeological context can be used to establish chronologies, to understand the history of architecture and to reconstruct cultural landscapes and natural vegetation in the past. In this study, we use the xylological identifications and radiocarbon dating results of five timber fragments recently discovered in three palaces or palace-like structures in Fars (SW Iran) dating back to the period of the Sasanian Empire (224-651 CE). We show that Qal'a-ye Dokhtar, a fortified palace to the north of Firuzabad, was constructed during the power transition from the Parthian to the Sasanian period. On the other hand, the so-called Palace of Ardashir I besides Firuzabad, was accomplished after the power takeover by the Sasanians and the political stabilisation of SW Iran under the reign of Ardashir I (224-240 CE) and his son Shapur I (240-270 CE). We also demonstrate that the 'Palace of Sarvistan' was mainly used right after the fall of the Sasanian Empire during the first centuries of Islamic domination over Iran. The discovery of timber in stone-dominated Sasanian architecture adds information on timber use in the Late Antique Near East. Mediterranean cypress (Cupressus sempervirens L.) was the only timber found in Sasanian palatial architecture, and its use suggests that the tree was one of the major cultivated elements in ancient 'Persis' most probably for its shade, beauty and building timber, but possibly also for its symbolic significance and sacred status to the Zoroastrians. Cypress trees may have played a major role in Persian gardens since antiquity, along with plane trees.

Introduction

The Sasanian Empire (224-651 CE) was the most powerful political and economic rival of the Roman Empire for about half a millennium. The Empire was established by Ardashir I (224-240 CE) in 224 CE, after he had defeated the last Parthian King Ardawan (Artabanus IV) in the plains of Hormozgan in Southern Iran [START_REF] Frye | The political History of Iran under the Sasanians: the rise of the Sasanian[END_REF][START_REF] Wiesehöfer | Ancient Persia: From 550 BC to 650 AD[END_REF][START_REF] Daryaee | Sasanian Persia, the Rise and Fall of an Empire[END_REF]. Similar to the Teispids and the Achaemenids (550-330 BCE) centuries before, the Sasanians also rose to power and established their empire first in 'Persis', the modern province of Fars in Southwestern Iran [START_REF] Alram | Sylloge Nummorum Sasanidarum Paris-Berlin-Wien. Band I, Ardashir I.-Shapur I[END_REF][START_REF] Weber | Pāgab, der Vater Ardašīrs I., in der historiographischen Überlieferung[END_REF]. Already in the Early Sasanian period, Persis stood out due to a large-scale urbanisation project and extensive agricultural production, which was guaranteed by the development of sophisticated irrigation systems [START_REF] Daryaee | The Persian Gulf Trade in Late Antiquity[END_REF][START_REF] Mousavi | The Sasanian Empire: an archaeological survey, c.220-AD 640[END_REF]. Over time, many cities were built in the province, these hosted the immigrant populations from the countryside and deported people [START_REF] Wiesehöfer | Ancient Persia: From 550 BC to 650 AD[END_REF]. The Sasanian 'Kings of Kings' also constructed a number of impressive palaces, fortifications, and Zoroastrian fire temples in Fars [START_REF] Huff | Architecture: III: Sasanian period[END_REF][START_REF] Mousavi | The Sasanian Empire: an archaeological survey, c.220-AD 640[END_REF]. Among the most famous of these buildings are the fortified complex of 'Qal'a-ye Dokhtar' and the 'Palace of Ardashir I' (or the so-called 'Ateshkadeh' or 'Ātaškada') near modern Firuzabad, both generally dated to Ardashir I's reign, and the 'Palace of Sarvistan' (hereafter 'Sarvestan'), near the modern town of Sarvestan, a building of uncertain function dated to the Late Sasanian-Early Islamic period (Fig. 1) [START_REF] Bier | Sarvistan: a study in Early Iranian architecture[END_REF][START_REF] Huff | Fürsten, Festungen und Feuertempel. Das sassanidische Imperium (224-651[END_REF][START_REF] Askari Chaverdi | Archaeological excavation in the so-called "Palace of Sāsān[END_REF]. Particular architectural features of these structures are the chahartaq (also chāhārtāq), a dome built on squinches above a square hall. These can be considered as the major Sasanian architectural innovation contributing into the later Middle Eastern architecture. Furthermore, the ayvan (also ayvān and īwān) also developed in the area, this is a large vaulted hall walled on three sides and open at the front [START_REF] Huff | Architecture: III: Sasanian period[END_REF][START_REF] Huff | Čahāṭāq[END_REF][START_REF] Callieri | Architecture et Représentations dans l'Iran Sassanide[END_REF]. Sasanian architecture certainly influenced the Early Islamic palatial architecture and urban design in Southern Iran, in Iraq and to some extent also in Syria; however the nature and degree of this influence is still a matter of debate [START_REF] Fontana | Art in Iran XII. Iranian pre-Islamic elements in Islamic art. Encyclopedia Iranica 11[END_REF][START_REF] Huff | Architecture: III: Sasanian period[END_REF][START_REF] Bier | The Sassanian Palaces and their influence in Early Islam[END_REF].

Although a chronology has been suggested for Ardashir's palaces, which was mainly based on the coins discovered during the excavations at the sites [START_REF] Huff | Ausgrabungen auf Qal'a-ye Dukhtar bei Firuzabad[END_REF], there has been, until now, no absolute dating available to confirm the suggested chronologies. The situation is even more complicated for the Palace of Sarvestan, for which a date encompassing the Sasanian and Early Islamic period has been attributed by archaeologists [START_REF] Bier | Sarvistan: a study in Early Iranian architecture[END_REF][START_REF] Askari Chaverdi | Archaeological excavations in the building known as the Palace of Sarvestan[END_REF].

Another problem concerning the aforementioned three buildings and Sasanian architecture in general is that, to date, there have been no archaeobotanical data available on the use of timber in the palaces and other Sasanian structures. Fortunately, wood samples have recently been discovered in the excavated materials and wall structures of some Sasanian monuments of Fars [START_REF] Askari Chaverdi | Archaeological excavations in the building known as the Palace of Sarvestan[END_REF]. These samples are of great importance in order to understand not only the building techniques of Sasanian architecture but also the arboricultural history and the use of timber in ancient Persia. Timber fragments in archaeological excavations can provide a wealth of information on the chronology of archaeological sites, the history of architecture and technology, agricultural practices, and the past distribution of plant species in the landscape (e.g. [START_REF] Rival | La charpenterie navale romaine: matériaux, méthodes, moyens[END_REF].

In this study, five timber fragments coming from the aforementioned Sasanian palaces and palace-like structures were analysed and radiocarbon-dated. Our main objectives were: i. to present a material source useful to provide an absolute chronology of the buildings and consequently also to give information on the 'history' of the construction and use of the buildings; ii. to provide information on the type of timber used in the construction of Sasanian buildings and; iii. to use the datings and wood identifications to shed new light on the history of tree cultivation in ancient Iran.

Material and methods

Five wood samples were collected from three monuments: the Palace of Ardashir I to the north of Firuzabad (28°53″ 51.44″N, 52°32″ 20.46″E, 1364 m), Qal'a-ye Dokhtar on the road from Shiraz to Firuzabad (28°55″ 14.79″N, 52°31″47.10″E, 1490 m), and the Palace of Sarvestan situated to the east of Lake Maharlu (also Maharlou) in the plain of Shiraz (29°1 1″ 44.20″N, 53°13″ 51.85″E, 1547 m). Fig. 2 displays the exact position of each wood fragment in the architectural plans of the buildings. Here are the descriptions of the studied wood fragments: Dokht-1 Wood fragment cut from a timber in the western part of southern wall of room 16 of Qal'a-ye Dokhtar (Fig. 2a).

Ard-1

Wood fragment from the archaeological debris in the northwestern corner of ayvan A of the Palace of Ardashir I (Fig. 2b). Sarv-1 Timber section from the wall located in the southwestern corner of room 9 of the Palace of Sarvestan (Fig. 2c; Fig. 4a).

Sarv-2

Wood fragment from the debris in the northeastern corner of room 1 supporting the large dome at the same place (Fig. 2c).

Sarv-3

Timber section from the wall in the northeastern corner of room 10 supporting the small semi-dome at the same place (Fig. 2c; Fig. 4b).

All five samples were AMS-radiocarbon dated in the Poznan Radiocarbon Laboratory (Table 1). For samples Sarv-1 and Sarv-3, which display almost complete series of tree-rings (Fig. 4a andb), only the outermost ring representing the date of the tree felling was subsampled for dating.

Results and discussion

Chronology of construction and use of Sasanian royal palaces

As reported in Table 1, the radiocarbon datings provide ages ranging from the 2nd to the 3rd century CE (Qal'a-ye Dokhtar and the Palace of Ardashir I) at the very beginning of the Sasanian Empire, and from the 7 th to the 9 th century CE corresponding to the centuries of domination of the Arab conquerors and the subsequent Iranian dynasties in the Islamic period (Fig. 3). These datings provide the first absolute ages from the Sasanian and early post-Sasanian palaces of Fars and help to shed some light on the history of the construction and use of the buildings in the Sasanian and post-Sasanian period.

Qal'a-ye Dokhtar and the Palace of Ardashir I

As illustrated in Fig. 3 (also reported in Table 1), the probability distribution curve of radiocarbon age for the Qal'a-ye Dokhtar sample points to a definitive absolute age older than 246 CE for this fortresspalace structure. This strongly suggests that it was constructed in the transitional period from the Parthian to the Sasanian rule over Persis. It seems that this absolute dating is in line with the historical events and the available archaeological remains. The Sasanian campaigns against the Parthians started in Persis at about 205-206 CE, when Pabag of the 'house' of Sasan, Ardashir I's father, dethroned the local ruler of the city of Istakhr [START_REF] Weber | Pāgab, der Vater Ardašīrs I., in der historiographischen Überlieferung[END_REF]. After succeeding his father, Ardashir I continued to conquer Parthian territories and finally defeated the last Parthian king on 28 th April 224 CE. In 226 CE, he was crowned in Ctesiphon as "King of Kings", marking the beginning of Sasanian rule [START_REF] Wiesehöfer | Ancient Persia: From 550 BC to 650 AD[END_REF][START_REF] Daryaee | Ardaxšīr and the Sasanian's rise to power[END_REF]. Qal"a-ye Dokhtar, a fortified palace that controls the access to the Firuzabad plain, was probably constructed in the Parthian-Sasanian transition period (2 nd to 3 rd century CE), most probably during King Ardashir I's reign, when he was fighting for supremacy in Persia [START_REF] Huff | Qal'a-ye Dokhtar[END_REF]. The fortress-palace would have been built with a defensive function before the final victory of Ardashir I over Artabanus IV in 224 CE; a fact that is also confirmed by pre-Sasanian coins found at the site [START_REF] Huff | Ausgrabungen auf Qal'a-ye Dukhtar bei Firuzabad[END_REF]. According to [START_REF] Huff | Ausgrabungen auf Qal'a-ye Dukhtar bei Firuzabad[END_REF], Qal'a-ye Dokhtar lost its importance after Ardashir's victory and the construction of the new palace (Palace of Ardashir I) in the plain of Firuzabad (Fig. 1). It only regained some military importance under King Yazdgerd III (632-651 CE) in his last efforts to organise resistance to the Arab conquest of Persia.

Absolute radiocarbon dating for the Palace of Ardashir I, near modern Firuzabad, suggests an age older than 257 CE (Fig. 3; Table 1), a time potentially corresponding to Shapur I's reign (240-270 CE). The comparison of the probability curves of samples from the two palaces 1). Plans are re-drawn after [START_REF] Hugi | Ghala-Dokhtar-Atechkade: monuments sasanides près de Firouzabad, Iran. Rapport préliminaire des projets de consolidation et restauration[END_REF] for the Palace of Ardashir I and Qal'a-ye Dokhtar and from [START_REF] Bier | Sarvistan: a study in Early Iranian architecture[END_REF] for the Palace of Sarvestan. See Materials and methods for more details. also suggests that the Palace of Ardashir I could also have been completed after Qal'a-ye Dokhtar, although, the radiocarbon dates provide no absolute certainty as to which of the two was built first (Fig. 3). These data are in accordance with the chronology suggested by historians and archaeologists [START_REF] Huff | Ausgrabungen auf Qal'a-ye Dukhtar bei Firuzabad[END_REF]. The strategic position and the fortifications of Qal'a-ye Dokhtar indicate the military use of this site and/or the necessity for protection of the royal family and the administrative centre of the young expanding kingdom of Ardashir I; while the absence of fortifications in the Palace of Ardashir in the Firuzabad plain shows that a degree of stability was already attained by the completion of the latter palace. Only after the Sasanian supremacy was definitively established, Ardashir I would have decided to build a residential palace in the well-watered plain of Firuzabad, a perfect starting place for his programme of revitalisation of land use and agriculture [START_REF] Huff | Ausgrabungen auf Qal'a-ye Dukhtar bei Firuzabad[END_REF].

Palace of Sarvestan

All radiocarbon datings for wood samples of the Palace of Sarvestan gave post-Sasanian ages (Table 1; Fig. 3). The highest age probability for Sarv-1 falls from the middle 7 th to the early 8 th centuries CE (655-724 CE). Sarv-2 sample gives the highest age probabilities for the middle 7 th to middle 8 th centuries CE (664-773 CE). Finally, sample Sarv-3 gives the highest probabilities for the middle 8 th to the end of the 9 th centuries CE (762-887 CE). The chronology of Sarvestan monuments has long been a matter of debate [START_REF] Askari Chaverdi | Archaeological excavations in the building known as the Palace of Sarvestan[END_REF]). An initial attribution of the palace to the Achaemenid period was soon rejected when in the late 19 th and early 20 th centuries several western travellers and archaeologists such as Flandin, Coste, Dieulafoy, De Morgan and Reuther attributed the ruins of Sarvestan to the Sasanian period. The Sasanian King Wahram V (420-438 CE) was traditionally identified as the builder of the palace in the 5 th century CE. Later, a Sasanian date was also put into question by several archaeologists, after a more detailed examination of the architecture of the building complex [START_REF] Bier | Sarvistan: a study in Early Iranian architecture[END_REF]. Especially the complex vault systems seem to have approached the structure of the Early Islamic architecture (see references in [START_REF] Huff | Architecture: III: Sasanian period[END_REF]. Some authors went even further by calling into question not only the date but also the function of the building complex. For instance, [START_REF] Bier | Sarvistan: a study in Early Iranian architecture[END_REF] suggested that the Palace of Sarvestan served as a fire temple during the first three centuries of the Islamic period. This is a very important observation from a historical perspective as it also suggests that the Arab conquerors arrived at a compromise with the Zoroastrian religious authorities and that Zoroastrianism would have been practiced for centuries after the establishment of the Islamic domination over Persia [START_REF] Askari Chaverdi | Archaeological excavations in the building known as the Palace of Sarvestan[END_REF]. Continuity in the settlement structures and their use has also been observed at several other sites in Iran in the transitional Sasanian-Islamic period [START_REF] Morony | The Islamic conquest of Sasanian Iran[END_REF]. New archaeological excavations that aimed at attributing a robust chronology to the Sarvestan complex were conducted at the site by an Iranian team of archaeologists directed by Askari Chaverdi in 2002. The excavations revealed that the building complex was mainly used during the Early Islamic period with the most active occupation phase dating to the 10 th century CE.

Our dating results confirm the recent interpretations of the Sarvestan monuments by the aforementioned archaeologists [START_REF] Bier | Sarvistan: a study in Early Iranian architecture[END_REF][START_REF] Askari Chaverdi | Archaeological excavations in the building known as the Palace of Sarvestan[END_REF][START_REF] Askari Chaverdi | Archaeological excavation in the so-called "Palace of Sāsān[END_REF]. The Sarvestan complex must be dated to the transitional period from the Sasanian to the Islamic domination over Persis and clearly represents an example of Sasanian architectural heritage in post-Sasanian times. Bier's identification of the building complex as a Zoroastrian fire temple used for several centuries after the invasion of the Arab Muslims may thus stand, especially given that palatial architecture was closely linked to religious architecture in Sasanian Persia [START_REF] Callieri | Architecture et Représentations dans l'Iran Sassanide[END_REF].

Archaeobotanical implications of cypress wood finds in Sasanian palaces

All examined wood fragments belonged to Mediterranean cypress Cupressus cf. sempervirens. Fig. 4a andb illustrate two sections of Sarv-1 and Sarv-3 timbers coming from rooms 9 and 1, respectively (Fig. 2c). Fig. 4a-1 and 4b-1 illustrate the transversal and radial sections of the two timbers. Although the two timbers have almost the same diameter (11.0 and 8.5 cm), their cutting ages are different. Sarv-3 was cut in a much younger age than Sarv-1, with 12 years against 36 years for Sarv-1. Tree-ring series show different patterns, evidencing probable different growth conditions. Indeed, tree-rings of Sarv-1 timber show an initial period of low radial growth rate indicating adverse site conditions perhaps due to high competition (Fig. 4a). Sarv-3 timber shows the absence of latewood for its last ring suggesting the felling of the tree during the growing season. This tree has benefited from less stressed growth conditions as is seen in its wider earlywood rings (Fig. 4b).

The discovery of cypress wood fragments in the Sasanian and post-Sasanian palatial architecture provides new evidence not only on timber use in Sasanian and post-Sasanian architecture but also on the cultivation history of this tree in the inland plateaus of the ancient Near East. Contrary to walnut, olive tree, grapevine and plane tree, whose arboricultural history has been partly reconstructed in Persis using their pollen records (Djamali et al., 2011a;[START_REF] Djamali | Olive cultivation in the heart of the Persian Achaemenid Empire: new insights to agricultural practices and environmental changes reflected in a late Holocene pollen record from Lake Parishan, SW Iran[END_REF], the history of cypress cultivation cannot be easily reconstructed because of the difficulties in distinguishing its pollen from other members of the Cupressaceae family. Indeed, Cupressus pollen is grouped under Junipertype or Cupressaceae and cannot be easily distinguished from the pollen of Juniperus and Thuja (e.g. [START_REF] Beug | Leitfaden der Pollenbestimmung für Mitteleuropa und angrenzende Gebiete[END_REF]. Juniperus is still surviving in some isolated populations in the central and southern Zagros [START_REF] Browicz | Chorology of trees and shrubs in South-West Asia and adjacent regions[END_REF] but was more widespread in prehistory [START_REF] Miller | Economy and environment of Malyan, a third millennium B.C. urban center in Southern Iran[END_REF]. Our archaeobotanical findings are thus of great importance to complete our picture of the arboricultural practices and wood use in the ancient Near East. Below, the implications of our finding of the cypress wood in the Sasanian archaeological sites are discussed in more detail.

Ecology and distribution of C. sempervirens

C. sempervirens is a Mediterranean conifer with very small and disjunct populations mostly constricted to the Eastern Mediterranean [START_REF] Quézel | Biogéographie et écologie des conifères sur le pourtour méditerranéen[END_REF]. Outside the Mediterranean biogeographical region, in the continental Middle East, C. sempervirens is mainly present in a few localities in the valleys cutting the Alborz Mountains in Northern Iran [START_REF] Riedl | Cupressaceae[END_REF][START_REF] Klein | La végétation altitudinale de l'Alborz Centrale (Iran)[END_REF]. In these valleys, C. sempervirens populations benefit from a mild Mediterranean Pluviseasonal Oceanic bioclimate restricted to a few N-S oriented valleys in natural geographic defiles (Djamali et al., 2011b). The climate of these restricted geographical zones is, however, in clear contrast to the rest of the vast Iranian and Eastern Anatolian plateaus which have a very continental climate with a long dry season (Djamali et al., 2011b[START_REF] Djamali | Climatic determinism in phytogeographic regionalization: a test from the Irano-Turanian region, SW and Central Asia[END_REF]. [START_REF] Browicz | Chorology of trees and shrubs in South-West Asia and adjacent regions[END_REF] suggests that the isolated populations of C. sempervirens in southern Iran including a cypress population (N100 trees) to the south of Shiraz might also be native to the region, thus forming relict populations. The presence of natural relict populations of Mediterranean elements in the southern Zagros Mountains has been known for a long time, as many isolated populations of Myrtus communis and some other Mediterranean trees can be found here and there [START_REF] Migliore | From Mediterranean shores to central Saharan mountains: key phylogeographical insights from the genus Myrtus[END_REF][START_REF] Akhani | First observations on the flora and vegetation of three islands in the NW Persian Gulf (Iran)[END_REF]. Although cypress is represented by only one species in Iran (C. sempervirens), Iranian botanists distinguish three cypress varieties including i) pyramidal/fastigiate cypress: C. sempervirens var. pyramidalis (O.Targ.Tozz.) Nyman called 'sarv-e Shirazi' and 'sarv-e Kashi' in Persian, ii) cereiform cypress or C. sempervirens L. cv. cereiformis Rehd. called 'sarv-e nāz' in Persian, and iii) horizontal cypress or C. sempervirens var. horizontalis (Mill.) Loudon called 'zarbin' in Persian [START_REF] Alam | Cypress. Encyclopedia Iranica[END_REF]. While the first two varieties are commonly planted as ornamental species in Iranian gardens, the third variety grows naturally in the Alborz Mountains. Today, all varieties and subspecies of C. sempervirens are, however, simply synonymous to Cupressus sempervirens L. according to The Plant List version 1.1. (http://www.theplantlist.org).

C. sempervirens, like many other Mediterranean elements, is sensitive to winter freezing temperatures and long dry seasons preventing its natural development under the Mediterranean xeric-and desertic continental bioclimates that dominate the Middle Eastern inlands (Djamali et al., 2011b). However, it is among the most resistant Mediterranean trees to winter temperature minima, whose leaves and cambium resist temperatures as low as -16 and -29 °C, respectively [START_REF] Quézel | Ecologie et biogéographie des forêts du bassin méditerranéen[END_REF]. This is why it can grow under irrigation in most parts of the continental Middle East. The relative indifference of C. sempervirens to the soil type also allows this tree to grow easily in different parts of Iran including the Fars region (e.g. [START_REF] Herbert | Some Years Travels Into Divers Parts of Asia and Afrique[END_REF][START_REF] Alam | Cypress. Encyclopedia Iranica[END_REF].

Cypress tree in the Iranian culture

The cypress tree has often been associated with Iranian traditional beliefs. One of the pioneers of western studies on oriental religions, [START_REF] Lajard | Recherches sur le culte du cyprès pyramidal chez les peuples civilisés de l'antiquité[END_REF], already linked the cypress tree to Zoroastrianism, the official religion of the Sasanian Empire, and quoted a great number of Persian and Arabic texts and short passages by Western orientalists that corroborated his statement. According to the tradition, the Prophet Zoroaster introduced the cypress cult in Persia by planting several cypress trees in the cities of ancient Khorasan in Northeastern Iran. In particular, the cypress planted by Zoroaster (or by his patron King Vishtaspa/Goshtasp) in front of the fire temple of the city of Kashmar in Khorasan (the so-called 'Cypress of Zoroaster') is of special significance for Zoroastrianism [START_REF] Stausberg | Faszination Zarathushtra. Zoroaster und die Europäische Religionsgeschichte der Frühen Neuzeit[END_REF][START_REF] Williams | The Zoroastrian Flame. Exploring Religion, History and Tradition[END_REF]. The origin of this mythical tree is uncertain, but according to a Persian tradition Zoroaster received the cypress directly from heaven. King Vishtaspa, who was one of the earliest followers and supporters of Zoroaster, ordered all governors of his empire to come to the feet of the cypress tree to listen to Zoroaster, to adopt his message and to abandon the cult of the idols of Turan and China.

These legendary tales suggest that the cypress tree assumed an important symbolic and cosmic significance in the course of history of Zoroastrianism. The pyramidal form of cypress also reminds us of the flames which go from earth to heaven, thus perfectly symbolizing the Zoroastrian doctrine itself. Some scholars go even further and argue that the 'paisley' motif ('boteh' in Persian), the twisted teardrop-shaped motif frequently used on Persian and Central Asian textiles, might be a stylized form of a cypress tree bent under wind [START_REF] Lajard | Recherches sur le culte du cyprès pyramidal chez les peuples civilisés de l'antiquité[END_REF]. This interpretation might be confirmed by the presence of this motif on old silk tissues dating back to the Sasanian period [START_REF] Eduljee | Boteh. Zoroastrian Heritage[END_REF]. Most probably, cypress trees had a strong symbolic and cosmic meaning also in Sasanian Iran and at the Sasanian court, where the kings made large use of Zoroastrian symbolism in their political propaganda (e.g., in their coinage) [START_REF] Schindel | Sasanian coinage[END_REF][START_REF] Daryaee | Sasanian Persia, the Rise and Fall of an Empire[END_REF]. Thus, it can be easily hypothesised that the Sasanian kings, who are known to have promoted agriculture and arboriculture [START_REF] Pope | A Sasanian garden palace[END_REF][START_REF] Wiesehöfer | Ancient Persia: From 550 BC to 650 AD[END_REF], particularly fostered the plantation of cypress trees in gardens and parks throughout their empire. The literary representation of this scenario is offered by a passage in the epic poem Shahnameh ('The Book of Kings'), in which Ferdowsi (ca. 940-1020 CE) tells the story of the mythical Iranian king Fereydun who decorated the world with cypresses (and roses) in order to turn it into an earthly paradise [START_REF] Khaleghi Motlagh | Abu'l Qasem Ferdowsi: The Shahnameh (The Book of Kings) 1[END_REF]. Cypress is frequently mentioned as a typical garden tree in later Persian poetry and continues to be an essential element of Persian gardens even today (Mahmoudi Farahani et al., 2016). It is very common to find centennial cypress trees in famous gardens in Iran such as the Fin Garden of Kashan, and shoots of the 'Cypress of Zoroaster' are still worshipped in some places such as Cham and Abarkuh [START_REF] Langer | Pīrān und Zeyāratgāh: Schreine und Wallfahrtsstätten der Zarathustrier im neuzeitlichen Iran[END_REF]. Today, the tree is considered a symbol of immortality [START_REF] Moynihan | Paradise as a Garden: In Persia and Mughal India[END_REF] because of its evergreen character that distinguishes it from many other deciduous trees planted in gardens.

Interestingly, the Middle and New Persian word for cypress is 'sarv' and the word 'Sarvestan' thus signifies the 'Garden of Cypress'. This suggests the presence of large-scale cypress plantations in the Sarvestan area at least since the foundation of the town of Sarvestan and probably also in earlier times.

Cypress wood use in ancient architecture

Cypress wood was used in a variety of ways in antiquity. In Mesopotamia, cypress (Akkadian: šurmenu) was considered a valuable wood, second only to cedar (i.e. Cedrus cf. libani). It was mainly imported as timber for palatial and temple architecture and shipbuilding since the 3 rd millennium BCE [START_REF] Moorey | Ancient Mesopotamian Materials and Industries: The Archaeological Evidence[END_REF]. The royal inscriptions from Mesopotamia attest that cypress also started to be cultivated as ornamental tree in the gardens of the Neo-Assyrian and Neo-Babylonian palaces in the 1 st millennium BCE. Cypress wood continued, however, to serve as constructional material especially for roofing and tall doors. Timber of cypress was probably also used in Achaemenid palatial architecture. When describing the palace of Ecbatana (modern Hamadan), Polybius (2011) reports that all the woodwork was made of cedar and cypress (Polyb. 10.27). Therefore this species was probably already planted in the Persian gardens and parks by the Achaemenids. The conifer, represented as scene-divider on the reliefs of Persepolis, is traditionally believed to be cypress [START_REF] Farrar | Gardens and Gardeners of the Ancient World: History, Myth and Archaeology[END_REF], and the Greco-Roman authors also mention cypress among all sorts of trees planted in the Persian paradeisoi. [START_REF] Strabo ; Stuiver | Extended 14 C database and revised CALIB radiocarbon calibration program[END_REF], for example, refers to the presence of cypress in the gardens and parks of Babylon in the Late Achaemenid period (Strab. 16.1.11), and Plutarch (1926) records a mixed forest of cypress and pine in a park in Northern Media (Plut. Artax. 25.1). The cypress trees of the Babylonian parks are said to have been cut by Alexander the Great to build part of his fleet (Strab. 16.1.11 Strabo, 1930).

Cypress wood was particularly appreciated in the Greco-Roman world for its durability, insect-repelling properties, and resistivity to humidity and seawater [START_REF] Lieutaghi | Le livre des arbres, arbustes et arbrisseaux[END_REF]. These properties combined with the fragrant scent of the tree are among the most important characteristics of cypress wood as timber (The Wood Database available at: http://www.wood-database.com/mediterranen-cypress/). Theophrastus (1916) praises cypress (Greek: kupárittos) for its durability and because it retains a fine and durable polish. He suggests that cypress wood was frequently used in house building and to construct the doors of the Greek temples (Theoph. Hist. pl. 5.7.4, 5.4.2 Theophrastus, 1916;cf. also Plin. HN 16.79 Pliny, 1945). The doors of the temple of Artemis at Ephesus, one of the Seven Wonders of the Ancient World, were made of cypress wood. Cypress trees often also adorned the surroundings of the temples of Artemis and of other gods in ancient Greece (Paus. 2.11.6, 2.13.3, 2.15.2, 3.22.9, 8.41.4, 10.38.9 Pausanias, 1918[START_REF] Pausanias | Description of Greece. I: Books 1-2[END_REF][START_REF] Pausanias | Description of Greece. I: Books 1-2[END_REF].

The Romans also used cypress (Latin: Cupressus) to construct buildings, statues and temple doors (Vitr. De arch. 2.18. 12 Vitruvius, 1931;Plin. HN 16.79 Pliny, 1945;Liv. 27.37.11 Livy, 1943). Vitruvius stresses in particular, the importance of furring strips and ties made of durable cypress wood in the construction of vaults (Vitr. De arch. 7.3.1 Vitruvius, 1934). The insect-repelling character of cypress wood (Plin. HN 16.80-81 Pliny, 1945;[START_REF] Nardi Berti | La struttura anatomica del legno ed il riconoscimento dei legnami italiani di piú corrente impiego[END_REF] greatly enhances the lifetime of its timber, and this fact surely increased its value as building material in the eyes of the Romans. In Roman Italy, cypress trees were probably not only cultivated in gardens, for their beauty and shade (Tac. Hist. 2. 78 Tacitus, 1925), but also in timber plantations. Recent studies have shown that it was utilised in the construction of luxurious Roman villas, such as the Villa of Poppaea at Oplontis near Naples and the houses of Herculaneum [START_REF] Moser | Archaeobotany at Oplontis: woody remains from the Roman Villa of Poppaea (Naples, Italy)[END_REF][START_REF] Moser | Timber economy in the Roman Age: charcoal data from the key site of Herculaneum (Naples, Italy)[END_REF]. In Roman times, cypress wood was largely used in shipbuilding and other structures in long contact with water [START_REF] Allevato | Woodland exploitation and roman shipbuilding, first data from the Shipwreck Napoli C (Naples, Italy)[END_REF][START_REF] Colaianni | Well lining in ancient hydrological system of Fratte settlement (Salerno, Italy): cypress wood use in the Archaic period (6th-5th century BC)[END_REF][START_REF] Sadori | Archaeobotany in Italian ancient Roman harbours[END_REF]. Cupressus wood is known for its resistance to fungal activity in humid environments [START_REF] Okino | Technological characterization of Cupressus spp. wood. Floresta e Ambiente[END_REF] and its use in shipbuilding, well lining etc., indicates that this property was already recognised in ancient times. Several wrecks revealed the use of cypress wood in ancient naval architecture such as a plank and a wale in Fiumicino 1 [START_REF] Boetto | L'épave de l'Antiquité tardive Fiumicino 1: analyse de la structure et étude fonctionnelle[END_REF], the keel in the Bourse wreck in Marseilles [START_REF] Rival | La charpenterie navale romaine: matériaux, méthodes, moyens[END_REF], planks and wales in Pantano Longarini [START_REF] Kampbell | The Pantano Longarini shipwreck: a reanalysis[END_REF]. Cypress wood was also the most used species for the construction of the medieval ship Dor 2001/1 [START_REF] Liphschitz | Wood species used in ancient shipbuilding in Turkey: evidence from dendroarchaeological studies[END_REF].

Cypress wood in the Sasanian palaces

The use of cypress wood as building timber in the Sasanian palaces (and especially in vaulting) may have been recommended and/or practiced by the Roman engineers/workers, who were captured during the Parthian and Sasanian campaigns against the Roman Empire. We know from the ancient written sources that the Romans were well aware of the technological properties of cypress wood as timber and that in Roman architecture it was particularly appreciated in the construction of vaults (Vitr. De arch. 7.3.1). Moreover, the use of the technical skills and manpower of Roman prisoners of war can be seen in the construction of many bridges, dams and drainage/irrigation systems in Sasanian cities, which bear Roman architectural elements [START_REF] Huff | Architecture: III: Sasanian period[END_REF]. Timber was most probably used in the Sasanian bridges in horizontal structures resting on the stony piers [START_REF] Kleiss | [END_REF]. However, considering the old tradition of cypress use in ancient Mesopotamia and its cultic value in ancient Iran -with possible extensive plantations in Persis and in other regions of the Empire -and also the indigenous status of C. sempervirens in Northern Iran and possibly also in Fars, it is highly probable that the Persians appreciated this wood as timber long before their first contacts with the Romans.

The most used building materials in Sasanian architecture are rubble with gypsum mortar and also mud bricks. Timber seems to have been a minor component of the Sasanian buildings, as it has not so far been reported from archaeological excavations of Sasanian structures. Our results show, however, that Sasanian architects did use some timber in their constructions. Although it is still not clear what was the exact function of timber in the Sasanian buildings, we know that in Islamic architecture, wood was used to support the vaulting systems during their construction and later for their reinforcement [START_REF] Kleiss | [END_REF]. Such utilisation may have been inspired by Sasanian architecture, where timber would have been used both for supporting vault systems and domes under construction and in later restoration works, especially in the remaking of facades. Sasanian palatial architecture took, in turn, inspiration from the Late Parthian brick masonry [START_REF] Hauser | The Arsacids (Parthians)[END_REF][START_REF] Callieri | Architecture et Représentations dans l'Iran Sassanide[END_REF] and it may also have been influenced by earlier Near Eastern traditions, especially in the technology and use of locally available materials. The use of timber to reinforce mud brick and stone constructions, for roofing, doors and columns is attested to in the ancient Near East since the Neolithic period and characterises, in particular, the palatial and temple architecture of Mesopotamia and Elam in historical times [START_REF] Moorey | Ancient Mesopotamian Materials and Industries: The Archaeological Evidence[END_REF]. In particular, the wood most commonly used in construction was coniferous wood for temples and palaces, while poplar was used for common buildings.

The choice of cypress wood as construction material in Sasanian palaces may have been determined by: i) the availability of the tree as a commonly planted species in parks, gardens, and along watercourses, ii) the absence or rarity of alternative suitable timbers, iii) the rectitude of its trunks, and iv) the well-known properties of the tree such as the resistance to insect attack, its high bending strength and moderate shrinkage values [START_REF] Bektaş | Principal mechanical properties of cypress wood (Cupressus sempervirens L.) naturally grown in Kahramanmaraş, Eastern Mediterranean of Turkey[END_REF][START_REF] Dogu | Wood identification of wooden marine piles from the ancient Byzantine port of Eleutherius/Theodosius[END_REF]. The latter property is probably important in Sasanian architecture because, to avoid the formation and development of fissures and cracks in the gypsum mortar used in the vaults and domes, the wood should present the least variations in volume. It can also not be excluded that the exclusive choice of cypress wood for the construction of important buildings such as royal palaces and temples was favoured due to the symbolic and religious meaning of this tree for the ancient Iranians.

Conclusions

Wood preserved in the walls from mud brick-and stone-dominated Sasanian palaces and palace-like buildings of Fars have proven to be highly informative. On the one hand, it has shed more light on their history, on the other, it has also increased our knowledge of the use of wood in Sasanian architecture, as well as the history of arboriculture in ancient Persia. Our results suggest that:

1. the Qal'a-ye Dokhtar palace-fortification dates back to the period corresponding to Parthian-Sasanian power transition;

2. the Palace of Ardashir I near Firuzabad also dates back, at the very latest, to the very beginning of Sasanian rule; 3. the Palace of Sarvestan must be dated to the transitional period from the Sasanian to the Islamic domination over Iran (7th century CE) and was actively used during the Early Islamic period for a long time, which corresponds to several archaeological phases. It is a Sasanian heritage in the post-Sasanian period; 4. cypress was among the major cultivated trees in ancient Iran. Its plantation was most probably not only due to its shade and beauty but also to its timber; 5. the sacred status of cypress in Zoroastrianism may have been one of the reasons for its cultivation in gardens and possibly also in the vicinity of temples in ancient Persia.

Our initial results showed that wood fragments in the archaeological materials from palaces and other residential or religious monuments constitute an important source of information to understand the history of architecture and agricultural practices in the Antiquity. Future investigations of archaeological monuments dating to this period should thus include a systematic study of wood fragments in parallel with other archaeological techniques.

Fig. 1 .

 1 Fig. 1. a. The Sasanian Empire during the reign of Shapur I (240-270 CE). b. Position of the study sites: QD: Qal'a-ye Dokhtar, AP: Palace of Ardashir I, SP: Palace of Sarvestan.

Fig. 2 .

 2 Fig. 2. Architectural plans of a) the Qal'a-ye Dokhtar, b) the Palace of Ardashir I, and c) the Palace of Sarvestan with the exact position of the studied timber fragments (star). The sample codes are those of the Poznan Radiocarbon Laboratory (see Table1). Plans are re-drawn after[START_REF] Hugi | Ghala-Dokhtar-Atechkade: monuments sasanides près de Firouzabad, Iran. Rapport préliminaire des projets de consolidation et restauration[END_REF] for the Palace of Ardashir I and Qal'a-ye Dokhtar and from[START_REF] Bier | Sarvistan: a study in Early Iranian architecture[END_REF] for the Palace of Sarvestan. See Materials and methods for more details.

Fig. 3 .

 3 Fig.3. Radiocarbon age probability distribution curves for timber wood fragments of Sasanian monuments. Note that both Qal'a-ye Dokhtar and the Palace of Ardashir I date to the very beginning of the Sasanian Empire and that the end of the latter palace may postdate the former one. Also note that despite their bimodal distribution, the probability curves for Sarvestan all fall in the post-Sasanian and Early Islamic period.

Fig. 4 .

 4 Fig. 4. Wood samples from the Palace of Sarvestan. a. Sarv-1 from Room 9. b. Sarv-3 from Room 1. a-1 and b-1 show the wood anatomies of Sarv-1 and Sarv-3 samples in radial crosssection. Absence of horizontal tracheids in the rays, particular bordered pits on radial faces of vertical tracheids and cupressoid type pits in parenchyma cells are typical features of Cupressus.

Table 1 14

 1 C age report for the studied wood samples. Calibrated radiocarbon ages are reported in CE (Common Era) ages. Radiocarbon age ranges with highest probabilities are highlighted in bold. Calibrations were performed in Calib 7.1[START_REF] Strabo ; Stuiver | Extended 14 C database and revised CALIB radiocarbon calibration program[END_REF] based on Intcal13 calibration dataset[START_REF] Reimer | IntCal13 and MARINE13 radiocarbon age calibration curves 0-50000 years cal BP[END_REF].

						Calibrated 14 C age ranges (CE)
	Sample	Provenance	Context	14 C lab code	14 C age (yr BP)	1σ range	2σ range
	Sarv1	Palace of Sarvestan	Room9	Poz-80705	1315 ± 30	660-694 (0.72)	655-724 (0.74)
						704-706 (0.01)	739-767 (0.26)
						746-763 (0.27)
	Sarv2		Room10	Poz-80706	1280 ± 30	682-720 (0.59)	664-773 (1)
						741-767 (0.41)
	Sarv3		Room1	Poz-80707	1220 ± 30	726-738 (0.1)	692-747 (0.22)
						768-779 (0.11)	762-887 (0.78)
						789-870 (0.78)
	Dokht1	Qal'a-ye Dokhtar		Poz-80708	1835 ± 30	134-216 (1)	86-110 (0.05)
							115-246 (0.95)
	Ard1	Palace of Ardashir I		Poz-80709	1820 ± 30	139-197 (0.66)	90-99 (0.01)
						207-234 (0.34)	124-257 (0.95)
							296-320 (0.04)

Acknowledgment

This study was funded by the Franco-German ANR-DFG project entitled "PALEO-PERSEPOLIS" (ANR-14-CE35-0026-01). We wish to thank the fruitful comments of Prof. Naomi F. Miller from the University of Pennsylvania and Prof. Josef Wiesehöfer from the University of Kiel. We are also grateful to Dr. Belinda Gambin (University of Malta) for her comments and the English-editing of the text.