Comment fonctionne l'interpolation Kriging?

10

Je travaille sur un problème dans lequel je dois utiliser Kriging pour prédire la valeur de certaines variables en fonction de certaines variables environnantes. Je veux implémenter son code par moi-même. Donc, j'ai parcouru trop de documents pour comprendre comment cela fonctionne, mais j'étais tellement confus. En général, je comprends que c'est une moyenne pondérée, mais je ne pouvais pas comprendre complètement le processus de calcul du poids, puis prédire la valeur d'une variable.

Quelqu'un peut-il m'expliquer en termes simples les aspects mathématiques de ces méthodes d'interpolation et comment cela fonctionne?

Dania
la source
3
L'implémentation de code est un excellent outil d'apprentissage mais ne peut pas être recommandée pour travailler sur des problèmes réels. Au moment où vous obtenez le code écrit, débogué et testé, vous découvrirez qu'il a besoin d'un ordre de grandeur plus d'efforts pour fournir des outils supplémentaires pour l'analyse des données d'exploration spatiale, la variographie, la validation croisée du variogramme, la recherche de quartier et la post- traitement des résultats krigés. Un compromis raisonnable et efficace serait de commencer avec du code de travail, tel que GSLib ou GeoRGLM , et de le modifier.
whuber
Merci beaucoup, c'est une excellente idée, mais je veux aussi comprendre l'aspect mathématique du Kriging, avez-vous une ressource qui l'explique clairement en termes simples? Je vous remercie.
Dania

Réponses:

15

Cette réponse consiste en une section d'introduction que j'ai écrite récemment pour un article décrivant une extension spatio-temporelle (modeste) de "Universal Kriging" (UK), qui est elle-même une modeste généralisation de "Ordinary Kriging". Il comprend trois sous-sections: la théorie donne un modèle statistique et des hypothèses; L'estimation passe brièvement en revue l'estimation des paramètres des moindres carrés; et Prediction montre comment le krigeage s'intègre dans le cadre GLS (Generalized Least Squares). J'ai fait un effort pour adopter une notation familière aux statisticiens, en particulier les visiteurs de ce site, et pour utiliser des concepts qui sont bien expliqués ici.

Pour résumer, le krigeage est la meilleure prédiction linéaire sans biais (BLUP) d'un champ aléatoire. Cela signifie que la valeur prédite à tout emplacement non échantillonné est obtenue sous la forme d'une combinaison linéaire des valeurs et covariables observées aux emplacements échantillonnés. La valeur (inconnue, aléatoire) présente une corrélation supposée avec les valeurs d'échantillon (et les valeurs d'échantillon sont corrélées entre elles). Cette information de corrélation est facilement traduite en variance de la prédiction. On choisit des coefficients dans la combinaison linéaire (les "poids de krigeage") qui rendent cette variance aussi petite que possible, sous réserve d'une condition de biais nul dans la prédiction. Les détails suivent.


Théorie

Le Royaume-Uni comprend deux procédures - l'une d'estimation et l'autre de prédiction - réalisées dans le contexte d'un modèle GLS pour une zone d'étude. Le modèle GLS suppose que les données d'échantillon sont le résultat d'écarts aléatoires autour d'une tendance et que ces écarts sont corrélés. Une tendance est conçue dans le sens général d'une valeur qui peut être déterminée par une combinaison linéaire de coefficients inconnus (paramètres) . (Tout au long de ce post, le prime désigne la transposition matricielle et tous les vecteurs sont considérés comme des vecteurs de colonne.)p β = ( β 1 , β 2 , ... , β p ) ' 'zi, (i=1,2,...,n)pβ=(β1,β2,,βp)

À n'importe quel endroit dans une zone d'étude, il existe un tuple d'attributs numériques appelés «variables indépendantes» ou «covariables». (En générale, est un «terme constant», et peuvent être des coordonnées spatiales, et les supplémentaires peuvent représenter des informations spatiales ainsi que d'autres informations auxiliaires disponibles à tous les emplacements de la zone d'étude, telles que la porosité d'un aquifère ou distance d'un puits de pompage.) À chaque emplacement de données , en plus de ses covariables , l'observation associéey 1 = 1 y 2 y 3 y i i y i = ( y i 1 , yy=(y1,y2,,yp)y1=1y2y3yiiz i Z i y jeyi=(yi1,yi2,,yip)ziest considérée comme une réalisation d'une variable aléatoire . En revanche, les sont considérés comme des valeurs déterminées par ou caractérisant les points ou petites régions représentées par les observations (les données «supportent»). Les ne sont pas considérés comme des réalisations de variables aléatoires et ne doivent être liés aux propriétés d'aucun des .ZiyiZ iyiZi

La combinaison linéaire exprime la valeur attendue de en termes de paramètres , qui est la valeur de la tendance à l'emplacement . Le processus d'estimation utilise les données pour trouver des valeurs qui représentent les paramètres inconnus , tandis que le processus de prédiction utilise les données aux emplacements pour calculer une valeur à un emplacement non échantillonné , qui est ici indexé comme . Les cibles d'estimation sont fixes ( c.-à-d.Z i βi β i β i i=1,2,...,ni=0 z 0 y

E[Zi]=yiβ=yi1β1+yi2β2++yipβp
Ziβiβ^iβii=1,2,,ni=0, non aléatoire) alors que la cible de prédiction est aléatoire, car la valeur inclut une fluctuation aléatoire autour de sa tendance . En règle générale, les prévisions sont faites pour plusieurs emplacements en utilisant les mêmes données en faisant varier l'emplacement . Par exemple, des prédictions sont souvent faites pour cartographier une surface le long d'une grille régulière de points adaptés au contour. z00y0β0

