Modélisation et Statistique Spatiale

Atelier RASMA – Université Gaston Berger

Saint Louis du Sénégal
29 novembre – 4 décembre 2010

Xavier Guyon – SAMM -- Université Paris 1

• Géostatistique, modèle du second ordre, krigeage

• Donnée sur un réseau :
Auto-Régression Spatiale (SAR et CAR, SARX)
Champ de Gibbs – Markov – Auto modèle de Besag
Simulation par chaîne de Markov (MCMC)
• Processus ponctuel
Atelier Spatial RASMA 1
Saint Louis du Sénégal
Atelier Spatial RASMA 2
Saint Louis du Sénégal
 Données géostatistiques X
(a) Cumul de pluies dans 100 stations météo suisse le jour du passage du
nuage de Tchernobyl : réseau irrégulier (sic.100 de geoR)
(b) Porosité d’un sol (soil250 de geoR) : réseau régulier
La dimension des symboles est proportionnelle à X

Atelier Spatial RASMA 3

Saint Louis du Sénégal
 Le logiciel R
Installation de R :
[Link]

• site miroir : i.e. Toulouse

• Deux fenêtres : R Console (RGui) et R
Graphics

Chargement du package geoR (données

géo-stat)
Atelier Spatial RASMA 4
Saint Louis du Sénégal

1 - Installer le package geoR (RGui et site miroir

 Données Porosité (soil250)
22 variables « chimiques » sur une grille régulière 10x25 points
espacés de 5 mètres (cf. soil250 dans la liste de geoR).
On sélectionne la coordonnées n°16, ctc ( catium exchange)

> data(soil250)
> ctc <- [Link](soil250, [Link]=16)
> plot(ctc)
4 graphiques

1 - les 4 quartiles (4 couleurs) de CTC

2 et 3 - les nuages (ctc(x,y), y) et (x,ctc(x,y))
4 - Histogramme de répartition des 250 valeurs de ctc

Conservation d’un graphique :

se placer dans la fenêtre graphique
historique
Ajouter (ou précédent, etc….)

Autre solution : placer la commande « > x11() » avant une commande graphique
conservera le graphique (aller dans fenêtre, les graphiques sont numérotés
séquentiellement)

Atelier Spatial RASMA 5

Saint Louis du Sénégal
 Données pluviométrie Suisse
> print(sic.100)
> points(sic.100,borders=[Link])

[Link] : fichier frontière

4 données pour chaque station : coordonnées (x,y), hauteur
de pluie, altitude
sic.100 : 100 stations choisies au hasard dans un réseau
de 367 stations
[Link] : toutes les 367 stations

> points([Link], borders=[Link])

> plot([Link])

Atelier Spatial RASMA 6

Saint Louis du Sénégal
 250

250
200

200
150

150
Y Coord
Y Coord
100

100
50

50
0

0
-50

-50
0 50 100 150 200 250 300 350 0 100 200 300 400 500 600
X Coord data
600

0.004
500
400

0.003
300

Density
data

0.002
200

0.001
100
0

0.000

0 50 100 150 200 250 300 350

X Coord 0 100 200 300 400 500 600
data
Atelier Spatial RASMA 7
Saint Louis du Sénégal
 Questions en Geostatistique

• Quelle structure de corrélation spatiale ?

- Stationnarité (covariance) , isotropie ?
- Non stationnarité (variogramme)
- Modèle avec covariables (données exogènes)

• Estimation (validation) de modèle

• Prédiction partout : carte de krigeage, simulation

conditionnelle

• Outil logiciel : geoR

Atelier Spatial RASMA 8

Saint Louis du Sénégal
 Données réelles sur un réseau discret
(a) % groupe sanguin A dans 26 comtés Irlande(eire, spdep)
(b) Image 256 x 256 de J. Lennon (193 niveaux de gris, lennon du
packagefields)

Packages : spdep, fields, …

Atelier Spatial RASMA 9

Saint Louis du Sénégal
Atelier Spatial RASMA 10
Saint Louis du Sénégal
Atelier Spatial RASMA 11
Saint Louis du Sénégal
Atelier Spatial RASMA 12
Saint Louis du Sénégal
 Questions

• Quel modèle ?
- voisinages d’influence pour chaque site
- SAR ou CAR
- stationnaire ou non
- avec variables exogènes (SARX)

• Estimation et validation de modèle

• Tests sur les paramètres ….

• Outils logiciel : spdep, fields, …

Atelier Spatial RASMA 13

Saint Louis du Sénégal
Atelier Spatial RASMA 14
Saint Louis du Sénégal
 Modèle de Gibbs - Markov

Ex : répartition spatiale d’une espèce végétale

présence / absence de la grande laîche
Modèle de Auto - Logistique {0,1}
Voisinage de dépendance? Estimation, validation, tests ? Simulation

Atelier Spatial RASMA 15

Saint Louis du Sénégal
 Données Ponctuelles

x = configuration spatiale de n points

3 exemples

(a) – 97 fourmilières : données ants de spatstat

(b) – 42 centres de cellules d’une coupe histologique (cells)

(c) -- 126 pins d’une forêt finlandaise (finpines)

Package spatstat

Atelier Spatial RASMA 16

Saint Louis du Sénégal
Atelier Spatial RASMA 17
Saint Louis du Sénégal
PP bivarié : 2 types de cellules de la rétine du lapin
(294 en tout)
> data(betacells) et > plot(betacells)

Atelier Spatial RASMA 18

Saint Louis du Sénégal
 Questions sur les Processus ponctuels

• Répartition spatiale au hasard ([Link] Poisson = CSR pour

Complete Spatial Randomness))
• Ou non :
- avec compétition (chaque centre de cellule développe une
zone d’influence)
- avec coopération (i.e. agrégats autour d’un père)

• Homogénéité spatiale ou non

• Quels modèles explicatifs ?
• Statistique

Atelier Spatial RASMA 19

Saint Louis du Sénégal
Atelier Spatial RASMA 20
Saint Louis du Sénégal
Atelier Spatial RASMA 21
Saint Louis du Sénégal
 Champ du second ordre X sur S (L**2)

• Domaine d’étude : sites s de S, sous ensemble de R**2

• Observation X(s) réelle et de variance finie : Var(X(s))<∞

• X caractérisé par ses lois finies dimensionnelles

Moyenne : m(s) = E(X(s))
Covariance : c(s,t) = cov (X(s), X(t))

• Le plus souvent, modèle gaussien (pas une nécessité)

Atelier Spatial RASMA 22

Saint Louis du Sénégal
 Différents Bruits Blancs (BB)

• BB fort : variables {e(s)} i.i.d.

• BB faible : variables centrées et de même variances

• BB gaussien : BB faible gaussien

• BB coloré : variables centrées même variances mais

corrélées

Atelier Spatial RASMA 23

Saint Louis du Sénégal
Caractérisation d’une covariance :
la semi définie positivité (sdp)

Atelier Spatial RASMA 24

Saint Louis du Sénégal
 Champ gaussien X

Si toute combinaison linéaire est gaussienne

X spécifié par sa moyenne m(.) et sa covariance c(.,.)

Atelier Spatial RASMA 25

Saint Louis du Sénégal
 Champ stationnaire
• Moyenne constante
• Covariance invariante par translation :

Champ isotrope
covariance invariante par isotropie :

Atelier Spatial RASMA 26

Saint Louis du Sénégal
 Propriétés d’une covariance stationnaire C

• C est semi-définie positive

• I C(h) I ≤ C(0)
• X(As) est stationnaire si s
As est linéaire
• Une somme pondérée à coefficients >0 de covariances est
encore une covariance
• Si C est continue en 0, alors C est uniformément continue
partout

Atelier Spatial RASMA 27

Saint Louis du Sénégal
 Quelques covariances isotropiques
Portée a >0 et Variance σ**2 >0

• Pépitique: C(0)= σ**2 et C(h)=0 sinon

• Exponentielle : C(h) = σ**2 exp(-a II h II)

• Sphérique si d≤3

C(h) = 0 sinon

• Gaussienne :

(cf. liste assez complète dans > [Link])

Atelier Spatial RASMA 28

Saint Louis du Sénégal
 Modèle(s) de Matern

• Plus un paramètre ν contrôle la régularité de C en 0

(K est la fonction de Bessel de première espèce)
• ν = ½
cov. exponentielle
• ν = ∞
cov. Gaussienne
• Plus ν augmente, plus C(h) est régulière en 0 et plus X est
régulier (en moyenne quadratique)

Atelier Spatial RASMA 29

Saint Louis du Sénégal
 Champ intrinsèque et variogramme

• Considérer le champ des h-accroissements :

