Pour modéliser le volume du sol urbain, on a recours aux observations archéologiques, géologiques et géotechniques disponibles à Tours. Ces observations sont représentées par des points au sein de la zone d’étude retenue.
Une base de données géographique a été constituée avec ArcGIS pour représenter ces observations à l’échelle de l’espace urbanisé ancien : BDdepot1. Les champs de cette couche sont :
- les coordonnées géographiques, X et Y
- un code identifiant
- la source de l’observation
- la cote altimétrique (en m NGF) du niveau actuel
- la cote altimétrique (en m NGF) du toit du niveau naturel (alluvions de la Loire)
- l’épaisseur du dépôt archéologique (en m)
Au total, 145 points ont été retenus dans l’espace étudié (Fig. 1). Le premier critère de sélection était d’obtenir des informations pour renseigner l’ensemble des champs. Le deuxième critère était de retenir que les points contenus dans la zone d’étude pour réaliser les cartes d’estimation du dépôt archéologique et du toit du niveau naturel.
On observe que ces données sont irrégulièrement réparties dans l’espace et que l’épaisseur varie entre 0,75 et 12 mètre. L’épaisseur moyenne est de 3,85 m, avec un écart-type de 1,9. Il faut noter que 80% des observations ont une épaisseur inférieure à 6 mètres.
Par ailleurs, le toit des alluvions varie entre 41,2 et 50,7 m NGF. La cote moyenne du toit des alluvions est de 46,14 m NGF avec un écart-type de 1,39. 80% des effectifs ont le toit des alluvions inférieur à 47 m NGF.
Un point important est que cette base contient 2 points dont les mesures ont été extrapolées à partir de données voisines (Fig. 1, en étoile) : ces données ont été construites afin de donner une valeur limite à l’espace étudié. C’est pourquoi, une autre base de donnée excluant ces points extrapolés a été établie et nommée BDdepot2. Cette deuxième base sera utilisée pour analyser la structuration spatiale des observations sur l’épaisseur et le toit des alluvions.

Pour construire une carte des épaisseurs de ce dépôt, il faut générer des surfaces à partir des  points par une interpolation*. C’est une opération qui consiste à déterminer, à partir de valeurs discontinues, de nouvelles valeurs correspondant à un caractère intermédiaire pour lequel aucune mesure n’a été effectuée. Le problème est la généralisation des observations ponctuelles à toute la ville. Il est donc nécessaire d’utiliser une méthode d’interpolation adaptée à la répartition spatiale des données et à la nature même de ces données.

Pour pouvoir représenter la variabilité spatiale de ces épaisseurs, l’interpolation spatiale est la procédure la mieux adaptée. Elle permet l’estimation de la valeur d’une variable en un point à partir des valeurs de la même variable mesurée en des points voisins (Arnaud, Emery 2000 ; Zaninetti 2005 : 213). Il existe deux approches pour l’interpolation :
- une approche déterministe où la pondération est définie à priori (fonction Spline, fonction inverse de la distance, polygones de Thiessen…)
- une approche probabiliste où la pondération est définie à partir des données (Krigeage).
Celle qui semble la mieux adaptée à l’étude est celle du krigeage. Cette méthode porte le nom de son précurseur, D.G. Krige, ingénieur minier, qui a développé des méthodes statistiques empiriques afin de déterminer la distribution spatiale de minerais à partir de forages à partir des années 1950 (Krige 1951). C’est G. Matheron (Matheron 1970) qui formalise l’approche « en utilisant les corrélations entre les forages pour en estimer la répartition spatiale » (Gratton 2002 : 1). Il met au point la théorie des variables régionalisées, connue sous le nom de géostatistique*. Elle est l'application des probabilités aux phénomènes naturels qui se développent dans l'espace et dans le temps. Le krigeage est un estimateur géostatistique qui devient alors une référence auprès des géographes, des sciences de la terre et des sciences humaines (Zaninetti 2005 : 233-234). Ce modèle probabiliste donne également un indicateur de la précision de l'estimation. « Cet indicateur, appelé variance d'estimation, est un outil précieux car il ouvre la voie à une possible gestion des incertitudes (traduites en termes de variance) »  (Sandjivy 2003).

