## lendo os dados wolf <- read.table("wolfcamp.txt", head=T, row.names=1) head(wolf) dim(wolf) colnames(wolf) rownames(wolf) dimnames(wolf) ## carregando o pacote geoR require(geoR) ## analise exploratoria wolfg <- as.geodata(wolf) args(as.geodata) args(as.geodata.default) summary(wolfg) plot(wolfg) plot(wolfg, low=T) plot(wolfg, low=T, trend="1st") points(wolfg) points(wolfg, trend="1st") points(wolfg, trend="1st", pt.div="quint") args(points.geodata) ##x11() ##points(wolfg, trend="1st", pt.div="quint", perm=T) wolfv <- variog(wolfg) plot(wolfv) wolfv <- variog(wolfg, max.dist=220) plot(wolfv) wolfv1 <- variog(wolfg, trend="1st", max.dist=220) plot(wolfv1) wolfv1.env <- variog.mc.env(wolfg, obj=wolfv1) lines(wolfv1.env) plot(wolfv1, env=wolfv1.env) wolfvfit <- eyefit(wolfv1) ## Max. Vero. wolfml0 <- likfit(wolfg, ini=c(200000,200)) wolfml0 wolfml1 <- likfit(wolfg, trend="1st", ini=c(3000,30), nug=1000) wolfml1 wolfml2 <- likfit(wolfg, trend="2nd", ini=c(3000,30), nug=1000) wolfml2 logLik(wolfml0) logLik(wolfml1) logLik(wolfml2) pchisq(8, df=3, low=F) summary(wolfml1) summary(wolfml2) wolfml1a <- likfit(wolfg, trend="1st", cov.model="sph", ini=c(3000,100), nug=1000) wolfml1a wolfml1b <- likfit(wolfg, trend="1st", cov.model="mat", kapp=1, ini=c(3000,100), nug=1000) wolfml1b wolfml1c <- likfit(wolfg, trend="1st", cov.model="mat", kap=1.5,ini=c(3000,100), nug=1000) wolfml1c ## cuidado -- a lina abaixo NAO PRECISA SER UM AJUSTE PERFEITO AOS PONTOS! lines(wolfml1, col=2) ## vmAOIS ESCOLHER O wolfml1a ## definindo borda arbitrária points(wolfg) bor <- locator(ty="o") bor <- matrix(unlist(bor), nc=2) bor <- rbind(bor,bor[1,]) polygon(bor) polygon(bor, col=2) wolfg$borders <- bor # interpolacao bbox <- apply(bor,2,range) bbox gr0 <- pred_grid(bor, by=5) points(wolfg) points(gr0, pch=19, cex=0.25, col=2) gr <- locations.inside(gr0, bor) points(gr, pch=19, cex=0.24,col=4) kc <- krige.conv(wolfg, loc=gr0, bor=bor, krige=krige.control(obj=wolfml1)) names(kc) kc$pred[1:10] image(kc) image(kc, col=gray(seq(1,0.2,len=21)), x.leg=c(-220,-50), y.leg=c(150,180)) kc2 <- krige.conv(wolfg, loc=gr0, bor=bor, krige=krige.control(obj=wolfml2)) image(kc2, col=gray(seq(1,0.2,len=21)), x.leg=c(-220,-50), y.leg=c(150,180)) plot(kc$pred, kc2$pred, asp=1) abline(0,1, col=2) summary(kc$pred - kc2$pred) ##