• X intrinsèque si ses h-accroissements sont

stationnaires

• Variogramme en h :

Atelier Spatial RASMA 30

Saint Louis du Sénégal
 Stationnaire ou intrinsèque ?
• Stationnaire
intrinsèque :

• Intrinsèque ≠> stationnaire :

Exemple : le mouvement brownien, γ(h) = I h I

• Un variogramme n’est pas toujours borné :

Exemple : vario puissance et auto-similarité (c = 1

pour le mouvement brownien).

Atelier Spatial RASMA 31

Saint Louis du Sénégal
Portée, palier, effet pépite d’un variogramme
(a) Les 3 caractéristiques d’un variogramme
(b) Variog. expo.(1), sphérique (2) et gaussien (3) : régularité
en 0 linéaire pour (1-2) et parabolique pour (3)

Atelier Spatial RASMA 32

Saint Louis du Sénégal
 Régularité du variogramme de Matern
• v = 1/2 donne le variogramme exponentiel
• v ↑, plus de régularité en 0
• v ≥ 2, C dérivable en 0 à dérivée nulle

Atelier Spatial RASMA 33

Saint Louis du Sénégal
[Link] de R
principales covariances spatiales
Exemple : expo., sphérique, gauss de même portée pratique.
variograms with equivalent "practical range"

1.0
exponential
spherical
gaussian

0.8
0.6
γ(h)

0.4
0.2
0.0

0.0 0.2 0.4 0.6 0.8 1.0

distance

Atelier Spatial RASMA 34

Saint Louis du Sénégal
 Anisotropies

Variogrammes différents suivant les directions

• Anisotropie géométrique :

• Anisotropie zonale:

Atelier Spatial RASMA 35

Saint Louis du Sénégal
 Anisotropies

(a) : géométrique (b) zonale

Atelier Spatial RASMA 36

Saint Louis du Sénégal
 La régulatité de la covariance C en 0 règle la
régularité en m.q. de X partout

• C continue en 0
X continu partout

• C’’ exite en 0
X dérivable partout

• Idem en remplaçant C par le variogramme γ

• Importance : régularité d’une carte de prédiction
(Krigeage) est fonction du choix de C (ou γ)

Atelier Spatial RASMA 37

Saint Louis du Sénégal
Atelier Spatial RASMA 38
Saint Louis du Sénégal
 Simulation d’un champ gaussien
le package RandomFields

GaussRF : simule un champ spatial ou spatio-temporel stationnaire

Il faut déclarer :
la fonction de covariance
la grille de simulation
la tendance si il y en a une
la méthode de simulation retenue
(cf. 1er exemple de champ stable pour différentes grilles)
CovarianceFct : donne liste des covariances/variogramme
spatiaux ou spatio-temporel
CondSimu : réalise la simulation conditionnelle d’un champ gaussien
en dehors des sites d’observation
ShowModels : démonstration interactive de simulation de modèles

Atelier Spatial RASMA 39

Saint Louis du Sénégal
Atelier Spatial RASMA 40
Saint Louis du Sénégal
Atelier Spatial RASMA 41
Saint Louis du Sénégal
 Prédiction à moyenne et covariance connue :
le Krigeage simple

• Vecteur des n observations X = (X(s(1)), X(s(2)), …, X(s(n)))

• Prédire X(s(0)) (carte de krigeage) partout sur S
• Choix d’une covariance connue (variogramme) C :
Σ = cov(X) et c = cov(X(s(0)),X)
• Reconstruction par MCO : minimiser

• si X de moyenne m connue, krigeage simple en s(0)

Prévision BLUP et Variance de prédiction :

Atelier Spatial RASMA 42

Saint Louis du Sénégal
Krigeage = Prédiction = Projection orthogonale
pour le produit scalaire de la covariance

Atelier Spatial RASMA 43

Saint Louis du Sénégal
 Le Krigeage universel : m inconnue

• X suit un modèle de régression : covariables Z,

paramètre δ inconnu, ε résidu de covariance Σ

• Le krigeage :
1. Estimer δ par MCG
2. Krigeage simple sur résidu ε

Atelier Spatial RASMA 44

Saint Louis du Sénégal
 Krigeage / prédiction avec geoR
• [Link]
Effectue la prédiction partout (en fait sur une grille à
définir) pour un modèle de variogramme donné (ou un
modèle estimé par variofit ou par likfit)
Option pour le krigeage simple, ordinaire opu universel.

• [Link]
Réalise des simulations conditionnelles; permet
d’évaluer la carte des probabilités de dépasser un seuil.

• Voir également les packages

Fields (Krig et [Link]) et
RandomFields (Kriging et CondSim)
Atelier Spatial RASMA 45
Saint Louis du Sénégal
Atelier Spatial RASMA 46
Saint Louis du Sénégal
 Krigeage des pluies Parana
Modèle : régression affine et résidu Gauss + Pépite

- à gauche carte des hauteurs de pluie (données en ●)

- à droite carte des écarts types des prédictions

Atelier Spatial RASMA 47

Saint Louis du Sénégal
 Propriété de la surface de Krigeage

1. Interpolateur universel : krigeage ≡ observation en

un point d’observation
2. Régularité de la surface de krigeage fonction de la
régularité du variogramme (covariance) en 0 :

• Pépitique : krigeage constant partout ≡ la moyenne arithmétique

des observations, discontinuité aux points d’observation.
• Linéaire pointu en 0 : surface continue mais non dérivable aux
points d’observations.
• Parabolique en 0 : continue et dérivable partout.

Atelier Spatial RASMA 48

Saint Louis du Sénégal
Régularité du Krigeage et régularité à l’origine de
la covariance spatiale

Exemple d’école : 2 observations en (x) et 3 reconstructions

Atelier Spatial RASMA 49

Saint Louis du Sénégal
Atelier Spatial RASMA 50
Saint Louis du Sénégal
 Krigeage ou simulation conditionnelle

(a) X simulé, gaussien, cov.

exponentielle
et
(b) 25 points échantillonnés

(c) reconstruction de X par

krigeage

(d) reconstruction de X par

simulation conditionnelle
aux 25 observations

Atelier Spatial RASMA 51

Saint Louis du Sénégal
Simulation conditionnelle à 100 points

X initial simulation cond.

Atelier Spatial RASMA 52

Saint Louis du Sénégal
Atelier Spatial RASMA 53
Saint Louis du Sénégal
 Nuée variographique (cas isotropique)

• Plus lissage avec bon noyau de convolution

• si X non isotropique, nuées dans les 4 directions cardinales
{E,NE,N,NW} avec une tolérance angulaire de ± 22.5 °

Atelier Spatial RASMA 54

Saint Louis du Sénégal
 Estimation empirique du variogramme
1 - N(h) : ensembles des couples r-voisins à ∆-près
(au moins 30 points dans chaque N(h) )
2 - Avantage : ne nécessite pas d’estimation préalable de la
moyenne (nécessaire pour estimer une covariance)
3 - Version robuste aux grandes valeurs (X(s(i))-X(s(j)))**2
(Cressie et Hawkins)

Atelier Spatial RASMA 55

Saint Louis du Sénégal
Estimation empirique du variogramme :
Choix de classes de distances et effectifs par classes

Atelier Spatial RASMA 56

Saint Louis du Sénégal
Nuée variographique et variogrammes lissés
(données de pluies au Parana)

(a) modèle isotropique (b) dans 4 directions

Atelier Spatial RASMA 57

Saint Louis du Sénégal
 variog : calcul du variogramme empirique

• isotrope ou non (variog4),

• Fixation ou non du nombre de classes (bin)
• estimation classique ou robuste

Retourne une estimation par classe (bin), le nuage

variographique, le variogramme lissé …

> variog(ca20)
> print (x)
> plot(x)

Atelier Spatial RASMA 58

Saint Louis du Sénégal
Données Parana : ajustement affine
Rouges : au dessus Noirs : au dessous

calcul des résidus puis variogramme des résidus

Atelier Spatial RASMA 59

Saint Louis du Sénégal
Variogrammes des résidus dans 4 directions pour:
(c) le modèle moyen affine (pluies au Parana);
(d) le modèle moyen quadratique ……………….

Atelier Spatial RASMA 60

Saint Louis du Sénégal
Estimation d’un modèle paramétrique :

Choisir k classes de distances rendant

identifiable le modèle de variogramme :

Atelier Spatial RASMA 61

Saint Louis du Sénégal
Estimations d’un modèle de variogramme

1 – Moindres carrés ordinaires (MCO)

2 – Moindres carrés pondérés (MCP)
3 – Moindres carrés généralisés (MCG) et MCQG
4 – Maximum de Vraisemblance (MV)