Le krigeage se différencie par les autres méthodes déterministes (Spline, distance-inverse, méthode des plus proches voisins, polygone de Thiessen…) du fait qu’il fait précéder l’interpolation d’une analyse variographique, qui consiste à élaborer un modèle de covariance ou de variogramme à partir de données et d’informations disponibles (Arnaud, Emery 2000 : 179 ; Zaninetti 2005 : 234).
En effet, l’intérêt de la variographie est de pouvoir analyser et modéliser le degré de ressemblance (autocorrélation) entre la mesure d’un phénomène obtenu entre deux points. En particulier, la variographie conduit à identifier la distance maximale entre deux points de mesures au delà de laquelle, il n’y a plus de corrélation.
En fait, pour les différentes approches d’interpolation spatiale, l’estimation d’une valeur est calculée suivant une fonction F où chaque valeur est pondérée. Le problème consiste à déterminer la pondération de chacun des points environnants. Dans le cas du krigeage, le poids correspond au degré de ressemblance obtenu par l’analyse variographique (Gratton 2002 : 1).
 
Pour mettre en œuvre le krigeage, il faut répondre à différentes conditions (Arnaud, Emery 2000) :
- la moyenne et la variance de la fonction F doivent être stationnaires, c’est-à-dire qu’elles ne dépendent pas de la position des points, seulement la distance entre les points
- les données voisines du point à estimer, appelé voisinage, doivent être en nombre suffisant pour conduire à une estimation précise et peu sensible à l’ajout ou au retrait d’une donnée périphérique.
En posant l’hypothèse de stationnarité, la méthode qui convient pour l’étude est le Krigeage ordinaire. 
 

Pour analyser la structuration spatiale de l’épaisseur du dépôt archéologique, il faut en tracer un semi-variogramme qui est un graphe qui présente le demi-écart quadratique moyen (fonction d’ajustement choisie ici) de couples de points en fonction de leur distance selon la relation suivante :

où zi et zi + h sont les positions des mesures et Qd(zi) est la valeur de résistance de pointe correspondant au point considéré. La variable h correspond à la distance entre deux mesures. La sommation est étendue à tous les couples de points (zi, zi + h) distants de h, N(h) représentant l’effectif de tels couples de points (Fig. 2).

Le but est de déterminer la distance de corrélation spatiale entre les points d’observations,  soit la portée a du semi-variogramme. 
Pour réaliser le variogramme expérimental, il faut définir plusieurs paramètres :
- la portée limite qui définit la limite de la recherche (moitié de la fenêtre d’étude) ;
- le nombre de pas de recherche ;
- la distance entre les pas de recherche ou distance de voisinage ;
- la direction du variogramme, omnidirectionnel ou directionnel ;
- le champ du cône d’orientation de l’orientation retenue.
 
Deux logiciels ont été utilisés afin d’élaborer le variogramme de l’épaisseur du dépôt archéologique : l’extension Geostastical analyst d’ARGIS© (ESRI) et le logiciel VARIOWIN© (Pannnatier 1996, Anselin 2003)
La portée limite a été définie à 1200 m. Pour déterminer la distance de voisinage, nous avons estimé la densité de points théorique. On obtient 1 point tous les 128 mètres (143 points sur 235,27 ha étudié ; la zone étant un polygone représentatif de l’espace urbanisé ancien). Nous avons alors choisi 100 m de distance de voisinage qui correspond à la tolérance selon J.-M. Zaninetti (Zaninetti 2005 : 237). Le nombre de pas utilisé est 12.
 
Dans un premier temps, on obtient un cercle variographique qui permet de déterminer la direction du variogramme (Fig. 3). A priori, on observe que le variogramme des épaisseurs est anisotrope.

Le variogramme de BDdepot2 est tracé suivant les directions (70°/160° avec un angle de tolérance de 10°), le pas d’étude (100 mètres) et le nombre de pas (12) (Fig.  4). Pour chaque pas de recherche, on possède le nombre de paires de points[1] qui ont contribué à calculer la variance correspondante.

[1] Il est conseillé d’avoir au minimum 20 paires de points de mesure ; à partir de 50 paires, les résultats sont robustes.

