Logistica nivelamento 02

Prévia do material em texto

Geostatística
Prof. Pedro Campos
Instituto Ciberespacial
Agredecimentos ao Prof. Joaquim Barbosa que
gentilmente cedeu parte de seu curso de
Geoestatística para uso nessa disciplina.
Prof. Pedro Campos
Estimação Local
• Semivariograma permite verificar e modelar a dependência
espacial de uma variável.
• Modelo de dependência espacial estimação de um atributo
em localizações não amostradas.
• Krigagem: interpolador que utiliza o semivariograma em sua
modelagem. Homenagem ao matemático sul-africano Daniel Krige
(mineração).
• Algumas Áreas de Aplicações:
• mapeamento geológico (Verly et al., 1984)
• mapeamento solo (Burgess e Webster, 1980)
• mapeamento hidrológico (Kitanidis et. al., 1983)
• mapeamento atmosférico (Lajaunie, 1984)
1. Métodos de estimação:
Tipos de interpoladores
• Envolvem combinação linear dos n valores em pontos vizinhos. 
x1 x2
x3
x4
x0
Inverso do Quadrado
da Distância
KrigagemMédia Local
n
1=
∑=
=
i
ii
*
λ
ZλZ
n
1i0
xx
2. Estimador da KRIGAGEM
x1 x2
x3
x4
x0
• Os pesos são extraídos de uma análise de correlação espacial
baseada no Semivariograma
Análise de Correlação Espacial
baseada no Semivariograma
1
Ajuste do Semivariograma
Experimental
2
4
Estimador de
KRIGAGEM
Modelo de ajuste do semivariograma
3
( ) ( )[ ]hhhh SphCCa2
1
a2
3CC 1010
3
+=-+=
ú
ú
ú
ú
û
ù
ê
ê
ê
ê
ë
é
÷
÷
÷
÷
ø
ö
ç
ç
ç
ç
è
æ
÷
÷
÷
ø
ö
ç
ç
ç
è
æ
γ
Estimador da krigagem
Aplicação da krigagem: 
1. Assume-se conhecidas as realizações
2. Semivariograma da variável já foi calculado 
3. Interesse em estimar um valor z* na localização u.
{ }nαz α ,...,1=),(u
.
• Krigagem : nome genérico para uma família de algoritmos de
regressão de mínimos quadrados generalizados.
[ ]∑
)(
1=
* )(-)()(=)(-)(
u
uuuuu
n
α
ααα mZλmZ
.
)(uαλ )( αz u
)( αZ u
: peso atribuído ao dado interpretado como uma
realização da VA
)(um )( αm ue )(uZ )( auZ: médias das VAs e
(3.0)
• O número de dados envolvidos na estimação, assim como
seus pesos, pode variar de uma localização a outra. Na
prática, somente )(un dados mais próximos da localização u sendo
estimada são utilizados.
•Três modelos de krigagem podem ser distinguidos de acordo
com o modelo considerado para a tendência )(um
1. Krigagem Simples: )(um = m ,conhecida e constante ˅ u € A
)(u'm
2.Krigagem Ordinária: considera flutuações locais da média pela
limitação do domínio de estacionaridade da média à vizinhança local
W(u): = constante, mas desconhecida ˅ u’ € W(u) .
3.Krigagem com modelo de Tendência ou Krigagem Universal: a
média local desconhecida )(u'm
vizinhança local W(u). A componente de tendência é modelada como
uma combinação linear de funções fk(u) de coordenadas:
varia suavemente dentro de cada
∑
0=
)()(=)(
K
k
kk fam u'u'u' kk aa ≈)(u' constante mas desconhecido ˅ u’ €W(u) 
PROPRIEDADES DA KRIGAGEM COMO INTERPOLADOR 
1. É um interpolador exato. No ponto ui 
e a variância é zero.
2. O pesos da krigagem são calculados com auxílio do
variograma.
3. O valor máximo dos pesos é 1. Entretanto, eles podem ser
negativos. Portanto, não é verdadeira a hipótese
)(=)(* ii ZZ uu
[ ] [ ])(min≤)(≤)(max * iii ZZZ uuu
4. Os pesos da krigagem não influenciam nos valores medidos. Se a
mesma configuração aparece em duas localizações diferentes os
pesos da krigagem serão os mesmos, independentemente dos valores
medidos. Os valores medidos influenciam o variograma que á a base
para os pesos da krigagem.
5. Os pesos da krigagem mostram o efeito screening. Pontos
distantes têm peso menor do que pontos mais próximos.
2.1.1. Krigagem Ordinária (KO)
- Efeito proporcional: a média local pode variar significativamente na
área de estudo.
2.1. Principais tipos de Krigagem
- KO: considera tal variação local da média limitando o domínio de
estacionaridade da média à vizinhança local W(u) centrado na
localização u sendo estimada.
-O estimador linear: combinação linear das n(u) VAs Z(uα) mais a
média constante local m(u):
[ ]∑
)(
1=
* )(-)()(=)(-)(
u
uuuuu
n
α
ααα mZλmZ
[ ] )(+)(-)()(=)( ∑
)(
1=
* uuuuu
u
mmZλZ
n
α
ααα
)()(-1)()(=)( _
)(
1=
)(
1=
* ∑ ∑ uuuuu
u u
mλZλZ
n
α
n
α
ααα
)(+)()(-)()(=)( ∑ ∑
)(
1=
)(
1=
* uuuuuu
u u
mmλZλZ
n
α
n
α
αααα
∑
)(
1=
* )()(=)(
u
uuu
n
α
α
OK
αOK ZλZ
A melhor estimativa de é obtida quando:)(* uOKZ
a) O estimador é não tendencioso
{ } 0=)(-)(* uu ZZE OK
para ∑
)(
1=
1=)(
u
u
n
α
OK
αλ
b) A variância da estimativa é mínima
{ } [ ]{ } mínimaZZEZZVar OKOK =)(-)(=)(-)( 2** uuuu
SISTEMA DE EQUAÇÕES DA KRIGAGEM. 
• Hipótese de estacionaridade de segunda ordem (ou hipótese
intrínseca):
mZE =)]([ u para todo u ∊A
• Estimador não tendencioso: erro médio é igual a zero.
{ }
0=)(-)(=
)(-)()(=)(-)( ∑
)(
1=
*
uu
uuuuu
u
mm
mmλZZE
n
α
OK
αOK
1=)(∑
)(
1=
u
u
n
α
OK
αλ
{ } =)(-)(=)(-)(=)(
2(
1=
**2 ∑ uuuuu
u
ZZλEZZVarσ
n
α
OKαOK
)(+)()(2-)()(= ∑ ∑∑
)(
1=
2
)(
1=
)(
1=
u uu
uuuuu
n
α
n
β
αα
n
β
βαβα ZZZλZZλλE
Sabe-se que:
Var{Z(u)} = E[Z2(u) – m2 ]=s2 = C(0) para todo u Є A
Cov{Z(u),Z(u+h)} = E[Z(u)Z(u + h) – m2]= C(h) para todo u Є A
+)0(++),(2-+),(=∑ ∑∑
)(
1=
22
)(
1=
)(
1=
2
u uu
uuuu
n
α
n
β
αα
n
β
βαβα mCmCλmCλλ
)0(+),(2-),(=)( ∑ ∑∑
)(
1=
)(
1=
)(
1=
2
u uu
uuuuu
n
α
n
β
αα
n
β
βαβα CCλCλλσ
A equação acima deve ser minimizada sob a restrição de que 1=)(∑
)(
1=
u
u
n
α
αλ
Este processo de minimização é feito através de técnicas de Lagrange 
Para satisfazer a condição de variância minima, é preciso que as N
derivadas parciais
)(∂
1-)()(2-)(∂ ∑
)(
1=
2
u
uuu
u
α
n
α
α
λ
λμσ
sejam igualadas a 0 (zero); é um multiplicador de Lagrange. μ
Fazendo a derivada parcial em relação a :
0=)(2-)(2-)()(2∑
)(
1=
uu,uu,uu
u
μCCλ α
n
β
βαβ a = 1 ,..., n(u) 
)(ual
Fazendo a derivada parcial em relação a :µ(u)
0=1-)(∑
)(
1=
u
u
n
α
αλ
Cancelando o fator 2, rearranjando tem-se o sistema da krigagem 
expresso em termos de função covariância: 
1=
)(=)(-)()(
∑
∑
)(
1=
)(
1=
u
u
u,uuu,uu
n
β
β
α
n
β
βαβ
λ
CμCλ
α= 1 ,..., n(u) 
Rearrajando as primeiras N equações no sistema de krigagem:
)(+)(=)()(∑
)(
1=
u,uuu,uu
u
α
n
β
βαβ CμCλ
)0(+),(2-),(=)( ∑ ∑∑
)(
1=
)(
1=
)(
1=
2
u uu
uuuuu
n
α
n
β
αα
n
β
βαβα CCλCλλσ
Substituindo a equação acima na equação da variância 
)0(+),(2-)),(+((=)( ∑ ∑∑
)(
1=
)(
1=
)(
1=
2
u uu
uuuuu)u
n
α
n
β
αα
n
β
αα CCλCμλσ
)0(+),(2-)),(+((=)( ∑ ∑
)(
1=
)(
1=
2
u u
uuuuu)u
n
α
n
β
αααα CCλCμλσ
∑∑
)(
1=
)(
1=
2 ),(2-),(+(+)0(=)(
uu
uuuuu)u
n
β
αα
n
β
αα CλCλμCσ
∑
)(
1=
2 ),(-(+)0(=)(
u
uuu)u
n
β
ααCλμCσ
Expressão da variância mínima da estimativa 
Considerando a relação o sistema de 
krigagem e da variância da estimativa, podem ser expressos 
em termos de semivariograma como
γ(h) = C(0) – C(h)
[ ]
1=
)(-)0(=)(-)(-)0()(
∑
∑
)(
1=
)(
1=
u
u
u,uuu,uu
n
β
β
α
n
β
βαβ
λ
γCμγCλ
Graças à condição de não-tendenciosidade, ,
o termo da variância C(0) é cancelado,
1=)(∑
)(
1=
u
u
n
α
αλ
1=
)(=)(+)()(
∑
∑
)(
1=
)(
1=
u
u
u,uuu,uu
n
β
β
α
n
β
βαβ
λ
γμγλ
A variância da estimativa, fica
[ ]∑
)(
1=
2 ),(-)0(-(+)0(=)(
u
uuu)u
n
β
αα γCλμCσ
∑
)(
1=
2 ),(+(=)(
u
uuu)u
n
β
ααγλμσ
Sistema matricial 
O sistema de krigagem pode ser escrito em notação matricial
[ ][ ] [ ] b C =λ
Função de Covariância
[ ][ ] [ ]b=λγ
Função de Semivariograma
Soluções [ ][ ] [ ]bC =1λ [ ] [ ] [ ]b1= γλ
[ ]C
[ ]γ
[ ]λ
[ ]b
Matriz de Covariância
Matriz de Semivariância
Matriz dos pesos procurados
Lado direito das equações 
do sistema de krigagem
A variância da estimativa, podem ser expressas em notação 
matricial, como 
[ ] [ ] λC= σ T b0u -)()(2 [ ] [ ]bu Tλσ =)(2
Função de Covariância Função de SemivariogramaSuponha que N=2. Então a matriz (ou [C]) é uma matriz 3x3
e pode ser, explicitamente, escrita como
[ ]γ
[ ]
011
1
1
222
2
)u,u(γ)u,u(γ
)u,u(γ)u,u(γ
= 1
111
γ [ ]
μ
λ
λ
= 2
1
λ )u,(uγ
)u,u(γ
= 
1
1
02
0
[b]
Exemplo 3.1. Considere 3 pontos em uma linha reta. Sabe-se 
que os valores de Z(u1)=2 e Z(u2)=4 e deseja-se estimar o valor 
do terceiro ponto. A Figura abaixo mostra a configuração, onde 
u=0, u1=1 e u2= – 2. Supor que o variograma seja linear 
u2 u1u
hγ =)(h
Equações da krigagem
1=+3+0 21 μλλ
2=+0+3 21 μλλ
1=+ 21 λλ
Sistema matricial 
[ ]
011
103
130
=γ [ ]
μ
λ
λ
= 2
1
λ = 
1
2
1
[b]
[ ] [ ] [ ]b1= γλSolução [ ]
0
333,0
667,0
== 2
1
μ
λ
λ
λ
667,2=4*333,0+2*667,0=
)()(=)( ∑
)(
1=
*
u
uuu
n
α
α
OK
αOK ZλZ
[ ] [ ] 333,1==)(2 bu Tλσ
Exemplo 3.2. Supondo que a configuração seja alterada e u2
seja movido para o outro lado da origem, i.e., u2=2
u2u1u
1=+1+0 21 μλλ
2=+0+1 21 μλλ
1=+ 21 λλ
Equações da krigagem
Sistema matricial 
[ ]
011
101
110
=γ [ ]
μ
λ
λ
= 2
1
λ = 
1
2
1
[b]
[ ] [ ] [ ]b1= γλSolução [ ]
0,1
0
0,1
== 2
1
μ
λ
λ
λ
0,2=4*0+2*0,1=
)()(=)( ∑
)(
1=
*
u
uuu
n
α
α
OK
αOK ZλZ
[ ] [ ] 0,2==)(2 bu Tλσ
OBS: A configuração dos 
dados tem um papel 
importante na krigagem
Exemplo 3.3. Estimar o valor da variável regionalizada Z4
indicada para a configuração abaixo, onde Z1=600, Z2=400 e 
Z3=500
SOLUÇÃO:
Amostra Coord. X Coord. Y Var Z
1 10 10 600
2 40 10 400
3 10 50 500
4 30 40 ?
( ) 2222211 )-(+...+)-(+)-(= pp yxyxyxd yx,
Variograma: Modelo esferico -3
2
=)(
3
a
h
a
hC
hγ
Supor: C = 100
a = 40
[ ]
0111
1094,35100
194,35041,91
110041,910
=γ [ ]
0168,533285,02207,04509,0
3285,00148.00125,00022,0
2207,00125,00161,00036,0
4509,00022,00036,00058,0
1- =γ
[ ]
1
12,75
88,93
59,98
=b [ ] [ ] [ ]
824,36
616,0
002,0
0382,0
== 1- bγλ 1=)(∑
)(
1=
u
u
n
α
αλ
996,537=616,0*500+002,0*400+0382,0*600=
)()(=)( ∑
)(
1=
*
u
uuu
n
α
α
OK
αOK ZλZ
[ ] [ ] 952,120==)(2 bu Tλσ
Exemplo 3.4. Fazer a krigagem das variáveis Ni e Co
#Carregar pacotes
library(geoR)
library(gstat)
library(lattice)
library(sp)
setwd("C:/PPGCA/GeoestatisticaR") #Definir a area de trabalho
Geo <- read.table("jura.txt", head=T) #carregar arquivo de dados
## lendo arquivo com bordas da área
jbor <- read.table("jurabor.txt", head=T) 
#Resumo do arquivo Geo
summary(Geo)
### Selecionar as variaveisl Cobalto e Niquel######
Ni <- as.geodata(Geo, coords.col = 1:2, data.col = 9)
Ni #mostrar a variavel Niquel area de trabalho
Co <- as.geodata(Geo,coords.col = 1:2, data.col = 6)
Co #mostrar a variavel Cobalto na area de trabalho
### Analise semivariografica: Niquel
variogNi <- variog(Ni, uvec = seq(0, 2, length = 20), 
option = "bin")
variogNi
#Semivariograma unidirecional experimental
plot(variogNi, xlab = "Distância", ylab = 
"Semivariância",type="l")
title("Semivariograma Experimental: Niquel")
#Semivariograma unidirecional e modelo ajustado
plot(variogNi, xlab = "Distância", ylab = "Semivariância")
lines.variomodel(cov.model = "spherical" , cov.pars = 
rbind(c(9,0.12),c(64,1.4)), nugget = 11, max.dist = 2, 
lwd = 2, col = "red")
text(1,0.1,"MODELO: 11 + 9,0 Sph(h/0,12) + 64 Sph(h/1,4)")
title("Semivariograma Experimental e modelo ajustado: 
Niquel")
###Mapa de Krigagem do Niquel
locNi <- expand.grid(seq(0,5.5,l=100), seq(0,6,l=120))
kcNi <- ksline(Ni, locations=locNi, cov.model = "spherical" , 
cov.pars = rbind(c(9,0.12),c(64,1.4)), nugget = 11)
contour(kcNi, filled=TRUE, bor=jbor, col=rainbow(22))
title("Niquel")
### Analise semivariografica: Cobalto
#Semivariogramas experimentais nas duas direções
coordinates(Geo) = ~ x + y
vgmCo <- variogram(Co~1, Geo, alpha=c(67.5,157.5))
plot(vgmCo,type="b", lty=2, col="black",main="Cobalto")
# Modelo anisotropico ajustado
model.Co <- vgm(11.5,"Sph",1.5,2.0,anis=c(67.5,0.67))
plot(vgmCo,model=model.Co,type="p",col=c("black"),main="Cobalto")
###Mapa de Krigagem do Cobalto
locCo <- expand.grid(seq(0,5.5,l=100), seq(0,6,l=120))
kcCo <- ksline(Co, locations=locCo, cov.model = "spherical",cov. 
pars=rbind(c(11.5,1.5)),nugget=2.0, aniso.pars = c(67.5,1.49))
contour(kcCo, filled=TRUE,borders=jbor, col=rainbow(30))
title("Cobalto")