Atelier Spatial RASMA 62

Saint Louis du Sénégal
Estimation du variogramme en présence
d’une tendance linéaire

1. Estimer par MCO

2. En déduire les résidus des MCO
3. Variogramme de ces résidus MCO

Atelier Spatial RASMA 63

Saint Louis du Sénégal
Estimations modèles de pluies en Suisse

Atelier Spatial RASMA 64

Saint Louis du Sénégal
 Estimation du modèle de Matern
(a) MCO, MCP et MV du variogramme de Matern
(b) Estimations empiriques + enveloppes sup et inf à partir de m=40
simulations du modèle de Matern estimé par MCP
(intervalle de confiance à 1- 2/(m+1) = 95 % )

Atelier Spatial RASMA 65

Saint Louis du Sénégal
 variofit et likfit : estimation d’un
variogramme paramétrique

• variofit
estime les paramètres d’un modèle de covariance (variogramme)
par MCO ou par MCP, ceci à partir du variogramme empirique
(variog)

• likfit
réalise l’estimation par maximum de vraisemblance

• Exercice : sur un jeu de données (personnelles), estimer un modèle

paramétrique par MCO et par MV. Représentations simultanées des
estimations empiriques et paramétriques.

Atelier Spatial RASMA 66

Saint Louis du Sénégal
 Krigeage pour le modèle de Matern
(données suisses)
(a) carte des pluies et (b) de leurs écarts types

Atelier Spatial RASMA 67

Saint Louis du Sénégal
Exemple : densité de pêche de la coquille
Saint Jacques dans l’atlantique nord

Atelier Spatial RASMA 68

Saint Louis du Sénégal
Atelier Spatial RASMA 69
Saint Louis du Sénégal
Estimation du modèle « sphérique + pépitique »
( après transformation isotropiante )

Atelier Spatial RASMA 70

Saint Louis du Sénégal
Atelier Spatial RASMA 71
Saint Louis du Sénégal
 Convergence et normalité asymptotique
des estimateurs
• Choix de k classes de distances identifiant le modèle
paramétrique
• Régularité suffisante de la fonctionnelle d’estimation
• Faible dépendance (mélange) du champ géostatistique
• Extension des résultat au cas où il existe une tendance
paramétrique (paramètres d’ordre 1, moyenne) et 2
(variogramme des résidus).

Voir les résultats (Lahiri, Cressie, C. Gaetan – X. G.)

Atelier Spatial RASMA 72

Saint Louis du Sénégal
 Validation d’un modèle paramétrique

• Validation croisée

• Validation par Bootstrap paramétrique

Atelier Spatial RASMA 73

Saint Louis du Sénégal
 Validation croisée d’un variogramme

• Estimer le modèle paramétrique

• Éliminer à tour de rôle une observation x(i) et la prédire par
krigeage avec le modèle estimé
• Calculer l’Ecart Quadratique Normalisé (par la variance de
prédiction) Moyen sur tous les x(i) :

• Si l’EQM proche de 1, valider le modèle (à 95% de sécurité) :

Atelier Spatial RASMA 74

Saint Louis du Sénégal
Exemple : validation du modèle de Matern
(données suisses)

Atelier Spatial RASMA 75

Saint Louis du Sénégal
 Validation par bootstrap paramétrique

• Modèle paramétrique
estimation du paramètre

• m simulations des données sous le modèle estimé
• pour chaque simulation
variogramme empirique
• enveloppes inférieure + supérieure des variogrammes
empiriques

bande de confiance au niveau α = (1- 2/(m+1))
• Si le variogramme empirique des données initiales dans
la bande de confiance, le modèle est validé.

Atelier Spatial RASMA 76

Saint Louis du Sénégal
 Validation du modèle de Matern
(données Parana)

- m = 40 simulations sous le modèle de Matern estimé

- intervalle de confiance à 95% pour le variogramme aux 13 distances
- contient les variogrammes empiriques
modèle de Matern valide

Atelier Spatial RASMA 77

Saint Louis du Sénégal
 Bande de confiance, validation de modèle

• [Link]
Répétition de simulations d’un modèle (donné ou estimé)
aux sites d’observations (modèle gaussien)

variogrammes empiriques pour chaque simulation

bande de confiance pour le variogramme

si variogramme théorique dans la bande
modèle
valide

• xvalid
validation croisée via le krigeage : chaque observation
est comparée à sa prédiction à partir des autres)

Atelier Spatial RASMA 78

Saint Louis du Sénégal
 Une étude épidémiologique : prévalence
du paludisme chez l’enfant
(donnée gambia de geoR, Diggle et altri)

• Expliquer la prévalence du palud (% d’enfants malades)

dans un village s
Covariables :

Existence ou non d’un centre de santé primaire

Un indice de végétation (données satellitaires)

Moustiquaire ou non

Si oui, traitée ou non

• Question : existe-t-il un facteur de risque spatial ?

Atelier Spatial RASMA 79
Saint Louis du Sénégal
 Données malaria en Gambie
1 – localisations s des villages échantillonnés
2 – prévalence prop = y/n (y = nb de malades, n nb d’enfants)
3 – nuage (green,prop) pour indice de végétation green

Atelier Spatial RASMA 80

Saint Louis du Sénégal
Le modèle log-linéaire d’étude

Atelier Spatial RASMA 81

Saint Louis du Sénégal
Atelier Spatial RASMA 82
Saint Louis du Sénégal
Atelier Spatial RASMA 83
Saint Louis du Sénégal
Atelier Spatial RASMA 84
Saint Louis du Sénégal
 Modèle sur un réseau discret
SAR ou Champ de Markov

Deux différences entre ces 2 familles :

(a) Spécificité des états E
(b) linéarité ou non du modèle

• Pour AR spatial : E = R ou R**d + modèle « linéaire »

d’équations Simultanées (SAR gaussien ou non)

• Pour Champ de Markov : E général (i.e. {0,1}, fini, N, R,

R+, R**d etc…) et modèle conditionnel (en général) non
linéaire

Pour E = R, les CAR sont aussi Markov

Atelier Spatial RASMA 85
Saint Louis du Sénégal
(I) Champ de Markov

Atelier Spatial RASMA 86

Saint Louis du Sénégal
Atelier Spatial RASMA 87
Saint Louis du Sénégal
Exemple de réseau : eire, 3 réseaux réguliers,
et sids

Atelier Spatial RASMA 88

Saint Louis du Sénégal
 Champ de Markov sur S fini

• S réseau fini (régulier ou non)

• Loi de X caractérisée par ses conditionnelles
• Loi conditionnelle est « locale »
• Espace d’état E général : binaire, fini, N, R+, R, R**p,
mixte (gris x variété), etc

Atelier Spatial RASMA 89

Saint Louis du Sénégal
 Champ de Gibbs

1. Potentiels réels Φ={Φ(x(A))} définis sur une famille de

parties A de S
2. Φ
énergie U(Φ; Λ) sur Λ conditionnelle ∂Λ
3. S’assurer que exp{U(Φ ; Λ)} intégrable (potentiel
admissible)
4. Le champ de Gibbs est de log-densité conditionnelle à
∂Λ proportionnelle à exp{U(Φ ; Λ)}

Atelier Spatial RASMA 90

Saint Louis du Sénégal
Energie et loi conditionnelles d’un champ
de Gibbs

Soit Λ une partie de S, ∂Λ son voisinage

Atelier Spatial RASMA 91

Saint Louis du Sénégal
 Spécification de Gibbs

• Famille des parties A sur lesquelles sont définit les

potentiels
• Les potentiels Φ(A)
• Vérifier l’admissibilité de exp U(Φ)
• Exemple : famille exponentielle de potentiels de
paramètre θ

Atelier Spatial RASMA 92

Saint Louis du Sénégal
 Exemple de potentiels contribuant à une
loi conditionnelle (loi aux 4 ppv)

• Cliques à un point {x}

• À deux points ●----●

• Potentiels contribuants à π(Λ / ∂Λ)

Atelier Spatial RASMA 93

Saint Louis du Sénégal
 Modèle d’Ising : état E = {-1,+1} et sites
S = {1,2,….,n}**2

• Cliques : singletons et paires de ppv (plus proche voisin)

• Potentiels : Φ(x(i))=α.x(i) et Φ(x(i),x(j))=β.x(i).x(j)

• Energie :

• loi conditionnelle en i :

Atelier Spatial RASMA 94

Saint Louis du Sénégal
 Généralisations du modèle d’Ising

• Etats {0,1} : présence - absence (écologie); sain - malade