On ajuste alors un modèle théorique sur le variogramme expérimental (Fig.  5). Celui qui est retenu est un modèle exponentiel avec un effet de pépite de 0,4, une portée ade 300 m et une variance totale C de 4,6, suivant la relation γ(h) = C*[1- exp (-h/ a)], avec h, la distance,
soit γ(h) = 0,4 + 4,2*[1- exp (-h/300)] pour la direction de 70° et γ(h) = 0,4 + 4,2*[1- exp (-h/200)] pour la direction de 160°. L’opération est délicate puisqu’elle dépend de l’objectivité de l’utilisateur.

On observe que :
- d’une part, l’effet de pépite correspond à une micro-structuration des données que la variogramme ne peut mesurer. Il correspond au 1/10^e de la variance total donc il est plutôt négligeable ;
- d’autre part, la distance de corrélation entre les points de mesure est estimée à environ 300 mètres. Pour une étude globale de la ville, ce résultat est plutôt raisonnable. Par contre, à l’échelle locale, il est à prendre avec des précautions puisque l’interpolation ne tient pas compte des effets des composantes socio-spatiales. Par exemple, la présence d’une enceinte présente une rupture significative dans la production de sol. L’idéal serait de réaliser une interpolation au sein de chaque enceinte urbaine. Ceci n’a pas pu être réalisé car le nombre de point est insuffisant pour certaines d’entre elle (enceinte des 3^e- 4^e s, du 10^e s. et du 12^e s.).
 
Il faut noter que pour la ville de Bordeaux, la portée est estimée à 500 mètres  et le krigeage a été effectué sous contrainte d’inégalité qui est une « technique qui permet d’utiliser une variable qui n’est pas définie par une valeur exacte en un point, mais par un intervalle dont les bornes constituent les données d’inégalité » (Dominique, Marache, Régaldo-Saint Blancard soumis ; Langlais 1990). Cette procédure consiste à simuler, au droit du point pour lequel on dispose une donnée d’inégalité, une valeur qui  appartient « à l’intervalle défini par les données d’inégalité et les données exactes » (Dominique, Marache, Régaldo-Saint Blancard soumis: 10) puis à intégrer ces points simulés aux points mesurés. Ainsi, elle permet de mettre en œuvre  des observations qui n’ont pas atteint les alluvions du Quaternaire ou le substratum.
 
Ensuite, il faut  déterminer le voisinage : il peut être unique si on veut prendre en compte tous les points dans le calcul de l’estimation, ou bien, glissant si on choisit une distance limite de recherche. Dans ce cas, le voisinage est unique et correspond aux 145 points d’observations.
Pour finir, il faut réaliser une validation croisée du modèle choisi afin de s’assurer de sa robustesse. Elle consiste à estimer la valeur de l’épaisseur en chaque point appartenant à la base de données en les éliminant temporairement et en appliquant le variogramme ajusté. Les données sélectionnées pour calculer les estimations sont celles de BDdepot1 (Fig.  6). On étudie la relation entre une valeur mesurée et estimée à partir du modèle en chaque point.

Le coefficient de corrélation (r) entre les valeurs mesurées et estimées est égal à 0,67. Celui obtenu pour Bordeaux présente le même résultat (0,71). L’erreur moyenne est correcte puisqu’elle est proche de zéro (-0,017). Enfin, la moyenne des écarts entre les valeurs mesurées et estimées est de 0,96 m.

La carte d’interpolation spatiale de l’épaisseur du dépôt archéologique a été élaborée avec ArcGIS, à partir du modèle défini ci-dessus (Fig. 7). L’épaisseur estimée varie de 1 à 10,8 mètres. Elle est plus importante au sein de la zone comprise dans l’enveloppe de la ville du 14^e siècle (environ supérieure à 5 mètres). Ces résultats sont proches des hypothèses de départ : une intensité des activités humaines au sein de la ville antique et médiévale traduite par une épaisseur plus importante que dans la ville extra-muros.

La carte des erreurs standard montre que les estimations calculées au bord de la zone d’étude sont peu fiables car il y a peu de points d’observation (Fig. 8).

Il est important de préciser le résultat dépend de la résolution choisie pour la représentation des données.
Enfin, puisqu’il s’agit d’un modèle, l’objectif de cette carte est de poser des hypothèses sur l’épaisseur à l’échelle de la ville. A l’échelle du site, l’archéologue doit revenir aux points d’observations et travailler sur des coupes stratigraphiques restituées.

