Calculate variogram of raster data from NA to R

Question

Summary: I have a raster dataset containing NA values, and you want to compute a semivariogram ignoring NA. How can i do this?

I have an image that I uploaded to R using the readGDAL function saved as im . To make this reproducible, the dput result in the image can be found at https://gist.github.com/2780792 . I am trying to show a variogram of this data and am afraid. I will consider what I have tried so far:

I tried the gstat package but could not get a function call that would work. I realized that basically I need the data values ​​themselves ( im@data $band1 ) and coordinates ( coordinates(im) ). I tried various commands like:

 > variogram(locations=coordinates(im), y = im@data $band1) Error in is.list(object) : 'object' is missing 

and

 > variogram(coordinates(im), im@data $band1) Error in variogram.default(coordinates(im), im@data $band1) : argument object and locations should be lists 

What am I doing wrong here?

Since this did not work, I tried the geoR package, which I called with:

 > variog(coords=coordinates(im), data=im@data $band1) variog: computing omnidirectional variogram Error in FUN(X[[1L]], ...) : NA/NaN/Inf in foreign function call (arg 4) 

The error looks like it is related to data containing NA, so I tried to delete it with na.omit , but that leaves all the NS there. This type makes sense since the raster file must have something in each square of the grid. Is there a way to remove NA somehow or at least make the variog command ignore them?

Any help would be greatly appreciated.

Answer 1

If the data object passed to gstat:variogram is a spatial object (your data is SpatialGridDataFrame ), you do not need to specify the locations as they are contained in the data.

However, obviously, NA values ​​cause a problem, so if we make the grid object be SpatialPointsDataFrame , this will remove the NA values

im contains the data https://gist.github.com/2780792

 library(gstat) point_data <- as(im, 'SpatialPointsDataFrame') gstat_variogram <- variogram(band1 ~ 1, data = point_data) 

To use geoR

 library(geoR) geor_variogram <- variog(coords = coordinates(point_data), data = point_data@data $band1) 

or even simpler (since geoR works with objects of the geodata class and contains the function as.geodata to convert from a SpatialPointsDataFrame object to geodata

 geor_variogram <- variog(as.geodata(point_data))