(épidémiologie)
• Plus d’états : niveaux de gris, E fini (K variétés)
• Anisotropie : β(H) pour horizontal, β(V) pour vertical
• Non stationnaire : α(i) et β(i,j) suivant les sites i, j
• Élargissement du voisinage ∂i : i.e. aux 8 ppv
• Potentiels de triplets, quadruplets ….
• Simulations par MCMC (échantillonneur de Gibbs)

Atelier Spatial RASMA 95

Saint Louis du Sénégal
 Exemples : 3 textures binaires aux 8 – ppv

- Potentiels au plus de paires
5 paramètres

- Simulation par échantillonneur de Gibbs (AntsInfields)
- Arrêt après 3000 itérations

Atelier Spatial RASMA 96

Saint Louis du Sénégal
 Modèle de Gibbs à nombre d’états fini
(modèle de Potts)

• Etats : E = {a(1),a(2),…,a(K)}, K états

• Potentiels : singletons et paires de sites voisins

• Energie : n(k) = nb sites modalité k; n(k,l) = nb de sites

voisins de modalités (k,l).

Atelier Spatial RASMA 97

Saint Louis du Sénégal
 Modèle de Potts échangeable

Très utile en traitement d’image

• Tous les états ont un comportement analogue

c’est-à-dire : α(k) = α et β(k,l) = β pour tout k,l

• Dans ce cas il y a un seul paramètre de dépendance

spatiale β et la loi jointe vaut :

Atelier Spatial RASMA 98

Saint Louis du Sénégal
 Modèle échangeable à 3 états

β ↑ augmente la régularité géométrique des plages constantes

(a) β = 0.5 (b) β = 0.6 (c) β = 0.7
(Echantillonneur de Gibbs à 5000 itérations avec AntsInFields)

Atelier Spatial RASMA 99

Saint Louis du Sénégal
Un champ gaussien est un champ de Gibbs

- Soit X = N(µ , ∑) de moyenne µ et covariance ∑ (Q inverse ∑)

- Alors X est de Gibbs avec :

potentiels de singletons et de paires

énergie – U
- Obtention des lois conditionnelles comme champ de Gibbs

Atelier Spatial RASMA 100

Saint Louis du Sénégal
 Simulation
d’un champ de Markov

Atelier Spatial RASMA 101

Saint Louis du Sénégal
 Simulation par dynamique de chaîne de
Markov
(MCMC pour Monte Carlo Markov Chain)

• Objectif : simuler une loi Π sur E

• Principe MCMC : construire une chaîne de Markov (X(n),
n>0) sur E de transition P(x,•) t.q. :
1 – P est irréductible (tous les états communiquent),
2 – P est Π - invariante (ΠP= Π)
3 – P est apériodique

Propriété : Si (1-2-3), la loi de X(n) tend vers Π

Atelier Spatial RASMA 102

Saint Louis du Sénégal
 Irréductible, invariance et apériodicité

• Irréductible : la chaîne P fait communiquer tous les états de E entre

eux

• Π - invariante (ΠP = Π) : si la loi de X(n) est Π, celle de X(n+1) est

encore Π (Π est la loi invariante de la chaîne)

• Périodique : si il existe une partition {E(1),E(2),…,E(k)} de E t.q. la

chaîne circule successivement dans E(1)
E(2)
…E(k)
E(1) etc

• Apériodique : si non périodique

Atelier Spatial RASMA 103

Saint Louis du Sénégal
 Comment construire une telle chaîne ?

• Apériodicité et irréductibilité se vérifient au cas par cas

• Pas facile de trouver P qui soit Π - invariante : Π est

vecteur propre de P associée à la valeur propre 1 !

• Condition suffisante assurant la Π – invariance de P :

la Π – réversibilité :
pour tout x, y : Π(x)P(x,y) = Π(y)P(y,x)

Atelier Spatial RASMA 104

Saint Louis du Sénégal
Deux algorithmes markoviens classiques

• Echantillonneur de Gibbs sur un espace produit E**S :

S={1,2,…,n} ensemble des sites, état E en chaque site

• L’algorithme de Metropolis (espace E général)

Pour l’un et l’autre, la transition P, est, par construction, Π –

réversible, donc Π – invariante

simulateur MCMC de Π si on vérifie de plus que P est

apériodique et irréductible

Atelier Spatial RASMA 105

Saint Louis du Sénégal
 Simulation par échantillonneur de Gibbs
d’une loi π sur E**S
• Loi π sur espace produit E**S
• Connaître les lois conditionnelles en tout i :

• Relaxation « site par site» sur S suivant la conditionnelle :

• Un « balayage» de S : itérer de i = 1, n
la transition P

• Itérer les balayages
• Propriété : P est π – réversible.
Donc l’échantillonneur de Gibbs simule approximativement π
Après un grand nombre de balayages de S.
Atelier Spatial RASMA 106
Saint Louis du Sénégal
 Transition pour un balayage de S

• Visite séquentielle de S : 1
2
3
….
n

• Au i-ème pas, relaxation au site i

• Enchaînement sur un balayage donne la transition :

Atelier Spatial RASMA 107

Saint Louis du Sénégal
 La transition P de l’algorithme de
Métropolis
• Π loi sur espace d’état E général, x, y deux états

Construction de P en deux étapes :

1 – Proposition de changement x
y suivant une transition

Q(x,y) symétrique (Q est la proposition de changement)
2 – Acceptation du changement avec une probabilité a(x,y)

Propriété : si a(x,y) = min {1,Π(y)/ Π(x)}, alors P est Π –

réversible (donc Π - invariante)

Si de plus P est irréductibilité et apériodique, on a un autre

algorithme de simulation de Π
Atelier Spatial RASMA 108
Saint Louis du Sénégal
 Algorithme de Métropolis
1 - x état initial. Changement x
y suivant Q(x,y)
2 - Si Π(y) ≥ Π(x), garder y.
3 – Sinon, tirer U uniforme sur [0,1] :
(a) si U > p = Π(y)/Π(x), garder x.
(b) Sinon, garder y.
4 – Revenir en (1)

L’algorithme de Métropolis – Hastings correspond à une

proposition de changement Q non symétrique.

Remarque importante : il suffit de connaître Π à un facteur

près pour construire cet algorithme (le cas pour Π de Gibbs).

Atelier Spatial RASMA 109

Saint Louis du Sénégal
Exemple : simulation d’un Ising isotropique
aux 4-ppv

(I) - Par échantillonneur de Gibbs

Les lois conditionnelles en chaque site i sont explicites
en terme de v(i), la somme aux 4-ppv

Atelier Spatial RASMA 110

Saint Louis du Sénégal
 AntsInFields

« Boite noire » illustrant le livre de G. Winkler (Springer, 2002) :

« Image Analysis, random fields and dynamic MC Methods »

Largement buggée et non ouverte à la programmation.

Outil de démonstration, sur des thèmes d’analyse d’image et de statistique de
champ de Gibbs.

Exemple : § 3 – 2

simulation Ising isotropique aux 4-ppv par échantillonneur de Gibbs
h=0, b=0, 0.2 et 0.4. Voir petit à petit se former des plages

simulation d’un modèle de Potts à 4 niveaux de gris, paramètres :

nombres de classes, h et b et distance retenue entre les configurations
voisines

Atelier Spatial RASMA 111

Saint Louis du Sénégal
 Simulation d’un Ising (suite)

(II) - Métropolis par échange de spins

1 - on tire au hasard 2 sites i et j et on permute les spins x(i)

et x(j) : on passe ainsi de x
y avec une probabilité de
transition Q(x,y)
2 – le quotient Π(x)/Π(y) s’explicite facilement en fonction de
x(i), x(j), v(i) et v(j) (cf. poly)
3 – Mettre en œuvre Métropolis.

Atelier Spatial RASMA 112

Saint Louis du Sénégal
 Champ de Markov et champ de Gibbs

• S={1,2;3, …., n} et G graphe symétrique sur S

• <i,j> : i et j voisins pour G.
• ∂A = voisinage de A pour G
Clique de G : singletons + parties A t.q. les points de A
sont tous voisins
• C(G) = toutes les cliques de G
Champ de Markov :

Atelier Spatial RASMA 113

Saint Louis du Sénégal
 Le théorème de Hammersley – Clifford
Gibbs ≡ Markov

• π un champ de Markov pour graphe G

• Positivité : pour tout x, π(x)>0

• Propriété (H-C) : alors π est- un champ de Gibbs dont

les potentiels sont limités aux cliques de G

• Réciproque : tout champ de Gibbs est un champ de

Markov pour le graphe engendré par les potentiels de
Gibbs.

Atelier Spatial RASMA 114

