Il existe une grande variété de méthodes de simulation géostatistique. Pour les regrouper, on peut proposer la subdivision suivante :
Méthodes gaussiennes :
- Méthodes matricielles basées sur une décomposition de la matrice de covariance
- Méthode gaussienne séquentielle
- Moyennes mobiles
- Bandes tournantes
- Méthodes autorégressives
- Méthodes fréquentielles
- Etc.
Méthodes non-gaussiennes :
- Recuit simulé
- Simulation séquentielle d’indicatrices
- Méthodes utilisant des champs de probabilité (« p-fields »)
- Etc.
Dans la suite, nous examinerons quelques-unes de ces méthodes.
8.3 Méthode de simulation matricielle (LU, décomposition de Cholesky)¶
Cette méthode est simple à programmer et efficace pour simuler de petits champs. Elle permet de réaliser facilement des simulations conditionnelles comme non-conditionnelles.
8.3.1 Simulations non-conditionnelles¶
Soit ( n ) points à simuler avec une covariance ( C(h) ). On construit la matrice de covariance ( K ) de taille ( n \times n ) (identique à celle utilisée pour le krigeage simple).
On effectue la décomposition de Cholesky :
[
K = L U, \quad \text{où } L = U^{T}
]
On tire ensuite ( n ) valeurs indépendantes ( y_i \sim \mathcal{N}(0,1) ) pour ( i=1,\dots,n ), formant le vecteur ( \mathbf{y} ).
La réalisation simulée est alors :
[
\mathbf{z} = L \mathbf{y}
]
Cette réalisation ( \mathbf{z} ) a la covariance correcte, car :
[
\operatorname{Cov}[\mathbf{z}] = \operatorname{Cov}[L \mathbf{y}] = L \operatorname{Cov}[\mathbf{y}] L^T = L I L^T = L L^T = K
]
8.3.2 Simulations conditionnelles¶
Le principe est similaire, mais on distingue deux ensembles de points :
- ( N ) points conditionnants (observés),
- ( n ) points à simuler.
On construit la matrice de covariance ( K ) de taille ( (N + n) \times (N + n) ), et on la décompose en blocs via Cholesky :
[
K = L L^T =
]
- ( L_{11} ) correspond aux points conditionnants,
- ( L_{22} ) aux points à simuler,
- ( L_{21} ) relie les deux groupes.
On génère un vecteur aléatoire ( \mathbf{y} = \begin{bmatrix} \mathbf{y}_1 \ \mathbf{y}_2 \end{bmatrix} ), où ( \mathbf{y}_1 \in \mathbb{R}^N ) et ( \mathbf{y}_2 \in \mathbb{R}^n ) sont des vecteurs de variables indépendantes ( \mathcal{N}(0,1) ).
La simulation s’écrit :
[
= L \mathbf{y} =
]
Comme les valeurs ( \mathbf{z}1 = \mathbf{z}\text{obs} ) sont connues, on impose :
[
\mathbf{y}1 = L{11}^{-1} \mathbf{z}_1
]
On tire ensuite aléatoirement ( \mathbf{y}_2 \sim \mathcal{N}(0, I) ), et on calcule
[
\mathbf{z}2 = L{21} \mathbf{y}1 + L{22} \mathbf{y}_2
]
Notez que ( L_{21} \mathbf{y}_1 ) est constant pour toutes les réalisations conditionnelles et ne nécessite pas d’être recalculé.
Notes importantes¶
- Cette méthode est très efficace car une seule inversion de matrice (celle de ( L_{11} )) est nécessaire. Les simulations additionnelles se génèrent rapidement à faible coût.
- La principale limitation est la mémoire, car ( (N + n) ) ne peut excéder quelques milliers sans que la matrice ( K ) devienne trop volumineuse pour être stockée.
- La matrice ( K ) doit être non-singulière pour permettre la décomposition de Cholesky. Cela implique que les points à simuler ne doivent pas coïncider avec les points observés.
- D’autres décompositions sont possibles, par exemple via la décomposition en valeurs propres et vecteurs propres de ( K ).
8.4 Méthode de simulation séquentielle gaussienne (SGS)¶
La méthode SGS, basée sur les travaux de Paul Lévy (1937), est une méthode de simulation conditionnelle (le cas non-conditionnel est un cas particulier).
Principe général¶
- On dispose de ( N ) points conditionnants, avec covariance ( C(h) ).
- On souhaite simuler ( n ) points inconnus.
- On génère un ordre aléatoire (un chemin) pour visiter successivement les ( n ) points à simuler.
Workflow de la méthode SGS¶
Initialisation
- Ensemble des points connus : les ( N ) points conditionnants.
- Points à simuler : ( {x_1, x_2, \dots, x_n} ) dans un ordre aléatoire.
Pour chaque point ( x_i ) dans l’ordre aléatoire :
a. Calcul du krigeage simple sur ( x_i ) en utilisant :
- Tous les points conditionnants ( N ),
- Tous les points déjà simulés ( {x_1, \dots, x_{i-1}} ).
Le krigeage simple fournit :
- La valeur krigée (espérance conditionnelle) :
[ Z_i^* = E[Z(x_i) \mid \text{données connues}] ] - La variance de krigeage (variance conditionnelle) :
[ \sigma_k^2 = \operatorname{Var}[Z(x_i) \mid \text{données connues}] ]
b. Simulation de la valeur au point ( x_i ) :
On tire aléatoirement une valeur suivant la loi normale conditionnelle :
[ Z(x_i) \sim \mathcal{N}(Z_i^*, \sigma_k^2) ]c. Ajout du point simulé ( x_i ) avec sa valeur ( Z(x_i) ) à l’ensemble des points connus.
Répéter jusqu’à ce que tous les ( n ) points soient simulés.
Remarques importantes¶
- La taille des systèmes de krigeage augmente avec le nombre total de points connus (points conditionnants + points simulés).
- En pratique, pour limiter la complexité, on utilise une fenêtre locale pour ne prendre en compte que les voisins proches (effet d’écran).
- Pour un effet de pépite dans le modèle :
- En non-conditionnel, il est préférable d’ajouter cet effet après simulation sous forme d’une erreur indépendante de variance ( C_0 ).
- En conditionnel, ce n’est pas possible car les observations incluent déjà cet effet.
- Certaines covariances (ex. sphérique, gaussien) peuvent être plus difficiles à reproduire fidèlement.
- L’ordre de visite des points peut influencer les résultats :
- Éviter un parcours systématique selon une direction.
- Recommandation : une visite en deux étapes,
- Première visite aléatoire sur une maille large pour capturer la grande échelle,
- Deuxième visite aléatoire sur une maille plus fine pour affiner la grille.
- Cette stratégie peut être étendue à plusieurs passes successives.