La carte du toit du niveau naturel a été élaborée de la même façon que celle de l’épaisseur. L’étude de la structuration spatiale a été effectuée à partir des observations de BDdepot2 puis l’interpolation a été réalisée à partir ce ceux de BDdepot 1.

Les mêmes paramètres sont conservés : un pas de mesure de 100 mètre et un nombre de pas de 12. Ici, le champ sélectionné correspond au toit des alluvions.
Le cercle variographique montre que le variogramme du toit des alluvions est omnidirectionnel (Fig.  9).

Ensuite, un modèle de variogramme est ajusté à partir du variogramme expérimental (Fig.  10). Le modèle retenu est sphérique avec un effet de pépite de 0,6, une portée a de 175 mètres et une variance totale C de 1,9 où γ(h) = C * [(3(h)/2 a) +( h³/2 a ³ )] lorsque h < a et où γ(h) = C lorsque h  ≥ a,
Soit γ(h) = 0,6 + 1,3* [(3(h)/2*175) +( h³/2*175³ )] lorsque h < 175 et où γ(h) = 1,6 lorsque h  ≥ 175.

Ainsi, la distance de corrélation théorique entre les points de mesure est de 175 mètres. Elle plus courte que celle des épaisseurs. Dans ce cas, l’effet de pépite est plus important (30% de la variance totale). Il y a donc des variations à micro-échelle plus importantes que pour les épaisseurs.
Ici, le voisinage est glissant et a été choisi sur une distance de 175m.
 
Enfin, la validation croisée des données mesurées et estimées par le modèle, pour la base BDdepot1, présente un coefficient de corrélation plutôt moyen (0,50) avec une erreur moyenne satisfaisante de 0,03 et une moyenne des écarts de 0,89 m (Fig.  11).

La carte d’interpolation du toit des alluvions est réalisée à partir du modèle présenté ci-dessus et des points d’observations de BDdepot1 (Fig. 12). Les cotes estimées varient entre 43 et 49 m NGF. On observe que :
- le toit des alluvions est plus élevé dans la partie nord-est de l’espace urbanisé ancien (à l’emplacement de l’amphithéâtre) et les parties les plus hautes sont localisées sous l’espace qui constituait la ville au 14^e siècle. Cette zone peut être interpréte comme une montille ancienne de la Loire (Gay-Ovejero, Macaire, Seigne à paraître) ;
- une dépression apparaît en direction de la Loire qui correspond à l’espace gagné sur la Loire de l’Antiquité au 17^e siècle ;
- une autre dépression dans la partie sud que l’on peut expliquer par la présence de paléochenaux.

La carte des erreurs standard permet également de prendre des précautions sur la qualité des estimations  sur les bords de la zone d’étude et notamment, dans la zone située à l’est de l’enceinte du 14^e s. (Fig. 13).

On peut obtenir un état de connaissance archéologique d'une ville à partir de l'emprise et la profondeur des fouilles ainsi que des destructions portées au sous-sol. Ces données pour la ville de Tours ont été intégrées dans le SIG, intitulé ToToPI (Topographie de Tours Pré-Industriel) du Laboratoire Archéologie et Territoires (Fig. 14). Le croisement de ces données donne un aperçu de l'activité archéologique en ville qui est directement lié aux politiques d'aménagement et à la politique de prescriptions de sauvegarde patrimoniale des services de l'Etat (Rodier, Laurent 2007).

Globalement, on observe que la superficie fouillée (1,7 ha) correspond à 1,3% à de l'espace urbanisé ancien (130 ha). L'ensemble des surfaces détruites par des aménagements depuis l'après-guerre représente 21,8% de la zone d'étude. Les interventions archéologiques ne portent que sur 28% de l'emprise des destructions qui sont à leur origine.

La carte d'estimation de l'épaisseur du dépôt archéologique est une vision du potentiel théorique de la ville. La carte du potentiel archéologique réel prend en compte l'état de conservation du sous-sol. Cet outil d'aide à la décision est un document sur lequel peuvent s'appuyer les archéologues et les aménageurs lors d'un projet urbain. Il peut être représenté en plan (Fig. 15) ou en 3D (Fig. 16). Le modèle en 3D a été simplifié ici (Rodier 2007). Il est constitué :
- d'un plan parcellaire du 19^e siècle sur lequel le bâti a été extrudé d'une valeur arbitraire ;
- du volume des destructions du sous-sol qui ont la même valeur de profondeur;
- du volume des opérations archéologiques calculé à partir de profondeur maximale atteinte par les fouilles de chacune d'entre elle ;
- du MNT du toit des couches géologiques (alluvions).