Saint Louis du Sénégal
 Recollement de lois conditionnelles
• Objectif : définir un modèle à partir de ses spécifications locales
• En général, des spécifications locales « ne se recollent pas »
• E sous ensemble réel
• Les spécifications de la « famille exponentielle » ci-dessous se
recollent.
Résultat : si pour tout i

Atelier Spatial RASMA 115

Saint Louis du Sénégal
 Auto – modèle de Besag (1974)

• E sous ensemble de R
• Famille exponentielle du type précédent
Alors les lois conditionnelles se recollent en la loi de
Gibbs π :

Atelier Spatial RASMA 116

Saint Louis du Sénégal
 Auto-modèle de Markov :
auto - régression pour espace E général

• E = {0,1} : auto-logistique (états binaires)

• E = {0,1,2, …, K} : auto-binomial
• E = N : auto-poisson (comptage en épidémiologie)
• E = R+ : auto-exponentiel (Gamma), pluviométrie
• E = R (R**d) : auto-gaussien
• Possibilité d’ajouter des covariables explicatives

Atelier Spatial RASMA 117

Saint Louis du Sénégal
 Modèle Auto-Logistique : E = {0,1}

• paramètre θ(x(i)) et loi Logit conditionnel

• L’énergie jointe du champ de Gibbs est

• Généralisation à l’auto-binomial
Atelier Spatial RASMA 118
Saint Louis du Sénégal
Auto – modèle de Poisson : E = N
(comptage en épidémiologie, etc)

- la loi conditionnelle Poisson suit un MLG

- admissibilité : β <0 (compétition)
- coopération possible en bornant E à K<∞

Atelier Spatial RASMA 119

Saint Louis du Sénégal
 Auto – exponentiel : E = R+
(variable >0 : pluviométrie, etc)

- loi conditionnelle exponentielle suit un MLG

- admissibilité : β <0 (compétition)
- coopération possible en bornant E à à K<∞

Atelier Spatial RASMA 120

Saint Louis du Sénégal
 Auto – modèle avec covariables z

• Il y a trop de paramètres α et β si le modèle est non

stationnaire
• Modéliser les α, β à partir de covariables z :
Exemple :

Atelier Spatial RASMA 121

Saint Louis du Sénégal
 Estimation d’un champ de Markov

3 procédures

• Max de Vraisemblance : efficace, difficile (cste de

normalisation)
méthodes numériques MCMC

• Pseudo Vraisemblance Conditionnelle (PVC) : facile à

mettre en place, bonnes propriétés, proche MV si peu de
dépendance spatiale.

• Codage : facile, moins efficace, test de chi2 direct

Atelier Spatial RASMA 122

Saint Louis du Sénégal
 Maximum de vraisemblance sur D(n)

• L’estimation du MV est convergente si π appartient à

famille exponentielle invariante par translation (S = Z**2)
• Difficulté de calcul de la constante de normalisation Z

calcul de Z par MCMC ou par algorithme de score

Atelier Spatial RASMA 123

Saint Louis du Sénégal
 Pseudo – Vraisemblance Conditionnelle
(PVC – Besag : 1974)

• PVC = produit en chaque site i des probabilités conditionnelles

• Bonne fonctionnelle d’estimation (convergence, normalité sous des
hypothèses raisonnables)
• Pour un CAR gaussien ≡ MCO sur les résidus
• Obtention estimation via logiciel dédié aux MLG ……
• Attention : calcul spécifique de la variance d’estimation

Atelier Spatial RASMA 124

Saint Louis du Sénégal
Influence du taux de nitrate des eaux sur la
mortalité par cancer, Valence – Espagne
(Ferrandiz et al. Biometrics – 1995)

Atelier Spatial RASMA 125

Saint Louis du Sénégal
Influence du taux de nitrate dans les eaux sur
la mortalité par cancer (suite)

Atelier Spatial RASMA 126

Saint Louis du Sénégal
 Estimation et simulation de 3 textures binaires:
(a) cailloux; (b) liège, (c) rideau (Cross et Jain)

• 3 textures binaires réelles de

taille 64 x 64

Estimation par PVC de 3

modèle auto-logistique

• Simulation des textures

estimées (éch. de Gibbs)

• Bilan : bonne adaptation de:

(1) la modélisation Markov,
(2) l’estimation par PVC et
(3) la simulation par ech. de
Gibbs.
Atelier Spatial RASMA 127
Saint Louis du Sénégal
 C « Codage » de S

Définition : C est un codage de S si 2

sites s ≠ t de C ne sont jamais voisins

• Ex 1 : les x rouges pour la relation aux 4 ppv

• Ex 2 : les □ verts pour la relation aux 8 ppv

Atelier Spatial RASMA 128

Saint Louis du Sénégal
 Estimation par C - codage

• Soit C un ensemble de codage de S

• La propriété fondamentale : indépendance des X(s), s dans C,
conditionnellement aux x(S\C) extérieurs
• La vraisemblance sur C conditionnelle aux x(S\C) est (exactement)
le produit des lois conditionnelles

Atelier Spatial RASMA 129

Saint Louis du Sénégal
 Conséquences

• Même propriétés que pour un estimateur du MV de

variables indépendantes (i.n.i.d.)
• Normalité, test du chi2 d’une sous hypothèse
• Calcul de l’estimation et de la variance d’estimation avec
un logiciel dédié aux MLG

• Plusieurs choix d’ensemble de codage possible, mais les

estimateurs associés sont dépendants

Atelier Spatial RASMA 130

Saint Louis du Sénégal
Modèle de répartition spatiale d’une
espèce végétale (présence / absence)
Variété : la grande laîche

Atelier Spatial RASMA 131

Saint Louis du Sénégal
 2 modèles auto - logistiques :
(I) aux 4 – ppv (v(1)) - (II) aux 8 – ppv (v(2))

Atelier Spatial RASMA 132

Saint Louis du Sénégal
Résultats : estimations MV, PMV et codage

Atelier Spatial RASMA 133

Saint Louis du Sénégal
Atelier Spatial RASMA 134
Saint Louis du Sénégal
(II) – Modèle AR spatial

Atelier Spatial RASMA 135

Saint Louis du Sénégal
 (II) - Modèle AR spatial

• L’espace d’état est E = R ou R**p

• Le modèle est définit par un ensemble d’équations

« spatiales » simultanées (comme en économétrie) en
référence à un graphe d’influence

• Le bruit de modèle est blanc (SAR) ou coloré (CAR)

• Souvent le modèle est supposé gaussien

Atelier Spatial RASMA 136

Saint Louis du Sénégal
Ex 1 : la mort subite du nourisson
(données sids de spdep, Cressie et al)

Atelier Spatial RASMA 137

Saint Louis du Sénégal
 Ex. de système de voisinage

à gauche : cantons avec frontière commune

à droite : centres sont à moins de 30 miles

Atelier Spatial RASMA 138

Saint Louis du Sénégal
Atelier Spatial RASMA 139
Saint Louis du Sénégal
 [Link] dans le package spdep

• Voir le descriptif de [Link]

• Représentation des 100 comtés et des 21
variables
• Graphes de voisinage

Atelier Spatial RASMA 140

Saint Louis du Sénégal
Ex. 2 : Taux groupe sanguin A
36 contés de l’Irlande (données eire)

Atelier Spatial RASMA 141

Saint Louis du Sénégal
Graphe de voisinage, données eire
« avoir de la frontière commune »

Atelier Spatial RASMA 142

Saint Louis du Sénégal
 Réseau régulier (agronomie, télédétection) ou non
Choix du graphe de voisinage reste à faire
Exemples de graphes réguliers symétriques :
(a) Lattice carré : voisinage aux 4 plus proches voisins (p.p.v.)

(b) Lattice carré : voisinage aux 8 - ppv

(c) Lattice triangulaire : voisinage aux 6 - ppv

Atelier Spatial RASMA 143

Saint Louis du Sénégal
 SAR général sur S = {1,2, …., n}
• recentrage de X
• Graphe et poids d’influence W + Bruit
• n équations simultanées de paramètres A avec ε BB
• X existe si A inversible
• Covariance ∑ en termes des paramètres A

Atelier Spatial RASMA 144

Saint Louis du Sénégal
SAR général : graphe et poids

Atelier Spatial RASMA 145

Saint Louis du Sénégal
 Spécifier une SAR

• ∂s = voisins de s (symétrique ou non)

• Dépendance « locale » : X(s) = F ( X(∂s),θ) + ε(s)
• F linéaire via θ inconnu et W matrice de poids connus
W={W(t,s), t voisins de s}
• ε un BB (éventuellement gaussien)

