[R-sig-Geo] gstat cross validation for accuracy ordinary kriging vs IDW

Alaba Boluwade alaba.boluwade at mail.mcgill.ca
Thu May 22 18:10:38 CEST 2014


Hi Moshood,

You did not follow what Michele sent to you:

ordkrigcv<- krige.cv(yield ~ 1, canmod.sp, 
model=exp.mod, nmax=20)


 try:

?krige.cv from the gstat library to see what  you need to supply for the function.


I would also advise you try reading  Bivand et al (2008). Some of  these difficulties you are encountering are pretty much solved there.

Bivand, R.S., E.J. Pebesma and V. Gómez- Rubio. 2008. Applied Spatial Data Analysis with R. Springer: New York, NY

Good luck.

Best Regards,

Alaba
________________________________________
From: r-sig-geo-bounces at r-project.org [r-sig-geo-bounces at r-project.org] on behalf of Moshood Agba Bakare [bakare at ualberta.ca]
Sent: Thursday, May 22, 2014 10:52 AM
To: Michele Fiori
Cc: r-sig-geo at r-project.org
Subject: Re: [R-sig-Geo] gstat cross validation for accuracy ordinary kriging vs IDW

Hi Michele,
I tried the script as advised but get the error message

ordkrigcv<- krige.cv(yield ~ 1, canmod.sp, newdata = grid,
model=exp.mod,nmax=20)Error in gstat(g = NULL, id = "var1", formula =
formula, data = locations,  :
  unused argument(s) (newdata = <S4 object of class "SpatialPixels">)


Likewise the same for the idw. Remember the yield data is irregularly
space and I created a regular grid as


grid <- expand.grid(easting=seq(from = 299678, to = 301299,
by=10),northing=seq(from = 5737278, to = 5738129, by=10))


## convert the grid to SpatialPixel class

coordinates(grid)<-~easting+northing
proj4string(grid)<-CRS("+proj=utm +zone=12 +ellps=WGS84 +datum=WGS84
+units=m +no_defs +towgs84=0,0,0")

gridded(grid)<- TRUE


I can not understand what the error message means but I suspect
newdata=grid which has now observed value to compare with the
interpolated value. I may be wrong.


Please what's next? Comments please and suggestions.


Moshood



On Thu, May 22, 2014 at 2:59 AM, Michele Fiori <mfiori at arpa.sardegna.it>wrote:

> Hi Moshood,
> I think you can do this way
>
> #### Ordinary kriging to create kriging prediction orkrig <- krige(yield ~
> 1, canmod.sp, newdata = grid, model=exp.mod,nmax=20)
>
>                 OrdKrigecv  <- Krige.cv(yield ~ 1, canmod.sp,
> model=exp.mod,nmax=20)
>                 RMSE.ok <- sqrt(sum(OrdKrigecv
> $residual^2)/length(OrdKrigecv $residual))
>
> ## Inverse Distance Weighting (IDW) Interpolation method maxdist=16.5
>
>                 idw1cv <- krige.cv(yield~1, canmod.sp, nmax=20,idp=1))
>                 RMSE.id <- sqrt(sum(idw1cv $residual^2)/length(idw1cv
> $residual))
>
> -----Messaggio originale-----
> Da: r-sig-geo-bounces at r-project.org [mailto:
> r-sig-geo-bounces at r-project.org]
> Per conto di Moshood Agba Bakare
> Inviato: mercoledì 21 maggio 2014 21:16
> A: r-sig-geo at r-project.org
> Oggetto: [R-sig-Geo] gstat cross validation for accuracy ordinary kriging
> vs
> IDW
>
> Hi all,
>
> I have been having a couple of challenge with my analysis. I have
> irregularly space spatial yield monitor data over four years. Pooling this
> data together is not feasible because of misalignment. That is, the
> coordinates of data point varies from one year to the other.
> I  created a common regular interpolation grid for each year with the same
> grid size of 10 x 10 m. I am able to get interpolated value for each point
> using ordinary kriging and inverse distance weighting method (IDW). Please
> how I cross validate this two interpolation methods to know which one give
> me the best estimate.
>
> The problem I notice is that there is no observe value in each
> interpolation
> point to assess the prediction accuracy of these methods.
> Please what do I do? see my script below. I correlated the interpolated
> values from the two methods. They are highly correlated (r=0.98). How do I
> know which  method gave good prediction?
>
> grid <- expand.grid(easting=seq(from = 299678, to = 301299, by=10),
>                     northing=seq(from = 5737278, to = 5738129, by=10))
>
> ## convert the grid to SpatialPixel class to indicate gridded spatial data
> coordinates(grid)<-~easting+northing
> proj4string(grid)<-CRS("+proj=utm +zone=12 +ellps=WGS84 +datum=WGS84
> +units=m +no_defs +towgs84=0,0,0")
>
> gridded(grid)<- TRUE
>
>
> #### Ordinary kriging to create kriging prediction orkrig <- krige(yield ~
> 1, canmod.sp, newdata = grid, model=exp.mod,nmax=20)
>
>
>
> ## Inverse Distance Weighting (IDW) Interpolation method maxdist=16.5
> idw1 = idw(yield~1, canmod.sp, newdata=grid,nmax=20,idp=1)
>
> Thanks while looking forward to reading from you.
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>

        [[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list