Leandro Valter
leandro.vvalter@gmail.com
Universidade Estadual da Paraíba
Centro de Ciência e Tecnologia
Departamento de Estatística
Pacotes necessários
require(akima)
require(geoR)
require(MASS)
require(moments)Importacao de dados para o R
setwd("C:/Users/Win7/Desktop/Estatística Espacial/aula_espacial")
solo= read.table('castro.txt', head=TRUE)
solo [1:2,]
#> Lat Lon X0.10 X10.15 X15.20 X20.25 X25.30 X30.35 X35.40
#> 1 608039.6 7250137 2.06 2.55 1.83 1.54 1.12 0.84 0.80
#> 2 608063.1 7250146 2.03 1.58 2.21 2.36 2.61 2.44 1.76
#> Xcd Ycd
#> 1 -49.93059 -24.85994
#> 2 -49.93036 -24.85986Foi dado início as análises com descritivas com o intuito de explorar os dados. As variáveis são as diferentes profundidades do solo.
Tabela 1: Estatísticas descritivas média, mínimo, máximo, desvio padrão, coeficiente de variação, curtose e assimetria para as profundidades.
| Profundidade(cm) | N | Média | Min. | Max. | D.P. | C.V.(%) | Curtose | Assimetria |
|---|---|---|---|---|---|---|---|---|
| 0-10 | 358 | 2,09 | 0,56 | 4,92 | 0,62 | 29,70 | 2,94 | 1,35 |
| 10-15 | 358 | 1,82 | 0,54 | 3,90 | 0,50 | 27,77 | 2,62 | 1,14 |
| 15-20 | 358 | 1,71 | 0,36 | 3,98 | 0,51 | 29,95 | 2,74 | 1,21 |
| 20-25 | 358 | 1,71 | 0,50 | 3,86 | 0,55 | 32,37 | 1,25 | 1,05 |
| 25-30 | 358 | 1,78 | 0,49 | 3,99 | 0,59 | 33,04 | 1,15 | 0,96 |
| 30-35 | 358 | 1,91 | 0,46 | 3,93 | 0,62 | 32,56 | 0,44 | 0,64 |
| 35-40 | 358 | 1,92 | 0,40 | 4,79 | 0,77 | 40,33 | 0,81 | 0,75 |
O coeficiente de variação de todas as profundidades foram maiores que 20% representando uma baixa homogeneidade. Todas as profundidades apresentam uma assimetria à esquerda (positiva), no ententato as profundidades 0-10,10-15,15-20 e 20-25 apresentam uma forte assimetria à esquerda, enquanto o restante apresenta uma fraca assimetria positiva. Verificando a curtose, todas as profundidades possuem uma função de distribuição leptocúrtica, ou seja, é menos achatada do que a distribuição normal. Provavelmente será necessário fazer o uso de algum artifício para corrigir as baixas homogeneidades, as assimetrias positivas e o achatamento menor que a da distribuição normal. Será utilizada a transformação Box-Cox, caso seja necessário.
Análise descritiva e gráfica profundidade 0-10
geosolo0_10 <- as.geodata(solo, coords.col=c(1, 2), data.col=3)
plot(geosolo0_10, low=T)
Verificação dos pressupostos
boxcox(geosolo0_10)
plot(geosolo0_10, low=T,lam=0)
points(geosolo0_10, pt.div="quintile", xlab="leste", ylab="norte")
Criação do polígono
setwd("C:/Users/Win7/Desktop/Estatística Espacial/aula_espacial")
bord=read.table("bord4.txt")
geosolo0_10$borders <- with(bord, cbind(V1,V2))
points(geosolo0_10, pt.div="quintile", xlab="leste", ylab="norte")
plot(geosolo0_10, low=T,lam=0)
Análise descritiva e gráfica profundidade 10-15
geosolo10_15 <- as.geodata(solo, coords.col=c(1, 2), data.col=4)
plot(geosolo10_15, low=T)
Verificação dos pressupostos
boxcox(geosolo10_15)
plot(geosolo10_15, low=T,lam=0)
#points(geosolo10_15, pt.div="quintile", xlab="leste", ylab="norte")Criação do polígono
geosolo10_15$borders <- with(bord, cbind(V1,V2))
points(geosolo10_15, pt.div="quintile", xlab="leste", ylab="norte")
plot(geosolo10_15, low=T,lam=0)
Análise descritiva e gráfica profundidade 15-20
geosolo15_20 <- as.geodata(solo, coords.col=c(1, 2), data.col=5)
plot(geosolo15_20, low=T)
Verificação dos pressupostos
boxcox(geosolo15_20)
plot(geosolo15_20, low=T,lam=0)
points(geosolo15_20, pt.div="quintile", xlab="leste", ylab="norte")
Criação do polígono
geosolo15_20$borders <- with(bord, cbind(V1,V2))
points(geosolo15_20, pt.div="quintile", xlab="leste", ylab="norte")
plot(geosolo15_20, low=T,lam=0)
Análise descritiva e gráfica profundidade 20-25
geosolo20_25 <- as.geodata(solo, coords.col=c(1, 2), data.col=6)
plot(geosolo20_25, low=T)
Verificação dos pressupostos
boxcox(geosolo20_25)
plot(geosolo20_25, low=T,lam=0)
points(geosolo20_25, pt.div="quintile", xlab="leste", ylab="norte")
Criação do polígono
geosolo20_25$borders <- with(bord, cbind(V1,V2))
points(geosolo20_25, pt.div="quintile", xlab="leste", ylab="norte")
plot(geosolo20_25, low=T,lam=0)
Análise descritiva e gráfica profundidade 25-30
geosolo25_30 <- as.geodata(solo, coords.col=c(1, 2), data.col=7)
plot(geosolo25_30, low=T)
Verificação dos pressupostos
boxcox(geosolo25_30)
plot(geosolo25_30, low=T,lam=0)
points(geosolo25_30, pt.div="quintile", xlab="leste", ylab="norte")
Criação do polígono
geosolo25_30$borders <- with(bord, cbind(V1,V2))
points(geosolo25_30, pt.div="quintile", xlab="leste", ylab="norte")
plot(geosolo25_30, low=T,lam=0)
Análise descritiva e gráfica profundidade 30-35
geosolo30_35 <- as.geodata(solo, coords.col=c(1, 2), data.col=8)
plot(geosolo30_35, low=T)
Verificação dos pressupostos
boxcox(geosolo30_35)
plot(geosolo30_35, low=T,lam=0)
points(geosolo30_35, pt.div="quintile", xlab="leste", ylab="norte")
Criação do polígono
geosolo30_35$borders <- with(bord, cbind(V1,V2))
points(geosolo30_35, pt.div="quintile", xlab="leste", ylab="norte")
plot(geosolo30_35, low=T,lam=0)
Análise descritiva e gráfica profundidade 35-40
geosolo35_40 <- as.geodata(solo, coords.col=c(1, 2), data.col=9)
plot(geosolo35_40, low=T)
Verificação dos pressupostos
boxcox(geosolo35_40)
plot(geosolo35_40, low=T,lam=0)
points(geosolo35_40, pt.div="quintile", xlab="leste", ylab="norte")
Criação do polígono
geosolo35_40$borders <- with(bord, cbind(V1,V2))
points(geosolo35_40, pt.div="quintile", xlab="leste", ylab="norte")
plot(geosolo35_40, low=T,lam=0)
Semivariograma
O método geoestatístico para diagnosticar a presença de correlação entre as unidades amostradas é o semivariograma, tem a função de caracterizar a estrutura de continuidade espacial de certa característica, por exigir hipóteses de estacionaridade menos restritivas. O semivariograma representa uma função de semivariâncias em relação às repectivas distâncias, sendo semivariância como a metade da variância de diferenças entre as observações de uma variável aleatória \(Z\), separadas por uma distância \(h\).
Os principais modelos de correlação são: Gaussiano, Esférico e a família das funções Matérn. A função Matérn com \(kappa\) igual a 0,5 é a função exponencial. As expressões destes modelos de correlação são apresentadas abaixo.
- Gaussiano: \[\rho(h)=\exp{\left(-\frac{h}{\phi}\right)^2}\]
- Esférico: \[\rho(h)=\left \{ \begin{array}{c} 1-1,5\left(-\frac{h}{\phi}\right)+0,5\left(-\frac{h}{\phi}\right)^3\\ 0,\ se \ h>0. \\ \end{array} \right. \]
- Exponencial: \[\rho(h)=\exp{\left(-\frac{h}{\phi}\right)}\] -Matérn: \[\rho(h)=\{ 2^{k-1} \Gamma{(K)} \}^{-1}\left(-\frac{h}{\phi}\right)^{k} K_k\left(-\frac{h}{\phi}\right)\]
Abaixo tem-se alguns semivariogramas teóricos para as profundidades utilizando a função Matern
Profundidade 0-10
v1=variog(geosolo0_10,max.dist = 250,lambda = 0)
#> variog: computing omnidirectional variogram
plot(v1,ylab=expression("Semivariância ("*matern*")"),
xlab="Distancia (m)", main = "Profundidade 0-10 (cm)")
#ef1=eyefit(v1)
#summary(ef1)
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 matern 0.03 435.29 0.06 0.41 <NA> 1199.23264011829
lines.variomodel(cov.model="matern", cov.pars=c(.03,435.29), nug=.06, max.dist=250,kappa =.41)
Profundidade 10-15
v2=variog(geosolo10_15,max.dist = 250,lambda = 0)
#> variog: computing omnidirectional variogram
plot(v2,ylab=expression("Semivariância ("*matern*")"),
xlab="Distancia (m)", main = "Profundidade 10-15 (cm)")
#ef2=eyefit(v2)
#summary(ef2)
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 matern 0.03 441.79 0.06 0.6 <NA> 1428.23261316673
lines.variomodel(cov.model="matern", cov.pars=c(.03,441.79), nug=.06, max.dist=250,kappa =.6)
Profundidade 15-20
v3=variog(geosolo15_20,max.dist = 250,lambda = 0)
#> variog: computing omnidirectional variogram
plot(v3,ylab=expression("Semivariância ("*matern*")"),
xlab="Distancia (m)", main = "Profundidade 15-20 (cm)")
#ef3=eyefit(v3)
#summary(ef3)
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 matern 0.03 175.42 0.06 0.41 <NA> 483.285602183215
lines.variomodel(cov.model="matern", cov.pars=c(.03,175.42), nug=.06, max.dist=250,kappa =.41)
Profundidade 20-25
v4=variog(geosolo20_25,max.dist = 250,lambda = 0)
#> variog: computing omnidirectional variogram
plot(v4,ylab=expression("Semivariância ("*matern*")"),
xlab="Distancia (m)", main = "Profundidade 20-25 (cm)")
#ef4=eyefit(v4)
#summary(ef4)
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 matern 0.05 272.87 0.06 0.27 <NA> 626.77887827009
lines.variomodel(cov.model="matern", cov.pars=c(.05,272.87), nug=.06, max.dist=250,kappa =.27)
Profundidade 25-30
v5=variog(geosolo25_30,max.dist = 250,lambda = 0)
#> variog: computing omnidirectional variogram
plot(v5,ylab=expression("Semivariância ("*matern*")"),
xlab="Distancia (m)", main = "Profundidade 25-30 (cm)")
#ef5=eyefit(v5)
#summary(ef5)
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 matern 0.05 110.45 0.06 0.41 <NA> 304.291955150124
lines.variomodel(cov.model="matern", cov.pars=c(.05,110.45), nug=.06, max.dist=250,kappa =.41)
Profundidade 30-35
v6=variog(geosolo30_35,max.dist = 250,lambda = 0)
#> variog: computing omnidirectional variogram
plot(v6,ylab=expression("Semivariância ("*matern*")"),
xlab="Distancia (m)", main = "Profundidade 30-45 (cm)")
#ef6=eyefit(v6)
#summary(ef6)
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 matern 0.06 175.42 0.07 0.4 <NA> 478.232229327521
lines.variomodel(cov.model="matern", cov.pars=c(.06,175.42), nug=.07, max.dist=250,kappa =.4)
Profundidade 35-40
v7=variog(geosolo35_40,max.dist = 250,lambda = 0)
#> variog: computing omnidirectional variogram
plot(v7,ylab=expression("Semivariância ("*matern*")"),
xlab="Distancia (m)", main = "Profundidade 35-40 (cm)")
#ef7=eyefit(v7)
#summary(ef7)
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
# 1 matern 0.09 72 0.12 1.89 <NA> 377.167271100791
lines.variomodel(cov.model="matern", cov.pars=c(.09,72), nug=.12, max.dist=250,kappa =1.89)
### Comprovando gráficamente a dependência espacial para cada profundidade
Comprovando gráficamente a dependencia espacial 0-10
par(mfrow=c(2,2))
var1=variog(geosolo0_10,max.dist = 250,lambda = 0);plot(var1)
#> variog: computing omnidirectional variogram
var2=variog(geosolo0_10,option = "cloud",lambda = 0);plot(var2)
#> variog: computing omnidirectional variogram
var3=variog(geosolo0_10,uvec=seq(0,350,l=50),lambda = 0)
#> variog: computing omnidirectional variogram
plot(var3)
env.var=variog.mc.env(geosolo0_10,obj.variog = var3,nsim = 99)
#> variog.env: generating 99 simulations by permutating data values
#> variog.env: computing the empirical variogram for the 99 simulations
#> variog.env: computing the envelops
plot(var3,env=env.var)
Comprovando gráficamente a dependencia espacial 10-15
par(mfrow=c(2,2))
var1=variog(geosolo10_15,max.dist = 350,lambda = 0);plot(var1)
#> variog: computing omnidirectional variogram
var2=variog(geosolo10_15,option = "cloud",lambda = 0);plot(var2)
#> variog: computing omnidirectional variogram
var3=variog(geosolo10_15,uvec=seq(0,350,l=50),lambda = 0);plot(var3)
#> variog: computing omnidirectional variogram
env.var=variog.mc.env(geosolo10_15,obj.variog = var3,nsim = 99)
#> variog.env: generating 99 simulations by permutating data values
#> variog.env: computing the empirical variogram for the 99 simulations
#> variog.env: computing the envelops
plot(var3,env=env.var)
Comprovando gráficamente a dependencia espacial 15-20
par(mfrow=c(2,2))
var1=variog(geosolo15_20,max.dist = 350,lambda = 0);plot(var1)
#> variog: computing omnidirectional variogram
var2=variog(geosolo15_20,option = "cloud",lambda = 0);plot(var2)
#> variog: computing omnidirectional variogram
var3=variog(geosolo15_20,uvec=seq(0,350,l=50),lambda = 0);plot(var3)
#> variog: computing omnidirectional variogram
env.var=variog.mc.env(geosolo15_20,obj.variog = var3,nsim = 99)
#> variog.env: generating 99 simulations by permutating data values
#> variog.env: computing the empirical variogram for the 99 simulations
#> variog.env: computing the envelops
plot(var3,env=env.var)
Comprovando gráficamente a dependencia espacial 20-25
par(mfrow=c(2,2))
var1=variog(geosolo20_25,max.dist = 350,lambda = 0);plot(var1)
#> variog: computing omnidirectional variogram
var2=variog(geosolo20_25,option = "cloud",lambda = 0);plot(var2)
#> variog: computing omnidirectional variogram
var3=variog(geosolo20_25,uvec=seq(0,350,l=50),lambda = 0);plot(var3)
#> variog: computing omnidirectional variogram
env.var=variog.mc.env(geosolo20_25,obj.variog = var3,nsim = 99)
#> variog.env: generating 99 simulations by permutating data values
#> variog.env: computing the empirical variogram for the 99 simulations
#> variog.env: computing the envelops
plot(var3,env=env.var)
Comprovando gráficamente a dependencia espacial 25-30
par(mfrow=c(2,2))
var1=variog(geosolo25_30,max.dist = 350,lambda = 0);plot(var1)
#> variog: computing omnidirectional variogram
var2=variog(geosolo25_30,option = "cloud",lambda = 0);plot(var2)
#> variog: computing omnidirectional variogram
var3=variog(geosolo25_30,uvec=seq(0,350,l=50),lambda = 0);plot(var3)
#> variog: computing omnidirectional variogram
env.var=variog.mc.env(geosolo25_30,obj.variog = var3,nsim = 99)
#> variog.env: generating 99 simulations by permutating data values
#> variog.env: computing the empirical variogram for the 99 simulations
#> variog.env: computing the envelops
plot(var3,env=env.var)
Comprovando gráficamente a dependencia espacial 30-35
par(mfrow=c(2,2))
var1=variog(geosolo30_35,max.dist = 350,lambda = 0);plot(var1)
#> variog: computing omnidirectional variogram
var2=variog(geosolo30_35,option = "cloud",lambda = 0);plot(var2)
#> variog: computing omnidirectional variogram
var3=variog(geosolo30_35,uvec=seq(0,350,l=50),lambda = 0);plot(var3)
#> variog: computing omnidirectional variogram
env.var=variog.mc.env(geosolo30_35,obj.variog = var3,nsim = 99)
#> variog.env: generating 99 simulations by permutating data values
#> variog.env: computing the empirical variogram for the 99 simulations
#> variog.env: computing the envelops
plot(var3,env=env.var)
Comprovando gráficamente a dependencia espacial 35-40
par(mfrow=c(2,2))
var1=variog(geosolo35_40,max.dist = 350,lambda = 0);plot(var1)
#> variog: computing omnidirectional variogram
var2=variog(geosolo35_40,option = "cloud",lambda = 0);plot(var2)
#> variog: computing omnidirectional variogram
var3=variog(geosolo35_40,uvec=seq(0,350,l=50),lambda = 0);plot(var3)
#> variog: computing omnidirectional variogram
env.var=variog.mc.env(geosolo35_40,obj.variog = var3,nsim = 99)
#> variog.env: generating 99 simulations by permutating data values
#> variog.env: computing the empirical variogram for the 99 simulations
#> variog.env: computing the envelops
plot(var3,env=env.var)
Modelos para cada profundidade
Dando continuidade precisa-se selecionar o melhor modelo geoestatstico univariado para as variaveis em estudo, que são as diferentes profundidades.
Tabela 2: Estimativas de máxima verossimilhança dos parâmetros associados ao modelo com diferentes funções de correlação e média constante sobre a profundidade 0-10 de estudo.
| Modelo | \(\hat\beta_0\) | \(\hat\beta_1\) | \(\hat\tau^2\) | \(\hat\sigma^2\) | \(\hat\tau^2 + \hat\sigma^2\) | \(\hat\phi\) | AIC | BIC |
|---|---|---|---|---|---|---|---|---|
| Matern(Kappa=0.5) | -3573,2609 | 0,0005 | 0,0519 | 0,0257 | 0,0776 | 32,6138 | 590,2562 | 609,6588 |
| Matern(Kappa=1) | -3581,2987 | 0,0005 | 0,0582 | 0,0196 | 0,0778 | 39,1672 | 590,0510 | 609,4537 |
| Pure Nugget | -2569,4299 | 0,0004 | 0,076 | 0 | 0,076 | 0 | 602,7292 | 622,1318 |
| Esférico | -3430,4184 | 0,0005 | 0,0586 | 0,0191 | 0,0777 | 100,6181 | 588,9096 | 608,3123 |
| Exponencial Exponenciada | -3546,7604 | 0,0005 | 0,0614 | 0,0165 | 0,0779 | 51,1099 | 589,493 | 608,8957 |
| Circular | -3218,4641 | 0,0004 | 0,0587 | 0,0182 | 0,0769 | 77,5151 | 589,7123 | 609,115 |
| Gassiano | -3554,2849 | 0,0005 | 0,0628 | 0,0152 | 0,078 | 55,1947 | 589,2088 | 608,6115 |
Para profundidade 0-10(cm) os dados apresentaram tendência na longitude, por isso ganharam um \(\hat\beta\) a mais. O melhor modelo pelos critérios de seleção foi o esférico, apresentando um menor AIC e BIC.
Tabela 3: Estimativas de máxima verossimilhança dos parâmetros associados ao modelo com diferentes funções de correlação e média constante sobre a profundidade 10-15 de estudo.
| Modelo | \(\hat\beta_0\) | \(\hat\beta_1\) | \(\hat\tau^2\) | \(\hat\sigma^2\) | \(\hat\tau^2 + \hat\sigma^2\) | \(\hat\phi\) | AIC | BIC |
|---|---|---|---|---|---|---|---|---|
| Matern(Kappa=0.5) | -3582,9738 | 0,0005 | 0,0636 | 0,0063 | 0,0699 | 78,1634 | 469,1076 | 488,5102 |
| Matern(Kappa=1) | -3647,4574 | 0,0005 | 0,0649 | 0,0052 | 0,0701 | 65,2931 | 468,9586 | 488,3612 |
| Pure Nugget | -2918,0984 | 0,0004 | 0,0692 | 0 | 0,0692 | 0 | 473,1436 | 492,5463 |
| Esférico | -3798,342 | 0,0005 | 0,0646 | 0,0058 | 0,0704 | 282,1833 | 468,0374 | 487,4400 |
| Exponencial Exponenciada | -3701,6614 | 0,0005 | 0,065 | 0,0051 | 0,0701 | 121,393 | 468,7966 | 488,1992 |
| Circular | -3387,911 | 0,0005 | 0,0645 | 0,005 | 0,0695 | 142,0205 | 469,4089 | 488,8116 |
| Gassiano | -3888,1953 | 0,0005 | 0,0658 | 0,0051 | 0,0709 | 167,6857 | 468,3937 | 487,7964 |
Para profundidade 10-15(cm) os dados novamente apresentaram tendência na longitude, por isso ganhou um \(\hat\beta\). O melhor modelo pelos critérios de seleção foi o esférico, apresentando um menor AIC e BIC.
Tabela 4: Estimativas de máxima verossimilhança dos parâmetros associados ao modelo com diferentes funções de correlação e média constante sobre a profundidade 15-20 de estudo.
| Modelo | \(\hat\beta\) | \(\hat\tau^2\) | \(\hat\sigma^2\) | \(\hat\tau^2 + \hat\sigma^2\) | \(\hat\phi\) | AIC | BIC |
|---|---|---|---|---|---|---|---|
| Matern(Kappa=0.5) | 0,5624 | 0,0621 | 0,0343 | 0,0964 | 100,16 | 469,8359 | 485,358 |
| Matern(Kappa=1) | 0,5529 | 0,0659 | 0,029 | 0,0949 | 59,2571 | 469,8644 | 488,3612 |
| Pure Nugget | 0,4943 | 0,0861 | 0 | 0,0861 | 0 | 473,1436 | 492,5463 |
| Esférico | 0,5205 | 0,061 | 0,0285 | 0,0895 | 123,9855 | 470,3366 | 485,8587 |
| Exponencial Exponenciada | 0,5306 | 0,0673 | 0,0248 | 0,0921 | 76,3526 | 469,8147 | 485,3368 |
| Circular | 0,5209 | 0,0607 | 0,0297 | 0,0904 | 110,1741 | 469,2672 | 484,7893 |
| Gassiano | 0,5306 | 0,0673 | 0,0248 | 0,0921 | 76,3514 | 469,8147 | 485,3368 |
Para profundidade 15-20(cm) os dados não apresentaram tendência. O melhor modelo pelos critérios de seleção foi o circular, apresentando um menor AIC e BIC.
Tabela 5: Estimativas de máxima verossimilhança dos parâmetros associados ao modelo com diferentes funções de correlação e média constante sobre a profundidade 20-25 de estudo.
| Modelo | \(\hat\beta\) | \(\hat\tau^2\) | \(\hat\sigma^2\) | \(\hat\tau^2 + \hat\sigma^2\) | \(\hat\phi\) | AIC | BIC |
|---|---|---|---|---|---|---|---|
| Matern(Kappa=0.5) | 0,556 | 0,0542 | 0,0536 | 0,1078 | 63,9157 | 491,3096 | 506,8317 |
| Matern(Kappa=1) | 0,5478 | 0,0624 | 0,0441 | 0,1065 | 41,7891 | 491,9068 | 507,429 |
| Pure Nugget | 0,4895 | 0,0968 | 0 | 0,0968 | 0 | 538,4259 | 553,948 |
| Esférico | 0,5254 | 0,0579 | 0,0454 | 0,1033 | 114,9353 | 493,2808 | 508,8029 |
| Exponencial Exponenciada | 0,56 | 0,0453 | 0,0627 | 0,108 | 53,8456 | 490,987 | 506,5091 |
| Circular | 0,5576 | 0,0645 | 0,0518 | 0,1163 | 170,4582 | 495,7832 | 511,3053 |
| Gassiano | 0,5306 | 0,0686 | 0,035 | 0,1036 | 65,2539 | 493,9222 | 509,4443 |
Para profundidade 20-25(cm) os dados não apresentaram tendência. O melhor modelo pelos critérios de seleção foi o Exponencial Exponenciada, apresentando um menor AIC e BIC.
Tabela 6: Estimativas de máxima verossimilhança dos parâmetros associados ao modelo com diferentes funções de correlação e média constante sobre a profundidade 25-30 de estudo.
| Modelo | \(\hat\beta\) | \(\hat\tau^2\) | \(\hat\sigma^2\) | \(\hat\tau^2 + \hat\sigma^2\) | \(\hat\phi\) | AIC | BIC |
|---|---|---|---|---|---|---|---|
| Matern(Kappa=0.5) | 0,6234 | 0,0651 | 0,0488 | 0,1139 | 86,9752 | 539,1666 | 554,6887 |
| Matern(Kappa=1) | 0,6155 | 0,0724 | 0,0401 | 0,1125 | 57,0457 | 539,9216 | 555,4437 |
| Pure Nugget | 0,528 | 0,1039 | 0 | 0,1039 | 0 | 591,5088 | 607,0309 |
| Esférico | 0,6042 | 0,0703 | 0,0429 | 0,1132 | 194,3463 | 540,1954 | 508,8029 |
| Exponencial Exponenciada | 0,6216 | 0,0705 | 0,0436 | 0,1141 | 97,2441 | 539,637 | 555,1591 |
| Circular | 0,6101 | 0,0711 | 0,0448 | 0,1159 | 181,5245 | 539,958 | 555,4802 |
| Gassiano | 0,5999 | 0,0788 | 0,0309 | 0,1097 | 95,6383 | 542,6214 | 558,1435 |
Para profundidade 25-30(cm) os dados não apresentaram tendência. O melhor modelo pelos critérios de seleção foi o Matern com o Kappa fixado em 0.5, ou seja, também poderia ser utilizada a função exponencial, apresentando um menor AIC e BIC.
Tabela 7: Estimativas de máxima verossimilhança dos parâmetros associados ao modelo com diferentes funções de correlação e média constante sobre a profundidade 30-35 de estudo.
| Modelo | \(\hat\beta\) | \(\hat\tau^2\) | \(\hat\sigma^2\) | \(\hat\tau^2 + \hat\sigma^2\) | \(\hat\phi\) | AIC | BIC |
|---|---|---|---|---|---|---|---|
| Matern(Kappa=0.5) | 0,9429 | 0,1295 | 0,1162 | 0,2457 | 202,4701 | 588,2491 | 603,7712 |
| Matern(Kappa=1) | 0,9561 | 0,1387 | 0,115 | 0,2537 | 122,6755 | 586,8588 | 602,381 |
| Pure Nugget | 0,7274 | 0,1997 | 0 | 0,1997 | 0 | 659,3784 | 674,9006 |
| Esférico | 0,8421 | 0,1305 | 0,0755 | 0,206 | 215,2553 | 587,7314 | 603,2536 |
| Exponencial Exponenciada | 0,9448 | 0,1389 | 0,1087 | 0,2476 | 181,03 | 586,3591 | 601,8812 |
| Circular | 0,8457 | 0,1298 | 0,0798 | 0,2096 | 189,2876 | 586,6159 | 602,1381 |
| Gassiano | 0,8862 | 0,1437 | 0,0799 | 0,2236 | 132,9214 | 586,5714 | 602,0935 |
Para profundidade 30-35(cm) os dados não apresentaram tendência. O melhor modelo pelos critérios de seleção foi o Exponencial Exponenciada, apresentando um menor AIC e BIC.
Tabela 8: Estimativas de máxima verossimilhança dos parâmetros associados ao modelo com diferentes funções de correlação e média constante sobre a profundidade 35-40 de estudo.
| Modelo | \(\hat\beta\) | \(\hat\tau^2\) | \(\hat\sigma^2\) | \(\hat\tau^2 + \hat\sigma^2\) | \(\hat\phi\) | AIC | BIC |
|---|---|---|---|---|---|---|---|
| Matern(Kappa=0.5) | 0,9441 | 0,1511 | 0,2286 | 0,3797 | 135,4849 | 697,6118 | 713,134 |
| Matern(Kappa=1) | 0,9194 | 0,1755 | 0,191 | 0,3665 | 78,2376 | 698,8909 | 714,413 |
| Pure Nugget | 0,7177 | 0,7177 | 0 | 0,7177 | 0 | 807,6301 | 823,1522 |
| Esférico | 0,8682 | 0,1536 | 0,2067 | 0,3603 | 207,9633 | 698,3469 | 713,869 |
| Exponencial Exponenciada | 0,9841 | 0,146 | 0,2717 | 0,4177 | 175,4164 | 697,532 | 713,0541 |
| Circular | 0,8876 | 0,1549 | 0,2281 | 0,383 | 195,4817 | 698,4772 | 713,9993 |
| Gassiano | 0,8785 | 0,1935 | 0,1588 | 0,3523 | 113,0547 | 702,1604 | 717,6825 |
Para profundidade 35-40(cm) os dados não apresentaram tendência. O melhor modelo pelos critérios de seleção foi o Exponencial Exponenciada, apresentando um menor AIC e BIC.
Krigagem e representação espacial
Profundidade 0-10
#gr=pred_grid(geosolo0_10$borders,by=5)
#points(geosolo0_10)
#points(gr,col=2,pch=19,cex=.3)
#gr0=gr[.geoR_inout(gr,geosolo0_10$borders),]
#points(gr0,col=3,pch=19,cex=0.3)
#KC=krige.control(type="OK",obj.model=ml4.2)
#kc1=krige.conv(geosolo0_10,loc=gr,krige=KC)
#image(kc1,main="Mapa de interpolação",col=terrain.colors(200),x.leg = c(607500,607800),y.leg =c(7250150,7250200))
Profundidade 10-15
Profundidade 15-20
Profundidade 20-25
Profundidade 25-30
### Profundidade 30-35
Profundidade 35-40