Exemple : ρ = corrélation spatiale, (I- ρW) inversible

Atelier Spatial RASMA 146

Saint Louis du Sénégal
Atelier Spatial RASMA 147
Saint Louis du Sénégal
 Choix ad hoc de W

• Fonction des distances inter-centres, des (portions) de

frontières communes, des réseaux de communications
entre deux cellules, etc…
• Paramètres γ et τ préalablement calibrés

Atelier Spatial RASMA 148

Saint Louis du Sénégal
 SAR stationnaire sur Z**2

• bruit blanc : (gaussien ou non)

• Équations avec variables «spatialement retardées»

• X existe si P ≠ 0 sur le tore :

• Graphe R de voisinage orienté (ou non)

Atelier Spatial RASMA 149
Saint Louis du Sénégal
Exemples de SAR

Atelier Spatial RASMA 150

Saint Louis du Sénégal
 Précautions sur un SAR

• Un SAR bilatéral ≠> AR causal (pour l’ordre

lexicographique; cf. exemple dans le polycopié)

• Sans contrainte, un SAR non identifiable

• L’estimation des MCO est non convergente

• Avantage : SAR est parcimonieux en paramètres

Atelier Spatial RASMA 151

Saint Louis du Sénégal
 AR conditionnel général (CAR) sur S
• Écrire l’espérance conditionnelle de X(t) sur autres X :
Les résidus e sont corrélés entre eux, décorrélés des X :

• Notations : D diagonale des résidus, ∑ = Cov(X)

C paramètres CAR. On a l’identité :

 contraintes sur les paramètres du CAR

Atelier Spatial RASMA 152
Saint Louis du Sénégal
Contraintes sur les paramètres d’un CAR

Atelier Spatial RASMA 153

Saint Louis du Sénégal
 AR Conditionnelle (CAR) stationnaire

• L’espérance conditionnelle linéaire de X(s) sur les autres X(t)

• e(t) est décorrélé des X(s) pour s ≠ t
• Le graphe d’un CAR symétrique ainsi que les c.
• Le résidu conditionnel e est un bruit coloré (c-à-d corrélé) :

Atelier Spatial RASMA 154

Saint Louis du Sénégal
CAR : espérance conditionnelle linéaire
et résidu

Atelier Spatial RASMA 155

Saint Louis du Sénégal
 SAR ou CAR : pour résumer

• Un SAR est spécifié par n équations simultanées

à résidus BB

• Un CAR est spécifié par ses « espérances

conditionnelles linéaires ». Le bruit résiduel est
coloré

Atelier Spatial RASMA 156

Saint Louis du Sénégal
 CAR ou SAR ?

• Tout SAR est un CAR

• Si S fini, CAR ≡ SAR ( ≠ sinon)
• Écriture CAR intrinsèque, celle d’un SAR non
• Estimation MCO
d’un CAR convergente
d’un SAR non
• SAR : plus parcimonieux en nombre de paramètres
• CAR : contraintes sur les paramètres

Atelier Spatial RASMA 157

Saint Louis du Sénégal
 Correspondances des graphes :
R d’un SAR  et
G d’un CAR

• R graphe du SAR : orienté

• G du CAR : non orienté, le « double » de R :

Atelier Spatial RASMA 158

Saint Louis du Sénégal
Exemple de correspondance
(a) SAR et R (b) CAR associé et G
en pointillé : voisinage CAR de ●

Atelier Spatial RASMA 159

Saint Louis du Sénégal
 SAR aux 4 ppv et CAR associé

• SAR aux 4-ppv avec R = {(1,0),(-1,0),(0,1),(0,-1)}

CAR aux 12-ppv (cf. poly. pour les coeff. c(s)) avec
L+ = {(1,0),(2,0),(1,1),(0,1),(0,2),(1,-1)}

avec un gain de prédiction

κ**2 =

Atelier Spatial RASMA 160

Saint Louis du Sénégal
SAR factorisant aux 3 – ppv
et CAR aux 8-ppv associé
R = support du SAR et L = support du CAR (8 voisins)

Atelier Spatial RASMA 161

Saint Louis du Sénégal
 SARX avec exogènes
• X endogène, Z matrice des exogènes
• Matrice de retard W sur endogène et exogène

• 3 types de variables expliquent X(t) :

(1) endogène retardée WX,
(2) exogène Z et
(3) exogène retardée WZ.

• On obtient facilement E(X) et Cov(X)

Atelier Spatial RASMA 162

Saint Louis du Sénégal
 Deux modèles avec exogènes

• Modèle de Durbin spatial (X-Zβ ~ SAR(ρ,W) )

• Modèle à décalage spatial (γ = 0)

Atelier Spatial RASMA 163

Saint Louis du Sénégal
 Auto - corrélation de Moran

• X sur S={1,2,…,n} centré (adaptation si modèle de

régression sur E(X))
• W(i,j) matrice de poids i
j, i ≠ j (W(i,i)=0) donnée
• W – auto-corrélation de Moran :

Atelier Spatial RASMA 164

Saint Louis du Sénégal
 Test de non corrélation spatiale (H(0))

• I(M) petit de variance identifiée

• Si X gaussien, plus précisions sur E() et Var().

• En général, sous (H0), normalité :

Atelier Spatial RASMA 165

Saint Louis du Sénégal
 Indice de Geary

Mesure la dépendance spatiale comme le fait un variogramme :

I(G) est petit si les valeurs voisines sont proches

Atelier Spatial RASMA 166

Saint Louis du Sénégal
 Loi permutationnelle d’une statistique I(X)

• X={X(i),i=1,n} et I(X) une statistique réelle

• Distribution empirique des {I(X(σ)), σ permutation}
• Intervalle de confiance associé à la statistique d’ordre
• Mais n! est trop grand
le faire pour m permutations
choisies au hasard

Atelier Spatial RASMA 167

Saint Louis du Sénégal
 Application : test de permutation de
(H0) : indépendance des {X(i),i=1,n}
• Choix de I(X), l’indice de Moran : sous (H0),

• Calcul des I(σ,X) pour m permutations au hasard (i.e.

m=1000) et de l’intervalle de confiance empirique IC(1-α)
• Calcul de I(x) pour l’observation x
• Si I(x) n’est pas dans IC(1-α), rejet de (H0)

• Avantage : non – asymptotique, libre du modèle sur X

• Inconvénient : le niveau est approximatif

Atelier Spatial RASMA 168

Saint Louis du Sénégal
 Données eire : groupe sanguin A

1 - G = graphe de voisinage de contiguité des 26 contés

2 - w(i,j) = 1/(nb voisins de i) si j est voisin de i

Atelier Spatial RASMA 169

Saint Louis du Sénégal
 Indice de Moran (Geary) + indices réduits t(a)
probabilités p(a) de dépassement

1 – I = index, t(a) index réduit, asymptotique gaussienne

2 – p(a) = proba de dépassement de t(a)
3 – pmc(a) pour le test de permutation MC (m = 1000
permutations)

Atelier Spatial RASMA 170

Saint Louis du Sénégal
 Le package spdep :
étude sur l’exemple des données eire

• Données i.e. : eire groupe sanguin en Irlande

• Tester bloc par bloc l’exemple de traitement de ces
données
1 - représentation des données (valeur, graphe de
voisinage)
2 - indice de Moran et test de non corrélation spatiale
3 - régression sur town et pale, analyse des résidus
(indice de Moran, SAR sur les résidus)
• Autres données sur la consommation intérieure

Atelier Spatial RASMA 171

Saint Louis du Sénégal
 Quelques programmes de spdep

• moran, [Link] et [Link]

• lagsarlm :estimation du MV d’un SAR avec covariables
y = rho W y + X beta + e
• knearneigh : matrice des k-ppv pour un choix de distance
• [Link] : test de Moran pour l’auto-corrélation
spatiale des résidus d’un modèle linéaire
• [Link] : test exact de non corrélation
• [Link] : corrélogramme spatial pour l’indice de
moran
et d’autres programmes ….

Atelier Spatial RASMA 172

Saint Louis du Sénégal
 Estimation d’une régression spatiale

X = Zδ+ε où Cov(ε)=Σ

• Estimation MCO de δ :

• Sous bonnes conditions, consistance des MCO

MCO bonne estimation initiale dans une procédure

itérative type MCQG

Atelier Spatial RASMA 173

Saint Louis du Sénégal
 Moindres Carrés Généralisés (MCG)

• Si Σ = cov (ε) est connue, le BLUE vaut

• Si X gaussien, c’est l’EMV, efficace

• En général Σ inconnue
MC Quasi G (MCQG)

Atelier Spatial RASMA 174

