[R-sig-eco] Kriging interpolation and plot

Rubén Roa rroa at azti.es
Wed Jun 15 09:43:42 CEST 2011


I think you want something like Fig. 4 in this paper of mine:

Roa-Ureta, R., and Niklitschek, E. 2007. Biomass estimation from surveys with likelihood-based geostatistics.  -  ICES Journal of Marine
Science, 64: 1723-1734.

There I put side by side a plot of the data from a scallop bed, with sizes and color of points proportional to scallop density, and the predicted density field (i.e. kriging). I used geoR.

Regarding your errors with geoR, I suggest you create a proper geodata object instead of passing the coordinates and the density as separate vectors.
See

?as.geodata

The error message you are receiving is not from krige.conv() but from rep() because your coords do not have the right number of rows. This is all caused by you passing the data erroneously. Thus my suggestion to create a proper geodata object, that will solve your problem and will help you understand your spatial prediction.

If you are bent on using geoR (i.e. you want to do spatial modeling using pure-likelihood or Bayesian inference) I can send you a well-commented script with my approach to spatial modeling as shown in the paper and applications above. The script includes how to make my Fig. 4 which is what apparently you want to do with your data.

HTH

Rubén

____________________________________________________________________________________ 

Dr. Rubén Roa-Ureta
AZTI - Tecnalia / Marine Research Unit
Txatxarramendi Ugartea z/g
48395 Sukarrieta (Bizkaia)
SPAIN



> -----Mensaje original-----
> De: r-sig-ecology-bounces at r-project.org 
> [mailto:r-sig-ecology-bounces at r-project.org] En nombre de 
> Xurde Sanchez
> Enviado el: martes, 14 de junio de 2011 13:30
> Para: r-sig-ecology at r-project.org
> Asunto: [R-sig-eco] Kriging interpolation and plot
> 
> Dear all,
> 
>  I would like to perform a plot with kriging interpolation. 
> My datas are very simple. I have two distances (x,y) measured 
> in meters, which are distance to a river mouth and depth and 
> in those points we measured many parameter as Temperature, 
> Conductivity, ...
> 
> So our data are something like this:
> 
> x          y          Temp
> 
> 0          1          20.5
> 
> 0          2          20.5
> 
> 0          3          20.2
> 
> 0          4          19.5
> 
> 1          1          20.5
> 
> 1          2          20.5
> 
> 1          3          20.2
> 
> 1          4          19.5
> 
> 
> 
> We want to perform a plot with those data and with a kriging 
> interpolation, but we are not able to finish it properly. We 
> have tried *gstat* and *geoR*packages without success. In the 
> last case, geoR, we almost got it.
> 
> We used those scripts:
> 
> 
> 
> library(geoR)
> 
> 
> 
> data02<-read.table("Coord&Temp.txt",header=TRUE, sep="\t", 
> na.strings="NA", dec=".", strip.white=TRUE, row.names=1)
> 
> attach(data02)
> 
> names(data02)
> 
> summary(data02)
> 
> 
> 
> datacoord<-read.table("coords.txt", header=TRUE, sep="\t", 
> na.strings="NA", dec=".", strip.white=TRUE)
> 
> attach(datacoord)
> 
> names(datacoord)
> 
> summary(datacoord)
> 
> datacoord
> 
> datacoord<-c(datacoord)
> 
> datacoord
> 
> 
> 
> vdata<-c(Temp)
> 
> vdata
> 
> 
> 
> 
> 
> plot(xcoord,ycoord, xlab="X-Coordinate", ylab="Y-Coordinate")
> 
> 
> 
> j<-expand.grid (seq (0.00,35.50,1), seq(-16.000,0.000,1))
> 
> 
> 
> j
> 
> j[1:612,]
> 
> names(j)<-c("xcoord","ycoord")
> 
> j[1:612,]
> 
> 
> 
> 
> 
> plot(xcoord,ycoord, xlab="X-Coordinate", ylab="Y-Coordiante")
> 
> 
> 
> addgridtoplot(j) #it doesn't work
> 
> 
> 
> predgrid<-k
> 
> 
> 
> k<-krige.control(type.krige="ok", cov.model= "spherical",
> 
>             cov.pars=c(12920.1376,79505.96), nugget=704.3657)
> 
> k
> 
> 
> 
> k2<-krige.conv(vdata, locations=predgrid, borders=NULL, krige=k)
> 
> 
> 
> I must be doing something wrong , because  when I execute the 
> last command R provides me an error:
> 
> 
> 
> k2<-krige.conv(vdata, locations=predgrid, borders=NULL, krige=k)
> 
> krige.conv: model with constant mean
> 
> Error en rep(1, nrow(geodata$coords)) : argumento 'times' inválido
> 
> Además: Mensajes de aviso perdidos
> 
> In krige.conv(datacoord, vdata, locations = predgrid, krige = j) :
> 
>   krige.conv: coordinates provided as a vector, assuming one 
> spatial dimension
> 
> 
> 
> Which are my mistakes? I think it is relate with my data but 
> in geoR instructions we can read that if you don't have 
> geodata then  coords and data must be provided instead.¿?¿?¿?
> 
> 
> 
> Many thanks in advance for your help and for your answers!
> 
> 
> 
> Best regards,
> 
> 
> 
> 
> 
> Jorge R. Sanchez
> 
> 	[[alternative HTML version deleted]]
> 
> 



More information about the R-sig-ecology mailing list