Krigagem

Krigagem ou krigeagem (do inglês Kriging) é um método de regressão usado em geoestatística para aproximar ou interpolar dados. A teoria de Kriging foi desenvolvida a partir dos trabalhos do seu inventor, Daniel G. Krige, pelo matemático francês Georges Matheron, no começo dos anos sessenta. Na comunidade estatística, também é conhecido como “Processo Gaussiano de Regressão”. A estimação com base em apenas um atributo insere-se no âmbito da Krigagem; a estimação de um atributo à custa de outros atributos insere-se no âmbito da Co-krigagem.

Introdução Krigagem poderá ser entendido como uma predição linear ou uma forma da Inferência bayesiana. Parte do princípio que pontos próximos no espaço tendem a ter valores mais parecidos do que pontos mais afastados. A técnica de Krigagem assume que os dados recolhidos de uma determinada população se encontram correlacionados no espaço. Isto é, se num aterro de resíduos tóxicos e perigosos a concentração de Zinco num ponto p é x, é muito provável que se encontrem resultados muito próximos de x quanto mais próximos se estiver do ponto p (princípio da geoestatística). Porém, a partir de determinada distância de p, certamente não se encontrarão valores aproximados de x porque a correlação espacial pode deixar de existir.

Considera-se o método de Krigagem do tipo BLUE (Best Linear Unbiased Estimator - Melhor Estimador Linear não-Viciado): é linear porque as suas estimativas são combinações lineares ponderadas dos dados existentes; é não enviesada pois procura que a média dos erros (desvios entre o valor real e o valor estimado) seja nula; é a melhor porque os erros de estimação apresentam uma variância (variância de estimação) mínima. O termo Krigagem abrange um conjunto de métodos, sendo os mais usuais os seguintes:

Tipos de Krigagem

Krigagem Simples Assume que as médias locais são relativamente constantes e de valor muito semelhante à média da população que é conhecida. A média da população é utilizada para cada estimação local, em conjunto com os pontos vizinhos estabelecidos como necessários para a estimação.

Krigagem Normal As médias locais não são necessariamente próximas da média da população usando-se apenas os pontos vizinhos para a estimação. É o método mais usado em problemas ambientais.

Co-krigagem É uma extensão do anterior a situações em que duas ou mais variáveis são espacialmente dependentes e a variável que se quer estimar não está amostrada com a intensidade com que estão as outras variáveis dependentes, utilizando-se os valores destas e as suas dependências para estimar a variável requerida.

Conceitos matemáticos O Método de Krigagem serve-se de diversas teorias explanadas na estatística. No entanto, para deixarmos mais claras as teorias de estatística usadas e mais direcionadas ao escopo deste texto, explicaremos alguns conceitos. Semi-variância e semi-variograma thumb|300px|Variograma A semi-variância é a medida do grau de dependência espacial entre duas amostras. A magnitude da semi-variância entre dois pontos depende da distância entre eles, implicando em semi-variâncias menores para distâncias menores e semi-variâncias maiores para distâncias maiores. O gráfico das semi-variâncias em função da distância a um ponto é chamado de Semi-variograma. A partir de uma certa distância a semi-variância não mais aumentará com a distância e estabilizar-se-á num valor igual à variância média, dando a esta região o nome de silo ou patamar (sill). A distância entre o início do semi-variograma e o começo do silo recebe o nome de range ou amplitude ou alcance. Ao extrapolarmos a curva do semi-variograma para a distância zero, podemos chegar a um valor não-nulo de semi-variância. Este valor recebe o nome de Efeito Pepita (Nugget Effect).

Modelos de Variograma No Método de Krigagem normalmente são usados quatro tipos de variogramas. Neles, são usadas as seguintes variáveis: : variância : nugget : efeito pepita : variância assintótica : distância de separação Linear Este modelo não apresenta silo e é muito simples. Sua curva pode ser representada por:

Esférico A forma esférica é a mais utilizada e possui silo. Sua forma é definida por:

Exponencial A curva do variograma exponencial respeita a seguinte equação:

Gaussiano A forma gaussiana é dada por:

O Método de Krigagem Determinação do Semivariograma Toma-se como base a simulação de um sistema de duas dimensões (2D) que contém um número finito de pontos onde é possível a medição de uma grandeza qualquer. Após a aquisição destes dados, iniciar-se-á a interpolação por Kriging buscando alcançar uma maior resolução. O primeiro passo é construir um semivariograma experimental. Para tal, calcula-se a semivariância de cada ponto em relação aos demais e insere-se no gráfico da semivariância pela distância.

A partir deste gráfico estima-se o modelo de variograma que melhor se aproxima da curva obtida. O efeito pepita pode estar presente no semivariograma experimental e deve ser considerado. Determinado o modelo do semivariograma a ser usado, inicia-se a fase de cálculos. Sendo o semivariograma uma função que depende da direção, é natural que apresente valores diferentes conforme a direção, recebendo este fenómeno o nome de Anisotropia. Caso o semivariograma apresente uma forma semelhante em todas as direções do espaço, só dependendo de h, diz-se que a estrutura é Isotrópica, i. e., sem direções privilegiadas de variabilidade.

Cálculo dos Pesos Considere, para o cálculo da Krigagem, a seguinte fórmula:

onde é o número de amostras obtidas, é o valor obtido no ponto e é o peso designado ao ponto . A fim de obter os pesos de cada um dos pontos, para cada um deles é realizado um cálculo de . Tal procedimento depende do tipo de Krigagem que está a ser utilizado. Salienta-se a seguinte notação: : peso do j-ésimo ponto : valor da semi-variância de : variável temporária