Saint Louis du Sénégal
MCQG : Σ = Σ(θ), θ inconnu

Atelier Spatial RASMA 175

Saint Louis du Sénégal
 Régression Gaussienne : MV
• Régression à covariance non sphérique :

• Mardia-Marshall donnent le comportement limite de

l’EMV de (θ,∂) (cf. poly)

• Log-vraisemblance est explicite :

Atelier Spatial RASMA 176

Saint Louis du Sénégal
 Données eire : 2 modèles de régression
avec 2 covariables
- towns (densité urbaine) et
- pale (binaire, 1 si colonisation anglaise, 0 sinon)
(R1) : cste, towns et pale + résidus i.i.d.
(R2) : cste, pale + résidus SAR aux ppv : r = (ρ W) r + e

Atelier Spatial RASMA 177

Saint Louis du Sénégal
Atelier Spatial RASMA 178
Saint Louis du Sénégal
Atelier Spatial RASMA 179
Saint Louis du Sénégal
 Répartition de 97 fourmilières
(données ants du package spatstat)

Question : la répartition s’est elle faite au hasard ?

Atelier Spatial RASMA 180

Saint Louis du Sénégal
42 centres de cellules d’une coupe histologique
(données cells de spatstat)
1 - La répartition est – elle au hasard ?
2 - Sinon (+ de régularité), quel modèle proposer ?

Atelier Spatial RASMA 181

Saint Louis du Sénégal
 126 pins d’une forêt finlandaise
(données finpines de spatstat)

1 - La répartition est elle au hasard ?

2 - Sinon (des aglomérats ?), quel modèle proposer ?

Atelier Spatial RASMA 182

Saint Louis du Sénégal
 Modèle de Processus Ponctuel X (PP)

• Configuration x : ensemble fini de points de la fenêtre

d’observation S
• Configuration à n points : x = {x(1),x(2),…,x(n)}

• E(n) = espace des configurations à n points

• E = U E(n) : l’espace exponentiel de toutes les
configurations, réunion des E(n)

• N(A) : le nombre de points de X dans A

• Loi de X : loi jointe de toutes les variables de comptage
N(A), A partie de S

Atelier Spatial RASMA 183

Saint Louis du Sénégal
 PP de Poisson homogène d’intensité λ

PPP(λ) : répartition spatiale homogène et au hasard

(1) N(A) suit une loi de Poisson de paramètre λ IAI

(2) La répartition sur A est uniforme

(1-2) ≡ (1-2*) où
(2*) : si A et B sont disjoints, N(A) et N(B) sont indépendants

Atelier Spatial RASMA 184

Saint Louis du Sénégal
Deux répartitions spatiales homogènes à 70 points
• (a) une répartition de Poisson homogène
• (b) une répartition « à r - noyau dur» homogène
(couples de points à distance < r = 3.5 sont interdits)

Atelier Spatial RASMA 185

Saint Louis du Sénégal
 PPP inhomogène d’intensité λ(●)

Soit λ(●) une mesure sur la fenêtre d’observation S

X est un PPP(λ(●)) si :
(1) N(A) suit une loi de Poisson de paramètre λ(A)
(2) si A et B sont disjoints, N(A) et N(B) sont
indépendants

Atelier Spatial RASMA 186

Saint Louis du Sénégal
 Simulation d’un PPP(λ(●))

Supposons que pour tout x : λ(x) ≤ c < ∞

La méthode par effacement de points est la méthode de
simulationpar rejet suivante :

1. Simuler x*, un PPP homogène d’intensité c;

2. Effacer indépendamment un x(i) de x* avec la
probabilité p(x(i)) = {1- λ(x(i))/c} .

Atelier Spatial RASMA 187

Saint Louis du Sénégal
Simulation de 2 PPP inhomogènes
(λ plus grand
fond plus clair)

Atelier Spatial RASMA 188

Saint Louis du Sénégal
 Répartition plus régulière : modèle à
noyau dur (ou hardcore)

La règle : interdire les points trop proches

• Exemples :
- répartition spatiale d’animaux (compétition)
- boulangeries dans une ville
- centres de cellules
- arbres dans une forêt (??)
- en physique, centres d’«atomes impénétrables»
• Ces modèles vont être défini par leur densité de Gibbs

Atelier Spatial RASMA 189

Saint Louis du Sénégal
 Répartition moins régulière :
formation d’agrégats (clusters)
Exemple : le PP de Neymann – Scott

1. Un processus P « parent » : PPP homogène λ

2. Chaque parents P(i) engendre des enfants en nombre N
et en positions D centrées autour de P(i), aléatoires, N
et D indépendantes

Paramètres : λ, les lois N et D

Atelier Spatial RASMA 190

Saint Louis du Sénégal
 Simulation de PP spatiaux avec statspat

• owin : crée la fenêtre d’observation (si nécessaire)

• runifpoint : n points uniformes.
• runifpoint3 : idem mais dans R**3 (installer le package
scatterplot3d pour la représentation 3d).
> X = runifpoint3(5000)
> plot(X)

• rpoispp : simulation d’un PPP (homogène ou non)

• rNeymanScott : PP de N-S avec agrégat
• rThomas ….

• rmh : simulation d’un PP à partir de son modèle de densité

(Strauss, noyau dur, etc)

Atelier Spatial RASMA 191

Saint Louis du Sénégal
 Ajustement finpines sur un Neymann –
Scott
• (a) Données réelles : modèle de NS à 3 paramètres θ = (λ,µ,σ**2 )
1 - parents Poisson λ
2 - nombre de descendants d’un père Poisson µ
3 - répartition des fils autour d’un père Gaussienne sphérique σ**2

Ajustement par MCO (cf. poly) puis
• (b) Simulation du NS estimé (∆ parents et ● descendants)

Atelier Spatial RASMA 192

Saint Louis du Sénégal
 PP doublement Poissonien

• PPP à intensité aléatoire {Λ(s), s dans S}

• Exemple : PP de Cox log-gaussien X
- Λ suit le modèle log-linéaire à effet aléatoire

- ψ un champ Gaussien centré de covariance c

- c contrôle la corrélation spatiale de X
(Moller – Waagepertersen)

Atelier Spatial RASMA 193

Saint Louis du Sénégal
 Deux exemples de PP de Cox
- intensité Λ en fond grisé (Λ(s) grand, fond clair)
- Modèle log-Linéaire : β = z ≡ 1 partout
- Deux covariances c pour l’intensité aléatoire Λ

Atelier Spatial RASMA 194

Saint Louis du Sénégal
 PP marqué (PPM)

• Une marque m(x(i)) s’ajoute en chaque x(i) de X

• Exemples :
- marque « disque » centré en x(i) (diamètre arbre)
- rayon r du disque est fixé ou aléatoire R
- marques fibres curvilignes attachées à x(i) (système
racinaire d’une plante, segment de longueur et
orientation aléatoires)
- nombre fini K de marques : K états possibles, une
couleur est associée à chaque x(i) (i.e. K=2 deux états
« sain » ou « malade »)

Atelier Spatial RASMA 195

Saint Louis du Sénégal
 Deux exemples de PPM
(a) Marques binaires : localisation des 62 cas (●) de Leucémie d’un
canton et de 141 résidences (x) d’enfants sains (K=2, données
humberside)

(b) Marques continues : positions et tailles des 584 aiguilles de pin

d’un sous bois (données longleaf)
Question sur (a) : effet spatial influençant la maladie ?

Atelier Spatial RASMA 196

Saint Louis du Sénégal
 Densité f d’un PP

• f : densité de probabilité par rapport à un PPP(1)

• f : E
R, E = espace exponentiel de toute les
configurations x
• f admissible si intégrable, d’intégrale 1

• En général on définit f à une constante près :

f(x) = c g(x) où g explicite (mais pas c!)

• Inutile connaître c pour la simulation (Metropolis)

• Mais il faut connaître c = c(θ) pour l’estimation du MV de θ

Atelier Spatial RASMA 197

Saint Louis du Sénégal
 Exemple de PP à densité : PP de Gibbs

• U(x)  potentiel de Gibbs : U(x) = Σ Φ(A)(x)

• Admissibilité de exp {U(x)}

• Exemple : famille exponentielle

f(x)=c(θ) exp { < θ , T(x) > }

Atelier Spatial RASMA 198

Saint Louis du Sénégal
 PP de Strauss
• U(x) dérive de 2 statistiques issues de x :
1 - n(x) = nombre de points de x
2 - s(x) = nombre de couples de x à distances < r

• Energie : U(x) = a n(x) + b s(x) (ou a = log β et b = log γ)