On a démontré ici que l'emploi d'un SIG est utile pour l'évaluation du potentiel archéologique du sol en milieu urbain. La carte d'estimation de l'épaisseur du dépôt archéologique et les cartes du potentiel archéologique sont des outils d'aide à la décision pour l'aménagement urbain et pour la programmation de la recherche sur la ville. L'intérêt d'un SIG pour modéliser le sol urbain est de ne pas figer les informations mais de favoriser le renouvellement des connaissances.Le module Geostatistical Wizard d'ArcGIS permet d'éviter une pratique de "presse-bouton" pour réaliser des cartes d'interpolation à condition de définir objectivement les paramètres du modèle choisi.

Bibliographie :
 
Anselin L. - An Introduction to Variography using Variowin, www.sal.uiuc.edu/stuff/stuff-sum/pdf/variowin.pdf, University of Illinois, 17p, 2003.
 
Arnaud M., Emery X. - Estimation et interpolation spatiale. Méthodes déterministes et méthodes géostatistiques, Paris, Hermes Science Europe, 221 p., 2000.
 
Breysse D. et al. - Le pénétromètre et l'hétérogénéité des sols archéologiques urbains, Revue Française de Géotechnique, 100: 43-58, 2002.
 
Dominique S., Marache A., Régaldo-Saint Blancard P. - Morphologie naturelle du site de Bordeaux (France), Archéosciences, soumis.
               
Galinié H. et al. - Utilisation du pénétromètre dynamique de type PANDA en milieu urbain pour l'évaluation et la caractérisation du dépôt archéologique, Revue d'Archéométrie, 27: 15-26, 2003.
 
Gay-Ovejero I., Macaire J.-J., Seigne J. - Sondages carottés cave du Couvent des Dominicains, 12 rue du Général Meusnier, et au 5 rue Manceau, Tours (37).
 
Gratton Y. - Le krigeage : la méthode optimale d'interpolation spatiale, Les articles de l'Institut d'Analyse Géographique: www.iag.asso.fr, 2002.
 
Krige D.G. - A statistical approach to some basic mine valuation problems on the Witwatersrand, Journal of Chem., Metal. and Minning Soc. of South Africa, 52: 119-139, 1951.
Langlais V. - Estiamtion sous contraintes d'inégalités, Thèse de doctorat, E.N.S. des Mines de Paris, 1990.
 
Laurent A. - Evaluation du potentiel archéologique du sol en milieu urbain, mémoire de doctorat en histoire, spécialité archéologie (dir. H. Galinié), 4 volumes, Université de Tours, http://tel.archives-ouvertes.fr/, 2007.
 
Matheron G. - La théorie des variables régionalisées et ses applications, Les cahiers du centre de morphologie mathématique de Fontainebleau, Paris, 1970.
 
Pannatier Y. - VARIOWIN: Software for Spatial Data Analysis in 2D, New-York, Springer-Verlag, 1996.
 
Rodier X. - Modélisation du bâti ancien, des destructions et des fouilles, in : Galinié H. (dir.), Tours antique et médiéval., Lieux de vie, Temps de la ville, 40 ans d'archéologie urbaine; 30^e supplément à la RACF, n°spécial de la collection Recherches sur Tours, Tours, FERACF, 2007 : 43.
 
Rodier X., Laurent A. - La connaissance archéologique de la ville, in : Galinié H. (dir.), Tours antique et médiéval., Lieux de vie, Temps de la ville, 40 ans d'archéologie urbaine; 30^e supplément à la RACF, n°spécial de la collection Recherches sur Tours, Tours, FERACF, 2007 : 40 - 43.
 
Sandjivy L. - Procédé pour la détermination d'un indice de qualité spatiale de données regionalisées, Demande de Brevet Européen, publié le 03/09/2003, http://gauss.ffii.org/PatentView/EP1341094, n° EP1341094, 2003.
.
Zaninetti J.-M. - Statistique spatiale, méthodes et applications géomatiques, Paris, Hermes sciences, 321 p, 2005.