############################################################################# ## Análise dos teores de argila usando as variáveis de solo (MO e CTC) ## ############################################################################# ## Pacotes require(geoR) require(splancs) ### Lendo arquivo como geodata sensor <- read.geodata("sensor.txt", head=T, coords.col=c(1,2), data.col=3) biomassa <- read.geodata("amostras.txt", head=T, coords.col=c(1,2), data.col=3) argila <- read.geodata("amostras.txt", head=T, coords.col=c(1,2), data.col=4) #### adicionando contorno #### bor <- read.table("cont.txt", head=F) bor # apresenta os dados da borda plot(bor) polygon(bor) apply(bor,2,range) #Mostra o mínimo e máximo das coordenadas gr<-expand.grid(x=seq( 192065.2,192354.9,by=4), y=seq(7634138,7634613, by=4)) plot(gr) points(sensor, pt.div="equal") #monta o grid de interpolação gi<- polygrid(gr,bor=bor) points(gi, pch="+", col=2) #o novo grid considerando apenas a região limitada gi ###Gráficos Post Plot para avaliar a tendencia direcional #sensor points(sensor,l=1,pt.div="equal",col=gray(seq(1,0,l=5))) polygon(bor) #biomassa points(biomassa,l=1,pt.div="equal",col=gray(seq(1,0,l=5))) polygon(bor) #argila points(argila,l=1,pt.div="equal",col=gray(seq(1,0,l=5))) polygon(bor) ### variograma chute inicial e maxima verossimilhança em 3 modelos ### # sensor sensor.var = variog(sensor, max.dist=50) plot(sensor.var) env.sensor <-variog.mc.env(sensor, obj.v=sensor.var) plot(sensor.var,env=env.sensor,xlab="distância", ylab="semivariância",font.main = 3) exp.sensor = likfit(sensor, ini=c(0.0014,20), nug=0.001, cov.model="exp") exp.sensor sph.sensor = likfit(sensor, ini=c(0.0014,20), nug=0.001, cov.model="sph") sph.sensor # biomassa biomassa.var = variog(biomassa, max.dist=400) plot(biomassa.var) env.biomassa <-variog.mc.env(biomassa, obj.v=biomassa.var) plot(biomassa.var,env=env.biomassa,xlab="distância", ylab="semivariância",font.main = 3) exp.biomassa = likfit(biomassa, ini=c(13,100), nug=0,cov.model="exp") exp.biomassa sph.biomassa = likfit(biomassa, ini=c(13,100), nug=5,cov.model="sph") sph.biomassa # argila argila.var = variog(argila, max.dist=250) plot(argila.var) env.argila <-variog.mc.env(argila, obj.v=argila.var) plot(argila.var,env=env.argila,xlab="distância", ylab="semivariância",font.main = 3) exp.argila = likfit(argila, ini=c(50,150), nug=5,cov.model="exp") exp.argila sph.argila = likfit(argila, ini=c(50,150), nug=5,cov.model="sph") sph.argila #### Validação cruzada para modelos ajustados por MV #### ## sensor vc.sensor.exp=xvalid(sensor, model=exp.sensor) erro.sensor.exp=sum(abs(vc.sensor.exp$predicted-vc.sensor.exp$data)) erro.sensor.exp summary(vc.sensor.exp) vc.sensor.sph=xvalid(sensor, model=sph.sensor) erro.sensor.sph=sum(abs(vc.sensor.sph$predicted-vc.sensor.sph$data)) erro.sensor.sph summary(vc.sensor.sph) ## biomassa vc.biomassa.exp=xvalid(biomassa, model=exp.biomassa) # maxver não ajustou parâmetros erro.biomassa.exp=sum(abs(vc.biomassa.exp$predicted-vc.biomassa.exp$data)) erro.biomassa.exp summary(vc.biomassa.exp) vc.biomassa.sph=xvalid(biomassa, model=sph.biomassa) erro.biomassa.sph=sum(abs(vc.biomassa.sph$predicted-vc.biomassa.sph$data)) erro.biomassa.sph summary(vc.biomassa.sph) ## argila vc.argila.exp=xvalid(argila, model=exp.argila) erro.argila.exp=sum(abs(vc.argila.exp$predicted-vc.argila.exp$data)) erro.argila.exp summary(vc.argila.exp) vc.argila.sph=xvalid(argila, model=sph.argila) erro.argila.sph=sum(abs(vc.argila.sph$predicted-vc.argila.sph$data)) erro.argila.sph summary(vc.argila.sph) #### Mapas individuais ###### # pelos critério analisados na validação cruzada, os modelos que mais se ajustaram foram: # sensor: sph # biomassa: sph # argila: sph ## sensor c.sensor=krige.control(obj=sph.sensor) kg.sensor=krige.conv(sensor, loc=gi, krige=c.sensor) image(kg.sensor, loc=gr, border=bor, col=gray(seq(1,0,l=5)),x.leg=c(192140, 192440), y.leg=c(7634150, 7634200),zlim=range(kg.sensor$predict)) ## biomassa c.biomassa=krige.control(obj=sph.biomassa) kg.biomassa=krige.conv(biomassa, loc=gi, krige=c.biomassa) image(kg.biomassa, loc=gr, border=bor, col=gray(seq(1,0,l=5)),x.leg=c(192140, 192440), y.leg=c(7634150, 7634200),zlim=range(kg.biomassa$predict)) # obs: avaliando os valores krigados e o baixo ajuste da MV, verifica-se que # não foi possível moodelar eficientemente a variável "biomassa". ## argila c.argila=krige.control(obj=sph.argila) kg.argila=krige.conv(argila, loc=gi, krige=c.argila) image(kg.argila, loc=gr, border=bor, col=gray(seq(1,0,l=5)),x.leg=c(192140, 192440), y.leg=c(7634150, 7634200),zlim=range(kg.argila$predict))