β (ou a) règle l’intensité de x; γ règle la régularité spatiale :

- γ < 1 : d’autant plus régulière que γ petit
- γ = 1 : PP de Poisson homogène d’intensité b=log(γ)
- γ > 1 : formation d’agrégats

• PP à noyau dur : γ = 0, interdit les couples à distance < r

Atelier Spatial RASMA 199

Saint Louis du Sénégal
 Simulation Metropolis d’un PP de Gibbs
On circule dans les espaces E(n) en autorisant à une itération :

soit une naissance (proba ½)

soit une mort (proba ½)
suivant la règle suivante :

Atelier Spatial RASMA 200

Saint Louis du Sénégal
 Quelques outils statistiques

• Moments d’ordre 1 ou intensité

(modèle sur la moyenne)

• Moments d’ordre 2, corrélation repondérée

(indépendance spatiale ou non)

• Moment réduit K d’ordre 2 de Ripley

• Distances aux plus proches voisins

Atelier Spatial RASMA 201

Saint Louis du Sénégal
Moments (intensités) d’ordre 1 et 2 d’un PP

Atelier Spatial RASMA 202

Saint Louis du Sénégal
Moment K de Ripley (cas isotropique)

Atelier Spatial RASMA 203

Saint Louis du Sénégal
K et régularité spatiale ?

Atelier Spatial RASMA 204

Saint Louis du Sénégal
 Distances aux plus proche voisins (ppv)

1 - d’un point ● de X
2 – d’un point o de la fenêtre d’observation

Atelier Spatial RASMA 205

Saint Louis du Sénégal
Distances aux ppv et régularité spatiale ?

Atelier Spatial RASMA 206

Saint Louis du Sénégal
 Questions de base pour un PP

• La répartition spatiale X est-elle due uniquement au hasard ?

(CSR = Complete Spatial Randomness) ou non ?

• Si non, y a t’il compétition ? coopération ?

• Quel modèle de PP spatial proposer pour X ?

• Comment estimer le modèle ? Comment le valider ?

Méthodes paramétriques

ou méthodes de Monte Carlo

Atelier Spatial RASMA 207

Saint Louis du Sénégal
 Test de « CSR » : X est un PPP

Utilisation des distances aux PPV et fonctionnelle J

Estimation de J  estimations empiriques de G et F

Atelier Spatial RASMA 208

Saint Louis du Sénégal
Test de « CSR » : utilisation de K

Atelier Spatial RASMA 209

Saint Louis du Sénégal
 Test de « CSR » : comptages par quadrats
Chi 2 d’ajustement sur une loi de Poisson (finpines)
• On forme (par exemple) 4 x 4 = 16 quadrats même surface
• Comptage des effectifs N(i) pour chaque quadrat i
• Distribution empirique des N(i) (à droite)
• Distance à une distribution de Poisson : D = 46.7 >> χ²(15;5%)

Atelier Spatial RASMA 210

Saint Louis du Sénégal
 Estimation de l’intensité du PP

• Méthodes non paramétriques classiques : fenêtrage

et noyau de convolution

• Modèle paramétrique : définir le modèle et estimer

par MV en supposant CSR (c’est une PV si CSR non
vérifiée)

bonne propriétés asymptotiques si X est ergodique

Atelier Spatial RASMA 211

Saint Louis du Sénégal
 Estimation NP d’une densité d’un PP
(a) données érables (lansing) et
3 estimations (b-c-d) avec noyaux de + en + régularisant

Atelier Spatial RASMA 212

Saint Louis du Sénégal
Estimations NP des distances aux ppv et K
1 – à partir de statistiques empiriques de comptage
2 -- lissage postérieur éventuel

Moment K de Ripley : sur le disque B centré en 0 de rayon

h si X est de densité τ :

Distribution G aux ppv : x un point de X, h(i) les n

distances d(x,x(i)) pour les n points de X :

Atelier Spatial RASMA 213

Saint Louis du Sénégal
 Estimation d’un modèle d’intensité ρ(●,θ)

• Si X est PPP d’intensité ρ(●,θ), la log-vraisemblance de

{x(1),x(2),…., x(n)} sur A est donnée ci-dessous
si ρ(●,θ) suit un MLG, maximisation via un logiciel dédié

• Si X est PP de densité ρ(●), on maximise encore cette

Pseudo-Vraisemblance
Bonne propriété limite si X ergodique et A → R**2.

Atelier Spatial RASMA 214

Saint Louis du Sénégal
 Moindres carrés pour un modèle K( . ,θ )

1 – choix d’une famille H de distances identifiant θ
K( . ,θ )

2 – MCO sur une puissance c de K

Choix (Diggle) : c= 0.5 si X régulier, 0.25 si X avec agrégats

Atelier Spatial RASMA 215

Saint Louis du Sénégal
 Autres méthodes paramétriques

1 - Pour un PP, on sait définir une densité conditionnelle

comme pour un champ de Markov sur un réseau
(Jensen-Moller)
PVC pour un PP

2 - On maximise cette PVC : sous conditions d’ergodicité et

de faible dépendance, « bons résultats asymptotiques

3 - Estimation par MV : la difficulté est le calcul de la

constante de normalisation
l’approcher par Monte
Carlo

Atelier Spatial RASMA 216

Saint Louis du Sénégal
 Estimation d’un PP avec spatstat

Déclarer le modèle paramétrique : Poisson, Strauss,

Hardcore, StraussHard, Geyer, N-S, Thomas ….

• ppm : ajuste le modèle aux données via la pseudo

vraisemblance conditionnelle (si PPP, c’est la
vraisemblance); consulter les exemples, i.e :
> data(nztrees) puis > plot(nztrees)
> ppm(nztrees,~ x,Strauss(13),correction="periodic")

• rmh : simule le modèle estimé
procédure « visuelle » de

validation.

Atelier Spatial RASMA 217

Saint Louis du Sénégal
 Estimation d’un PP avec spatstat
(suite)

• fitin : donne l’interaction (0 pour un PPP) du PP estimé

• [Link] : ajuste un PP de Cox log-gaussien à
covariance exponentielle par MC sur la base de la
statistique K
• [Link] : idem pour un PP de Thomas
• [Link] : donne la log-vraisemblance de
l’ajustement pour un PPP; on en déduit le critère AIC.
• density : lissage de densité par convolution (noyau k
à choisir)

Atelier Spatial RASMA 218

Saint Louis du Sénégal
 Estimation d’un PP avec spatstat:
statistique et bande de confiance G, J et K
• Gest : estimation de la densité cumulée de la distance d’un point
de la configuration X à son ppv dans X
(F si distance d’un point courrant à son ppv dans X)
• Jest : J(r) = (1- G(r)) / (1 – F(r)).
J >1, =1 et <1 indique plus régulier, Poisson, moins régulier

• Kest : estime le moment d’ordre 2 K de Ripley.

• Kinhom : pour un PP inhomogène
• envelope : calcule les bandes de confiances par simulation des
statistiques de base (K, G, J …)

• [Link] : test du chi2 de l’hypothèse CSR

d’indépendance spatiale.

Atelier Spatial RASMA 219

Saint Louis du Sénégal
 Test de Monte Carlo de H(0) : « X est CSR »

Test basé sur un intervalle de confiance sur K

On observe x. Choisir L distances {h(1),h(2),h(L)}.

• Estimer (empiriquement) K par K*(h(l)) à ces L distances

• Sous H(0), estimation ρ* de l’intensité ρ d’un PPP
• Simulation de m réalisations x*(l), l=1,m, d’un PPP(ρ*) (i.e. m=20)
•
m estimations K*(i,h(l)) associées à chaque x*(i)
• Enveloppes inf et sup de ces m estimations
• Si K* se trouve entre ces 2 enveloppes, accepter H(0) (ici, niveau 10%)

On peut aussi comparer les K* à πh², valeur théorique de K pour un PPP.

Cette méthode décrit le principe du Bootstrap paramétrique.

Atelier Spatial RASMA 220

Saint Louis du Sénégal
Les 156 pins d’une forêt finlandaise
(finpines)

Atelier Spatial RASMA 221

Saint Louis du Sénégal
« Indépendance H(0) » contre « Neyman-Scott »
(données finpines)
• En continu, le graphe estimé de h
L(h) – h ( =0 si (H(0) ) pour x
• À gauche : bande de confiance pour m=40 simulations x* sous H(0)
• À droite : bande sous l’alternative N–S, avec en pointillé la courbe
h
K(h) théorique pour les paramètres de N–S estimés.

Atelier Spatial RASMA 222

Saint Louis du Sénégal
Atelier Spatial RASMA 223
Saint Louis du Sénégal