Krigagem Normal Neste caso é utilizada a média local dos pontos amostrados. Por conseguinte, deve-se normalizar a média dos pesos. Consequentemente, tem-se um resultado mais preciso do que a Krigagem Simples. Utilizar-se-ão as seguintes equações para a determinação dos valores dos pesos no p-ésimo ponto:

Krigagem Simples Para este caso, utiliza-se a média de todos os dados. Implica-se, portanto, em não se normalizar a média local dos pesos, como no caso anterior. Assim, teremos quase que a mesma equação, exceto pela exclusão de e pela última equação. A característica principal deste método é a geração de gráficos mais lisos e mais esteticamente suaves. Deve-se salientar que este caso é menos preciso que o caso anterior. Os valores dos pesos para o p-ésimo ponto serão dados por:

Obtendo o Ponto Interpolado Ao obtermos os valores de , calcula-se o valor de :

Desta maneira, calcula-se o valor interpolado para todos os pontos desejados. Ressalta-se que somente devem ser utilizados os valores adquiridos acima. Interpolando Outros Pontos A obtenção do valor interpolado em um outro ponto requer a repetição de todos os cálculos realizados a partir da obtenção do modelo de variograma. Desta forma, para aumentarmos a resolução que é pretendida, deve-se recorrer à métodos matemáticos para a resolução computacional. Diversos códigos foram desenvolvidos para esta resolução, mas um dos melhores algoritmos pode ser obtido no link abaixo. Ele fora desenvolvido inicialmente para a linguagem Fortran, porém foi recodificado para C com a ajuda da biblioteca fortran2c e apresenta-se totalmente em C: Kriging Interpolation Algorithm in C

Algoritmo da krigagem (Python)

A implementação bidimensional que se segue é feita na linguagem Python (versão 2.7.2) com recurso à biblioteca NumPy tendo por esse motivo o seguinte cabeçalho de importações (está considerado também a biblioteca matplotlib usada para visualização):

from __future__ import division import numpy as np import matplotlib.pyplot as plt

Repare-se que a primeira linha de importação refere-se à possibilidade de uma divisão entre inteiros poder ser lida como podendo ter um resultado decimal não sendo de maior consequência na transcrição do seguinte procedimento para qualquer outra linguagem (inclusive o próprio Python em versões mais recentes). A função que se segue faz, com efeito, a interpolação por krigagem simples com recurso a uma matriz de dados na qual a primeira coluna é a coluna da posição X, a segunda posição Y, e terceira valor da amostra. O último argumento trata-se de uma tupla ou lista com o valor da amplitude na direção X (posição 0) e direção Y (posição 1).

def KS(dados,blocos,variograma): resultado = np.zeros(blocos)

cov_angulos = np.zeros((dados.shape,dados.shape)) cov_distancias = np.zeros((dados.shape,dados.shape)) K = np.zeros((dados.shape,dados.shape)) for i in xrange(dados.shape-1): cov_angulos[i,i:]=np.arctan2((dados[i:,1]-dados[i,1]),(dados[i:,0]-dados[i,0])) cov_distancias[i,i:]=np.sqrt((dados[i:,0]-dados[i,0])**2+(dados[i:,1]-dados[i,1])**2) for i in xrange(dados.shape): for j in xrange(dados.shape): if cov_distancias[i,j]!=0: amp=np.sqrt((variograma*np.cos(cov_angulos[i,j]))**2+(variograma*np.sin(cov_angulos[i,j]))**2) K[i,j]=dados[:,2].var()*(1-np.e**(-3*cov_distancias[i,j]/amp)) K = K + K.T

for x in xrange(resultado.shape): for y in xrange(resultado.shape): distancias = np.sqrt((x-dados[:,0])**2+(y-dados[:,1])**2) angulos = np.arctan2(y-dados[:,1],x-dados[:,0]) amplitudes = np.sqrt((variograma*np.cos(angulos[:]))**2+(variograma*np.sin(angulos[:]))**2) M = dados[:,2].var()*(1-np.e**(-3*distancias[:]/amplitudes[:])) W = LA.solve(K,M) resultado[x,y] = np.sum(W*(dados[:,2]-dados[:,2].mean()))+dados[:,2].mean() return resultado

O algoritmo acima feito não usa qualquer motor de procura ou limite de amostras. Desta maneira todas as amostras são usadas no sistema de krigagem pelo que apenas é necessário calcular a matriz M, com o valor de variogram entre todas as amostras, uma vez. Isto é feito na seguinte parcela de código:

De resto apenas é necessário calcular a matriz M por cada nó da malha de resultados de maneira a resolver o sistema K*W=M (o objectivo é obter W, a matriz de pesos). Após o cálculo dos pesos basta multiplicar os mesmos aos dados e proceder à soma final (com atenção à forma de cálculo em krigagem simples onde: resultado = Soma(W*dados - média) + média). Repare-se que nesta krigagem simples poderia haver um outro argumento que é a média imposta pelo utilizador que posteriormente é utilizado no cálculo do ponto a estimar (nesta implementação é automáticamente considerado a média dos dados de entrada).

Para fazer krigagem normal as modificações são poucas. A matriz K tem mais uma coluna e mais uma linha e a matriz M tem mais um elemento também. A multiplicação dos pesos (sem considerar o último valor) pelos dados dá-nos o valor estimado em cada bloco.

def OK(dados,blocos,variograma): resultado = np.zeros(blocos)

Krigagem

Źródło i licencja