Re: creating a variogram with gstat in R

["Edzer J. Pebesma" <e.pebesma@geog.uu.nl>]

femke wrote:

>  
> Hello again,
>  
> Sorry to be dense about this but I'm still having trouble working out 
> how gstat operates within R in terms of data import (the data import 
> phase isn't in the examples in the demos).

OK -- maybe I aimed more at people familiar with R than
people familiar with gstat-standalone, sorry for that.

>  
> I import my data fine - and have been trying it with two approaches, 
> one with the original ascii grid (dataA) and the other with a csv list 
> of x,y, and value (dataB).
>  
> > dataA <- read.table(file="testData.txt", sep=" ", na.strings = 
> "-9999", skip = 6)

I wouldn't expect success here: asciigrid files only contain z values,
not x or y values (as they're implicitly defined by row,col, grid origin
and cell size values).

If you want to create a grid (without missing values), explore the
S function expand.grid, e.g.

expand.grid(x=1:10,y=1:10)

You could get the wanted table by:

predictionlocs = data.frame(expand.grid(x=..., y=...), dataA)

which glues x and y (with ... to be filled in correctly) to the
attribute values (and possible NA's) in dataA.

> > dataB <- read.table("testData2.csv", sep=",", header=TRUE,  row.names=1)
>  

I hade more success when using read.cvs, but maybe this works.
I don't think the row.names=1 does anything.

> In attempting to plot a basic variogram, I do the following and get 
> the an error message I can't interpret (nor is there much help on line 
> with this error message):
>  
> > v <- variogram(log(A1)~1,loc= ~x+y, dataB)
> Error in vector("double", length) : negative length vectors are not 
> allowed
> I'm also wondering how the variogram function knows which column is x 
> and which y (as is specified in the command files), or is the x and y 
> coded in?

This is in the manual pages. The loc=~x+y tells the variogram
function that locations are in columns named "x" and "y" in dataB.

>  
> With dataA (the asciigrid) I'm not even sure how to create a variogram 
> as it doesn't have a column header.
>  

Maybe you can export the grid from ArcGIS in a column format. In R,
it is certainly possible, e.g. using the expand.grid() function, but a
more direct way (a function reading the asciigrid header, calling
expand.grid directly with the right parameters) I wouldn't know.

You could post a question on this subject on R-sig-geo, see

https://www.stat.math.ethz.ch/mailman/listinfo/r-sig-geo

> Please help.
>  
> Thanks in advance, and apologies for the simplistic questions.
>  

No, thank you for the questions, from which we, developers,
learn from you, users, where you get stuck.
--
Edzer