Estimation

Le krigeage classique suppose que les fluctuations aléatoires ont des valeurs attendues de zéro et leurs covariances sont connues. Ecrire la covariance entre et comme . En utilisant cette covariance, l'estimation est effectuée à l'aide de GLS. Sa solution est la suivante: où est le vecteur des observations, (la «matrice de conception») est la matrice par dont les lignes sont les vecteursZ i Z j c i j β = H z , H = ( Y ' C -ZiZiZjcijz=(z1,z2(cij)nnpn

β^=Hz, H=(YC1Y)1YC1
n Y = ( y i j ) n p y i , 1z=(z1,z2,,zn)nY=(yij)npC =yi,1in et est la matrice de covariance par- qui est supposée être inversible (Draper et Smith (1981), section 2.11) . La matrice par , qui projette les données sur les estimations de paramètres , est appelée la «matrice chapeau». La formulation de comme application de la matrice de chapeau aux données montre explicitement comment les estimations des paramètres dépendent linéairement des données. Les covariancesC=(cij)nnpnz la ß la ß C = ( c i j )Hzβ^β^C=(cij) sont calculés de façon classique à l'aide d'un variogramme qui donne la covariance en termes de localisation des données, bien qu'il soit sans importance de savoir comment la covariance est réellement calculée.

Prédiction

Le Royaume-Uni prédit également au moyen d'une combinaison linéaire des données Les sont appelés les «poids de krigeage» pour la prédiction de . UK réalise cette prédiction de en répondant à deux critères. Tout d' abord, la prédiction doit être impartiale, qui est exprimé en exigeant que la combinaison linéaire des variables aléatoires est égal en moyenne: Cette attente est reprise sur le jointz 0 = λ 1 z 1 + λ 2 z 2 + + λ n z n = λ ' z . λ i z 0 z 0 Z i Z 0 0 = E [ Z 0 - Z 0 ] = E [ λ 1 Z 0 Z = ( Z 1 ,z0

z^0=λ1z1+λ2z2++λnzn=λz.
λiz0z0ZiZ0n+
0=E[Z^0Z0]=E[λZZ0].
n+1-distribution variable de et . La linéarité de l'attente associée à l'hypothèse de tendance (1) implique: Z00 = E [ λ Z - Z 0 ] = λ E [ Z ] - E [ Z 0 ] = λ ( Y β ) - y 0 βZ=(Z1,Z2,,Zn)
0=E[λZZ0]=λE[Z]E[Z0]=λ(Yβ)y0β=(λYy0)β=β(Yλy0)

peu importe ce que peut être. Ce sera le cas à condition queβ

Y^λ=y0.

Parmi toutes les solutions possibles de ce système d'équations sous-déterminé, UK choisit pour minimiser la variance de l'erreur de prédiction . En ce sens, le Royaume-Uni est «le meilleur» parmi tous les prédicteurs linéaires non biaisés. Parce que cette dernière relation implique que l'erreur de prédiction est nulle en moyenne, la variance est simplement l'attente de l'erreur de prédiction au carré: où est le vecteur de covariances entreZ 0 - Z 0 V a r ( Z 0 - Z 0 ) = E [ (λZ^0Z0...,c

Var(Z^0Z0)=E[(Z^0Z0)2]=E[(λZZ0)2]=c002λc0+λCλ
c0=(c01,c02,,c0n)Z0et le et est la variance de . Zi, i1c00Z0

Pour minimiser la variance, différencier par rapport à et introduire un vecteur de multiplicateurs de Lagrange à incorporer dans la contrainte . Cela donne un système d' équations linéaires , écrites sous forme de matrice de blocs comme où représente un parp μ Y ' λ = y 0 n + p ( C Y Y ' 0 ) ( λ μ ) = ( c 0 y 0 ) 0 p p 1 n n λ λ =λpμY^λ=y0n+p

(CYY0)(λμ)=(c0y0)
0ppmatrice de zéros. En écrivant pour la matrice d'identité par , la solution unique pour est donnée par 1nnλ
λ=Hy0+C1(1YH)c0.

(Les lecteurs familiers avec la régression multiple peuvent trouver instructif de comparer cette solution à la solution basée sur la covariance des équations normales des moindres carrés ordinaires , qui a presque exactement la même apparence, mais sans termes de multiplicateur de Lagrange.)

Cette relation présente les poids de krigeage, , comme la somme d'un terme dépendant uniquement de la matrice du chapeau et des covariables à l'emplacement de prédiction , plus un terme dépendant des covariances parmi les données et le prédictant, . La substitution dans le côté droit de l'équation de variance donne la variance de prédiction de krigeage, qui peut être utilisée pour construire des limites de prédiction autour de .λ Z 0 z 0[Hy0]Z0z^0

whuber
la source
1
Merci beaucoup whuber, c'est exactement ce que je recherche. Vous avez résolu ce problème pour moi, maintenant je comprends Kriging. J'apprécie vraiment votre aide, merci beaucoup.
Dania
Explication fantastique. Une question: que représente ? Comment est-il défini? Cela fait-il partie des données? Que signifie le premier? Cette variable est introduite sans être définie, donc je suis un peu confus quant à sa définition. Y^
DW
@DW Le premier indique la transposition tout au long de ce post. Ainsi, en prenant la transposition de la définition dans la réponse, nous pouvons décrire cette matrice comme " est la matrice par dont les colonnes sont les vecteurs . " Il encapsule ainsi l'ensemble de données des covariables. pn y i ,1inY=(yji)pnyi,1in
whuber