Dear Group, Happy New Year to Everybody.
I am a novice user of gstat and geostatistics. I have a data set (please see attached: data.csv) containing coordinates (East_x, North_y),Elev(Elevation), and Richness (Plant Species) collected from 125 plots of 10m x 10m plots. I want to interpolate plant species richness from the relationship between Richness and Elevation of the study area using Regression Krigging. I have written the following codes: #Variogram of residuals #=========================================================================== data<-read.csv("E://data.csv",header=TRUE) attach(data) names(data) m2=glm(Richness~Elev) summary(m2) residuals<-m2$residuals data$Residuals<-residuals library(maptools) coordinates(data) = ~ East_x+North_y library(gstat) var <- variogram(Residuals~1, data=data,alpha=c(0,22.5,45,67.5,90,112.5,135), tol.hor=30,cutoff=1000,width=67) plot(var) #Variogram shows that major axis is 45 & minor 112.5; Fitting variogram in #major axis var.45<-var$dir.hor==45 plot(var[var.45,]) eye.sph<-vgm(psill=10.5,model="Sph",range=705,nugget=2) plot(var[var.45,],eye.sph) #Fit anisotropic variogram fit.sph<-fit.variogram(var,model=vgm(psill=10.5,model="Sph",range=705, nugget=2,anis=c(45,0.72))) plot(var, model=fit.sph, as.table=TRUE) fit.exp<-fit.variogram(var, model=vgm(psill=10.5,model='Exp',range=705/3, nugget=2, anis=c(45, 0.72))) plot(var, model=fit.exp, as.table=TRUE) fit.gau<-fit.variogram(var,model=vgm(psill=10.5,model='Exp',range=705/sqrt(3), nugget=2, anis=c(45, 0.72))) plot(var, model=fit.gau, as.table=TRUE) #Cross Variogram cv.sph<-krige.cv(Residuals~1,data,fit.sph) View(cv.sph) mean(cv.sph$zscore) #-0.0001968425 var(cv.sph$zscore) #1.114187 cv.exp<-krige.cv(Residuals~1,data,fit.exp) View(cv.exp) mean(cv.exp$zscore) #-1.243423e-05 var(cv.exp$zscore) #1.116687 cv.gau<-krige.cv(Residuals~1,data,fit.gau) View(cv.gau) mean(cv.gau$zscore) #-1.235693e-05 var(cv.gau$zscore) #1.116684 #Cross variogram shows that spherical model shows the best fit. The sill varies in different direction of the varigram. So, we will have to fit a zonal #anisotropic variogram with spherical model (fit.sph). Variogram shows that #22.5 deg has the highest sill and 112.5 deg should have the lowest sill. We #will fit the zonal variogram to the 112.5 deg direction with an anisotropic #ratio of 293/(1600*10000) #Fitting zonal anisotropic variogram fit.sph1<-vgm(psill=2.5,"Sph",range=1600*10000,anis=c(122.5,1.825e-05), add.to=fit.sph) plot(var,fit.sph1,as.table=TRUE) #Cross Variogram of the final model cv.sph1<-krige.cv(Residuals~1,data,fit.sph1) View(cv.sph1) mean(cv.sph1$zscore) #-0.0009952294 var(cv.sph1$zscore) #1.071403 #Creating a grid file #============================================================================ grd <- expand.grid(x=seq(from=341960,to=343195,by=10), y=seq(from=2667253,to=2668873,by=10)) library(maptools) coordinates(grd) <- ~ x+y gridded(grd) <- TRUE So far,I think, things were going well but I found a error while running the following code #Regression kriggin #=========================================================================== prk.rich<-krige(Richness~Elev, data=data, newdata=grd,model=fit.sph1) Error in function (classes, fdef, mtable) : unable to find an inherited method for function krige for signature "formula", "missing" Can anyone help me out of this please? I shall be grateful. Best regards, --------------- Md. Abdul Halim Assistant Professor Department of Forestry and Environmental Science Shahjalal University of Science and Technology,Sylhet-3114, Bangladesh. Cell: +8801714078386. alt. e-mail: xo...@yahoo.com -- This message has been scanned for viruses and dangerous content by MailScanner, and is believed to be clean.
data.csv
Description: MS-Excel spreadsheet
_______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo