[R] Plot SpatialLinesDataFrame with xlim ylim

[Edzer J. Pebesma]

Michael,

The plot method for SpatialLinesDataFrame objects resides in package sp, 
and questions regarding it are easier noticed on the r-sig-geo mailing 
list.

The reason why they are plotted with aspect ratio 1 is that they are 
assumed to be spatial (geographical) data, and assume that 1 m north 
equals 1 m west -- think of a map. The exception is when the projection 
argument is set to longlat data (i.e. decimal degrees North/East), where 
the aspect ratio is computed differently, such that the argument above 
more or less holds.

You should be able to override the default aspect setting by explicitly 
passing the e.g. asp=0.5 argument to plot.

Here's the comment in the documentation of plot for Spatial objects 
(such as SpatialLinesDataFrame):

The default aspect for map plots is 1; if however data are not projected 
(coordinates are longlat), the aspect is by default set to 1/cos(My * 
pi)/180) with My the y coordinate of the middle of the map (the mean of 
ylim, which defaults to the y range of bounding box).

The argument |setParUsrBB| may be used to pass the logical value |TRUE| 
to functions within |plot.Spatial|. When set to |TRUE|, par(“usr”) will 
be overwritten with |c(xlim, ylim)|, which defaults to the bounding box 
of the spatial object. This is only needed in the particular context of 
graphic output to a specified device with given width and height, to be 
matched to the spatial object, when using par(“xaxs”) and par(“yaxs”) in 
addition to |par(mar=c(0,0,0,0))|.
--

Edzer

I'm running windows xp, R 2.3.1 with maptools 0.6-6, I guess. 
When plotting from a large SpatialLinesDataFrame and using xlim  ylim to 
reduce the area, the plot axes automatically have the same scale size, even if 
xlim and ylim ranges differ.  
E.g.:
tmp - readShapeLines(filepath)
plot(tmp,xlim=c(-126,-119),ylim=c(50,51))

The y-axis range is actually 47-54, same range as the x-axis.  What am I doing 
wrong?  Should I be using a different object for simple coastline  river data?
Thanks in advance!
Michael

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Lagrange multiplier for ordinary kriging in gstat in R??

[Edzer J. Pebesma]

Zia, 

this is not trivial, as gstat uses the kriging solution where 
Lagrange parameters are not involved. You can get the separate
components of it and work from there. Another angle would be the
difference between the ordinary and simple kriging variance, but
I'm not 100% sure that equals the Lagrange parameter. Perhaps need
to add a minus sign.

btw -- may I suggest to use r-sig-geo for these kind of questions.
--
Edzer



Does anybody know how to write the code for extracting value of  LAGRANGE
MULTIPLIER after ordinary prediction in gstat in R?

Thanks!

Zia


-- 
Zia Uddin Ahmed

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] [R-sig-Geo] plot methods in sp

[Edzer J. Pebesma]

Patrick, I can see this is confusing.

Besides methods(plot), try after library(sp) a

showMethods(plot)

to get an overview of the S4-style plot methods; methods(plot) only 
shows the S3-style plot methods. Then,

class?SpatialPolygons

gives methods available for this class, among which plot, and the html 
listing of plot methods has entries like

plot.SpatialPolygons,missing-method

that points to the same help page (which is rather brief).

getMethod(plot, c(SpatialPolygons, missing))

shows you the arguments this method takes, and shows that it actually 
calls plot.SpatialPolygons. This method is not exported from sp, so to 
view it you need to access it as

 sp:::plot.SpatialPolygons
 
the other sp functions it calls can also be accessed by preceding them 
with sp:::

I agree that more elaborate documentation along with examples would be 
useful.

Hope this helps,
--
Edzer

Patrick Giraudoux wrote:
 Dear listers,

 I am working since a while with the sp package and still wonder how the 
 plot methods are managed with sp spatial objects. For instance, 
 SpatialPolygonsDataFrame objects have obviously a plot method. However 
 it cannot be found in the list provided by methods(plot) . Furthermore 
 ?plot.SpatialPolygonsDataFrame, nor ?plot.SpatialPolygons, etc.. provide 
 a help, though the lattice function spplot is adequately documented.

 On the other hand, plot(myobject, border=grey), with myobject a 
 SpatialPolygonsDataframe is well interpreted and recalls the syntax of 
 plot.polylist of matools (though myobject is far from being a polylist...).

 Can anybody (especially the package's authors...) comment on this? Where 
 a help with the list of the plot function arguments can be found?

 Thanks for any hint,

 Patrick

 ___
 R-sig-Geo mailing list
 R-sig-Geo@stat.math.ethz.ch
 https://stat.ethz.ch/mailman/listinfo/r-sig-geo


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] [R-sig-Geo] plot() and Jpeg() increase font size and resolution

[Edzer J. Pebesma]

Herry, you need to set font size after opening a device; try after your 
example e.g.

png(out.png)
trellis.par.set(fontsize,fontsize)
plt1
dev.off()

Hth,
--
Edzer



[EMAIL PROTECTED] wrote:
 Hi Edzer,

 png also does not increase the file size. 

 Following the commands. The data are located at
 ftp://ftp.csiro.au/Herry, file qldproperty.Rdata.

 
 trellis.par.get(fontsize)-fontsize
  fontsize$default-16
  fontsize$points-16
  fontsize$text-18
   trellis.par.set(fontsize,fontsize)
 trellis.par.get(fontsize)
  
  
  
 plot(v1,plot.number=F, model=m2, ylim=c(4,9))

 plt1-plot(v1,plot.number=F, model=m1, ylim=c(4,9), col=black,
 pch=+, scales=list(cex=1.5,)) 

 jpeg(file=LOTPLAN_variogram_mod.png, bg=white, res=300,
 pointsize=16, width=800, height=800)
  plt1
 dev.off()
 XXX


 SessionInfo:
 R version 2.4.0 (2006-10-03)
 i686-pc-linux-gnu

 locale:
 LC_CTYPE=en_US;LC_NUMERIC=C;LC_TIME=en_US;LC_COLLATE=en_US;LC_MONETARY=e
 n_US;LC_MESSAGES=en_US;LC_PAPER=en_US;LC_NAME=C;LC_ADDRESS=C;LC_TELEPHON
 E=C;LC_MEASUREMENT=en_US;LC_IDENTIFICATION=C

 attached base packages:
 [1] methods   stats graphics  grDevices utils
 datasets
 [7] base

 other attached packages:
 geoRgstat MCMCpack MASS coda
 lattice
 1.6-11 0.9-34  0.7-4 7.2-29 0.10-7
 0.14-13
 mcmc RColorBrewer maps maptools   sp
 foreign
  0.5-1  0.2-3 2.0-32  0.6-3  0.9-4
 0.8-17
  spatial
 7.2-29


 Thanks
 Herry


 Dr Alexander Herr
 Spatial and statistical analyst
 CSIRO, Sustainable Ecosystems
 Davies Laboratory,
 University Drive, Douglas, QLD 4814 
 Private Mail Bag, Aitkenvale, QLD 4814
  
 Phone/www 
 (07) 4753 8510; 4753 8650(fax)
 Home: http://herry.ausbats.org.au
 Webadmin ABS: http://ausbats.org.au
 Sustainable Ecosystems: http://www.cse.csiro.au/
 


 -Original Message-
 From: Edzer J. Pebesma [mailto:[EMAIL PROTECTED] 
 Sent: Monday, November 27, 2006 5:34 PM
 To: Herr, Alexander Herr - Herry (CSE, Townsville)
 Cc: r-help@stat.math.ethz.ch; r-sig-geo@stat.math.ethz.ch
 Subject: Re: [R-sig-Geo] plot() and Jpeg() increase font size and
 resolution

 [EMAIL PROTECTED] wrote:
   
 Thanks to Edzer and Roger,
 I can now plot with increased font sizes. However, jpeg still does not
 

   
 reproduce these, nor does it show up in high quality. What I would 
 like to do is produce some highresolution jpegs.
   
 
 Please specify the commands that you used, preferrably such that we can
 reproduce them, and specify which worked and which didn't the way you
 expected. Also please give the output of sessionInfo(), and tell us the
 platform you work on. Have you tried producing png, did that work?
 --
 Edzer
   
 Any help would be appreciated

 Thanx
 Herry


 R2.4 on Mandriva 10.2 linux.

 Dr Alexander Herr
 Spatial and statistical analyst
 CSIRO, Sustainable Ecosystems
 Davies Laboratory,
 University Drive, Douglas, QLD 4814
 Private Mail Bag, Aitkenvale, QLD 4814
  
 Phone/www
 (07) 4753 8510; 4753 8650(fax)
 Home: http://herry.ausbats.org.au
 Webadmin ABS: http://ausbats.org.au
 Sustainable Ecosystems: http://www.cse.csiro.au/
 


 -Original Message-
 From: Edzer J. Pebesma [mailto:[EMAIL PROTECTED]
 Sent: Friday, November 24, 2006 5:37 PM
 To: Herr, Alexander Herr - Herry (CSE, Townsville)
 Cc: r-help@stat.math.ethz.ch; r-sig-geo@stat.math.ethz.ch
 Subject: Re: [R-sig-Geo] plot() and Jpeg() increase font size and 
 resolution

 plotting variograms in gstat is done through xyplot in lattice; you'll
 

   
 find where it gets it's defaults by

 library(lattice)
 trellis.par.get()
 ?trellis.par.set
 --
 Edzer

 [EMAIL PROTECTED] wrote:
   
 
  
 Dear list,

 I am having troubles increasing the fontize when plotting a 
 variogram{gstat} and its model (vgm) with plot and using jpeg(). Also
   

   
 the resolution in the jpeg call does not work. I am using R2.4 on 
 Mandriva 10.2 linux.

 I can change fontsize with cex.axis in a normal plot, so I presume it
   

   
 has to do with plotting the variogram model. Any help on how to 
 increase the font size and resolution would be appreciated?

 R-calls:
 v1-variogram(log(z)~x+y, loc=coordinates(shp1),data=shp1) 
 m1-vgm(0.0175, Gau, 20,0.052)

 jpeg(file=LOTPLAN_variogram_mod.jpg, bg=white, res=300, 
 pointsize=16, width=1200, height=1200, quality=100)

  plot(v1,plot.number=F, model=m1, ylim=c(0.04,0.08), col=black,
 cex.axis=1.5)

 dev.off()

 Thanks Herry

 Dr Alexander Herr
 Spatial and statistical analyst
 CSIRO, Sustainable Ecosystems
 Davies Laboratory,
 University Drive, Douglas, QLD 4814
 Private Mail Bag, Aitkenvale, QLD 4814
  
 Phone/www
 (07) 4753 8510; 4753 8650(fax)
 Home: http://herry.ausbats.org.au
 Webadmin ABS: http://ausbats.org.au
 Sustainable Ecosystems: http://www.cse.csiro.au

Re: [R] [R-sig-Geo] plot() and Jpeg() increase font size and resolution

[Edzer J. Pebesma]

[EMAIL PROTECTED] wrote:
 Thanks to Edzer and Roger,
 I can now plot with increased font sizes. However, jpeg still does not
 reproduce these, nor does it show up in high quality. What I would like
 to do is produce some highresolution jpegs.
   
Please specify the commands that you used, preferrably such that we can 
reproduce them, and specify which worked and which didn't the way you 
expected. Also please give the output of sessionInfo(), and tell us the 
platform you work on. Have you tried producing png, did that work?
--
Edzer
 Any help would be appreciated

 Thanx
 Herry


 R2.4 on Mandriva 10.2 linux.

 Dr Alexander Herr
 Spatial and statistical analyst
 CSIRO, Sustainable Ecosystems
 Davies Laboratory,
 University Drive, Douglas, QLD 4814 
 Private Mail Bag, Aitkenvale, QLD 4814
  
 Phone/www 
 (07) 4753 8510; 4753 8650(fax)
 Home: http://herry.ausbats.org.au
 Webadmin ABS: http://ausbats.org.au
 Sustainable Ecosystems: http://www.cse.csiro.au/
 


 -Original Message-
 From: Edzer J. Pebesma [mailto:[EMAIL PROTECTED] 
 Sent: Friday, November 24, 2006 5:37 PM
 To: Herr, Alexander Herr - Herry (CSE, Townsville)
 Cc: r-help@stat.math.ethz.ch; r-sig-geo@stat.math.ethz.ch
 Subject: Re: [R-sig-Geo] plot() and Jpeg() increase font size and
 resolution

 plotting variograms in gstat is done through xyplot in lattice; you'll
 find where it gets it's defaults by

 library(lattice)
 trellis.par.get()
 ?trellis.par.set
 --
 Edzer

 [EMAIL PROTECTED] wrote:
   
  
 Dear list,

 I am having troubles increasing the fontize when plotting a 
 variogram{gstat} and its model (vgm) with plot and using jpeg(). Also 
 the resolution in the jpeg call does not work. I am using R2.4 on 
 Mandriva 10.2 linux.

 I can change fontsize with cex.axis in a normal plot, so I presume it 
 has to do with plotting the variogram model. Any help on how to 
 increase the font size and resolution would be appreciated?

 R-calls:
 v1-variogram(log(z)~x+y, loc=coordinates(shp1),data=shp1) 
 m1-vgm(0.0175, Gau, 20,0.052)

 jpeg(file=LOTPLAN_variogram_mod.jpg, bg=white, res=300, 
 pointsize=16, width=1200, height=1200, quality=100)

  plot(v1,plot.number=F, model=m1, ylim=c(0.04,0.08), col=black,
 cex.axis=1.5)

 dev.off()

 Thanks Herry

 Dr Alexander Herr
 Spatial and statistical analyst
 CSIRO, Sustainable Ecosystems
 Davies Laboratory,
 University Drive, Douglas, QLD 4814
 Private Mail Bag, Aitkenvale, QLD 4814
  
 Phone/www
 (07) 4753 8510; 4753 8650(fax)
 Home: http://herry.ausbats.org.au
 Webadmin ABS: http://ausbats.org.au
 Sustainable Ecosystems: http://www.cse.csiro.au/

 ___
 R-sig-Geo mailing list
 R-sig-Geo@stat.math.ethz.ch
 https://stat.ethz.ch/mailman/listinfo/r-sig-geo

 

 ___
 R-sig-Geo mailing list
 R-sig-Geo@stat.math.ethz.ch
 https://stat.ethz.ch/mailman/listinfo/r-sig-geo


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] [R-sig-Geo] plot() and Jpeg() increase font size and resolution

[Edzer J. Pebesma]

plotting variograms in gstat is done through xyplot in lattice; you'll 
find where it gets it's defaults by

library(lattice)
trellis.par.get()
?trellis.par.set
--
Edzer

[EMAIL PROTECTED] wrote:
  
 Dear list,

 I am having troubles increasing the fontize when plotting a
 variogram{gstat} and its model (vgm) with plot and using jpeg(). Also
 the resolution in the jpeg call does not work. I am using R2.4 on
 Mandriva 10.2 linux.

 I can change fontsize with cex.axis in a normal plot, so I presume it
 has to do with plotting the variogram model. Any help on how to increase
 the font size and resolution would be appreciated?

 R-calls:
 v1-variogram(log(z)~x+y, loc=coordinates(shp1),data=shp1)
 m1-vgm(0.0175, Gau, 20,0.052)

 jpeg(file=LOTPLAN_variogram_mod.jpg, bg=white, res=300,
 pointsize=16, width=1200, height=1200, quality=100)

  plot(v1,plot.number=F, model=m1, ylim=c(0.04,0.08), col=black,
 cex.axis=1.5)

 dev.off()

 Thanks Herry

 Dr Alexander Herr
 Spatial and statistical analyst
 CSIRO, Sustainable Ecosystems
 Davies Laboratory,
 University Drive, Douglas, QLD 4814 
 Private Mail Bag, Aitkenvale, QLD 4814
  
 Phone/www 
 (07) 4753 8510; 4753 8650(fax)
 Home: http://herry.ausbats.org.au
 Webadmin ABS: http://ausbats.org.au
 Sustainable Ecosystems: http://www.cse.csiro.au/

 ___
 R-sig-Geo mailing list
 R-sig-Geo@stat.math.ethz.ch
 https://stat.ethz.ch/mailman/listinfo/r-sig-geo


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Interpolate univariate data on regular 3D grid to new 3D grid

[Edzer J. Pebesma]

Package gstat allows 3D interpolation. It doesn't accept
a 3D array directly; you'll have to provide it as

x1 y1 z1 obs1
x2 y2 z2 obs2

etc, e.g. by using expand.grid. Or, you may want to try
out the classes provided by package sp, which allow for
more than 2 dimensions in case of points and grids.
--
Edzer

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] [R-pkgs] new gstat version

[Edzer J. Pebesma]

Soon on CRAN a new version of package gstat will emerge, which
has a few minor changes and possible incompatibilities w.r.t. the
previous version(s).

The new gstat (0.9-23) now:
+ depends on sp, and uses internally with Spatial* classes from sp
   if data are provided in the old-fashoned way (as data.frame)
+ has a vignette to get you started with the classes in sp
+ defines krige as a generic; two typical uses are:

# old-style, using data.frame's:
library(gstat)
data(meuse)
data(meuse.grid)
v = vgm(0.6, Sph, 900, 0.1) # spherical variogram model
zn.kri1 = krige(log(zinc)~1, ~x+y, meuse, meuse.grid, v)
image(zn.kri1)

# new-style, using SpatialPointsDataFrame's:
coordinates(meuse) = c(x, y) # promote meuse to SpatialPointsDataFrame
coordinates(meuse.grid) = c(x, y) # promote to SpatialPointsDataFrame
gridded(meuse.grid) = TRUE # promote to SpatialPixelsDataFrame
zn.kri2 = krige(log(zinc)~1, meuse, meuse.grid, v)
spplot(zn.kri2[1])

+ provides a generic idw for inverse distance interpolation
(calling sequence identical to krige)

INCOMPATIBLE CHANGE:
  the krige generic is

 krige(formula, locations, ...)

  based on the class of argument 2 (locations) either of the two
  versions (see above) is called. This means that no named argument
  data should be passed in case locations is missing (i.e., when
  data is of class SpatialPointsDataFrame). Instead, pass it as
  argument 2 or as named argument locations.

--
Edzer

___
R-packages mailing list
[EMAIL PROTECTED]
https://stat.ethz.ch/mailman/listinfo/r-packages

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] How to use filled.contour(x,y,z) data for levelplot(z)?

[Edzer J. Pebesma]

Jan, 

libary(lattice)
rownames(zz) - xx
colnames(zz) - yy
levelplot(zz)

will not work -- levelplot needs a data frame with
zz, xx and yy values next to each other, so read
about data.frame(), rep() and as.vector(zz); you
will need the each= argument for one of xx or yy.

Sounds like you'd be helped by using package sp, 
see also http://r-spatial.sourceforge.net/ and look
for the graph gallery, function spplot.
--
Edzer

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] Variogram

[Edzer J. Pebesma]

Leaf, please note that r-help is not the appropriate place to ask
package-specific questions. We have r-sig-geo for questions related
to geographic data in R, and gstat has a mailing list on its own.

The answer is below.
--
Edzer


Leaf wrote:

Dear All,

Is there anybody has the experience in using variogram(gstat) ? Please kindly 
give me some hints about the results.


I used variogram() to build a semivariogram plot as:

tr.var=variogram(Incr~1,loc=~X+Y,data=TRI2TU,width=5)

then fir the variogram to get the parameters as:

 v.fit = fit.variogram(tr.var,vgm(0.5,Exp,300,1))

v.fit
  modelpsillrange
1   Nug 1.484879  0.0
2   Exp 3.476700 29.70914

This is the output of v.fit. Can anybody help me write the exponential formula 
for this variogram?  I have the problem in understanding the result.

BTW

The equation you're looking for is:

if h = 0, gamma(h) = 0
if h  0, gamma(h) = 1.484879 + 3.4767 (1 - exp(-h/29.70914))

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] Error in vector(double, length) : vector size specified is too large....VLDs

[Edzer J. Pebesma]

Tom,

please try to use the variogram function in package gstat;
it doesn't (try to) store all pairwise differences, but rather
accumulates them for distance intervals.

It will take a while to do this, and there is a chance that
you overflow the counter that keeps the number of point
pairs for each interval: 304000^2  2^32; it is stored as
a C long, so may work on a 64 bit architecture. Otherwise,
I'd suggest to sample your data set.

I'd be interested to hear whether you succeed (or not).
--
Edzer

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] adding text to the corner of a lattice plot

[Edzer J. Pebesma]

The following worked for me:

  xyplot(rnorm(100)~ rnorm(100), 
key=list(text=list(title=mytext),corner=c(0,0)))

but I'm sure there are better ways of doing this (esp. if you need
the key for the key!)
--
Edzer

Matt Wiener wrote:

I think this message from the help archives might address your question
(found using the query lattice mtext on R site search:
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/51605.html.

Hope this helps,

Matt Wiener

-Original Message-
From: r-help-bounces at stat.math.ethz.ch 
https://stat.ethz.ch/mailman/listinfo/r-help
[mailto:r-help-bounces at stat.math.ethz.ch 
https://stat.ethz.ch/mailman/listinfo/r-help] On Behalf Of RINNER Heinrich
Sent: Friday, September 09, 2005 4:16 AM
To: r-help at stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-help
Subject: [R] adding text to the corner of a lattice plot


Dear R community,

I am using R 2.1.1 on Windows XP, package lattice Version 0.12-5, and
want to add text (sort of a dat-stamp actually) to the lower left corner
of a lattice plot, prefarably _after_ the plot has been created.

Here is a simple example what I do in base graphics:

# base graphics:
 plot(rnorm(100), rnorm(100))
 mtext(as.character(Sys.Date()), side = 1,line = -2, outer = T, adj =
0, font = 1, cex = 0.7)

How can I get the same using lattice?
# lattice:
 require(lattice)
 xyplot(rnorm(100) ~ rnorm(100))
 ???

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] [R-pkgs] New CRAN package sp: classes and methods for spatial data

[Edzer J. Pebesma]

We're happy to announce the CRAN release of sp, an R package which
has new-style classes and methods for spatial data, version 0.7-9.

Spatial data types that sp implements are: points, grids, lines, and
polygons (i.e., rings) optionally with holes. Methods include

+ the usual print, summary, plot, [, [[, $, ...
+ coercion between types (e.g. points and grids, matrices, data.frames)
+ coordinates(x), which returns the spatial coordinates of x
+ bbox(x), returns the coordinates bounding box of x
+ overlay, to query the value of e.g. points in polygons or grid
   (essentially does a point-in-polygon or point-in-raster cell)
+ spsample, for random sampling methods over a spatial domain.

An additional package (spproj) provides coordinate reference system
transformation (projection and re-projection) using the PROJ.4 library
[2]. Others (will) provide interfaces to GRASS 6 and gdal.

A good deal of work has also gone into providing plotting methods using
base, grid and lattice graphics, through the spplot function, a front-end
to lattice plots for spatial data (see gallery [1]).

The home page of these packages is found at
http://r-spatial.sourceforge.net/
See also the Task View on Spatial Data Analysis, linked from CRAN.

The reason why we wrote this package is that we think R is an excellent
environment to deal with spatial data, but that it lacks a uniform way
to deal with spatial data. Compared to the handling of dates and times,
which can utilize base classes or those provided in the chron package,
spatial data handling is much more fragmented. As a consequence:

- various packages make their own assumptions about how spatial data
   are organized
- spatial data organized for a certain package cannot easily be used
   for another package
- few (or no) packages address the full range of spatial data types
   (points, grids, lines, polygons) and their interaction
- generic spatial functionality (e.g. I/O to GIS, plotting, projection)
   is scattered and often limited in functionality.

It also means that many different package authors have to use time writing
similar data handling code, rather than concentrating on analytical
functions. If the sp package achieves its goals, data I/O will become
many-to-one, and data access for analysis one-to-many, providing a shared
data object layer for which shared methods can be written.

Classes and methods for spatial data are only useful when the spatial
packages support them. The sp development team includes maintainers of
a number of spatial R packages, on which we will work, and we hope that
over time other spatial package maintainers will also provide support
for the classes provided by this package.

Although we are working towards a fixed set of classes at this moment
of the development of sp we cannot fully guarantee that the exact
representation of the sp classes will not change in the future. Therefore,
we advise users to keep the scripts with instructions of building the
data into sp classes; it is most likely that we will not change the
building functions and methods.

Work in progress currently involves:
- support of sp classes by several spatial statistics packages
- automatic determination of hole polygons from shapefiles
- plotting filled polygons with holes
- gdal read/write support
- GRASS 6.0 support

The development of this package is a joint effort of Virgilio Gomez-Rubio,
Barry Rowlinson, Roger Bivand and Edzer Pebesma, and followed from
discussions held at a pre-DSC2003 workshop [3], announcements on R-sig-geo
[4], and a meeting held last November in Lancaster [5]. A beta release
of sp was announced a while ago on R-sig-geo [6].

With best regards,
--
Roger Bivand and Edzer Pebesma


[1] http://r-spatial.sourceforge.net/
[2] http://www.remotesensing.org/proj/
[3] http://spatial.nhh.no/meetings/vienna/index.html
[4] e.g. https://stat.ethz.ch/pipermail/r-sig-geo/2003-October/28.html
[5] http://elearning.maths.lancs.ac.uk:8080/RSpatial/
[6] https://stat.ethz.ch/pipermail/r-sig-geo/2005-April/000378.html

___
R-packages mailing list
[EMAIL PROTECTED]
https://stat.ethz.ch/mailman/listinfo/r-packages

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] converting R objects to C types in .Call

[Edzer J. Pebesma]

Is there a specific reason why, instead of
   CHAR(STRING_ELT(chstr, 0));
The S-Plus compatible
  CHARACTER_POINTER(chstr)[0]
does not work in R?
--
Edzer
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] Re: Computation of space-time empirical crosscovariances

[Edzer J. Pebesma]

Giovanna,
space-time cross covariance has never seemed to me something
that can be done simple, but package gstat at least offers you:
1. to include time as a third dimension, and model 3D anisotropy,
2. to calculate cross variograms between different moments in
time, and proceed with cokriging.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] Cross-variograms

[Edzer J. Pebesma]

Jacques, provided that X and Y are colocated (i.e., have
exactly the same observation locations), you get the
cross variogram right; the definition of this cross
variogram is however:
gamma(h)= E[(X(s)-X(s+h))*(Y(s)-Y(s+h))]
also, where you select:
cv - v$gamma[1:14]
you may be better off using the more general
v$gamma[v$id == X.Y]
Best regards,
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] [R-pkgs] gstat 0.9-12: cokriging cross validation and class name incompatibilities

[Edzer J. Pebesma]

I uploaded gstat 0.9-12 to CRAN, which has a few important changes:
1. Cokriging cross validation
Cokriging cross validation is now possible with the function gstat.cv:
you simply pass a multivariable gstat object, and cross validation is done
for the first variable in the object. Optionally, secondary variable 
records
at locations coinciding with the validation locations are removed.

2. Class name changes
Both gstat and geoR used the name variogram for -- highly incompatible --
objects that contain information about a variogram. This led to errors when
calling plot() with a variogram object calculated from one, but plotted by
the other package. I changed the class name into gstatVariogram. The full
list of class name changes:
variogram - gstatVariogram
point.pairs - pointPairs
variogram.map - variogramMap
variogram.cloud - variogramCloud
You can still run scripts you have with the new packages, but e.g. plot(v)
when v is an old object of class variogram will not work with the new
package. As these are S3 classes, a simple solution would be to
reassign the class:
class(v) = c(gstatVariogram, data.frame)
3. Compatibility with package sp and variogram maps
On http://sourceforge.net/projects/r-spatial/ beta releases for the upcoming
package sp are available; sp provides classes and methods for spatial data.
Instead of:
 library(gstat)
 data(meuse)
 variogram(zinc~1, ~x+y, meuse)
you could, with package sp loaded, do:
 coordinates(meuse) = ~x+y # promotes meuse to SpatialDataFrame
 variogram(zinc~1, meuse) # no coordinates required.
In addition, sp provides classes for gridded data and polygon data. With sp
loaded, gstat can calculate variogram maps.
Any comments are welcome,
--
Edzer
___
R-packages mailing list
[EMAIL PROTECTED]
https://www.stat.math.ethz.ch/mailman/listinfo/r-packages
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] Simulating a landscape (matrix) in R

[Edzer J. Pebesma]

For creating larger landscapes (Gaussian random fields), you may
consider using package gstat, which uses the sequential Gaussian simulation
algorithm with local approximations. In the example below, only the nearest
20 neighbours are used to approximate each of the conditional distributions.
The example below which simulates  plots a 300x300 grid, took 1 minute
on a 2 GHz CPU.
Best regards,
--
Edzer
library(gstat)
g = gstat(NULL, id0, z~1, ~x+y, dummy=TRUE,nmax=20,
   model=vgm(1,Exp,500),beta=10)
locs=expand.grid(x=1:300,y=1:300)
out - predict.gstat(g, locs, nsim=1)
levelplot(sim1~x+y,out,col.regions=bpy.colors())
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] graphical interface -- gstat/lattice

[Edzer J. Pebesma]

Carlos, what you see are the lattice default settings;
use trellis.par.get/set to get and change them:
library(gstat) # requires lattice
data(meuse)
trellis.par.set(background, list(col=#ff))
trellis.par.set(plot.symbol, list(cex=1,col=#00,font=1,pch=+))
plot(variogram(zinc~1,~x+y,meuse))
for more help, try:
?Lattice
?trellis.par.set
Please use a more descriptive subject, next time.

Best regards,
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] Interpolation or Kriging

[Edzer J. Pebesma]

Sebastien,

Please have a look at the geostatistics packages refered to from

http://sal.agecon.uiuc.edu/csiss/Rgeo/

BTW -- is this page referenced from the R or CRAN home
pages? I couldn't find a link a few days ago.
Most of them may do what you want. I assume they can all
be fooled by providing one-dimensional data having a constant
second coordinate. Package gstat allows explicit specification
of one-dimensional data sets, but has no special provisions
for 1D data in terms of efficiency, compared to 2D data.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] Memory limitation in GeoR - Windows or R?

[Edzer J. Pebesma]

Usually, with that many observations, you use
kriging in a local neighbourhood, i.e. use only the
n nearest observations in a kriging system. If n is
pretty large, this is practically equivalent to kriging
with a global neighbourhood.
There are several other R packages that provide kriging,
some of which provide local kriging, e.g. gstat does
not require calculation of the full distance matrix.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] predict.lm with (logical) NA vector

[Edzer J. Pebesma]

I was surprised by the following (R 1.8.0):

R lm.fit = lm(y~x, data.frame(x=1:10, y=1:10))
R predict(lm.fit, data.frame(x = rep(NA, 10)))
1  2  3  4  5
-1.060998e-314 -1.060998e-314 -1.060998e-314 -1.060998e-314 -1.060998e-314
6  7  8  9 10
 0.00e+00  1.406440e-269  6.715118e-265  4.940656e-323  1.782528e-265
R predict(lm.fit, data.frame(x = as.numeric(rep(NA, 10
1  2  3  4  5  6  7  8  9 10
NA NA NA NA NA NA NA NA NA NA
shouldn't the first predict() call return NA's, or else issue an error 
message?
--
Edzer

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] Lattice: no grid name space

[Edzer J. Pebesma]

The following now occurs to me when I try to
load lattice (R 1.7.1, debian stable):
library(lattice)
Error in loadNamespace(i, c(lib.loc, .libPaths()), keep.source) :
   package `grid' does not have a name space
Error in library(lattice) : package/namespace load failed
I did an update.packages() as root recently.
Any idea what's wrong?
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] Indicator Kriging

[Edzer J. Pebesma]

Indicator kriging is nothing but ordinary or simple
kriging on binary (0/1) data; there's quite a list of packages
that can do kriging: gstat, geoR, sgeostat, fields,
spatial, ...
Be aware that indicator kriging outcomes only estimate
probabilities and may lie outside [0,1]
To my knowledge there's no package that ports GSLIB to R.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] More questions about R extension programming

[Edzer J. Pebesma]

I have been strugling to get character data in and out of c functions
using the .Call() interface, but I don't succeed getting them back
into R.
The example:

#include R.h
#include Rdefines.h
char *cp = xxx;
SEXP do_char(SEXP c) {
   int i;
   SEXP s;
   for (i = 0; i  LENGTH(c); i++)
   printf(%d:[%s]\n, i, CHARACTER_DATA(STRING_ELT(c,i)));
   /* CHARACTER_VALUE() does not work here! */
   s = NEW_CHARACTER(1);
   CHARACTER_POINTER(s) = cp;
   c = NEW_CHARACTER(3);
   SET_STRING_ELT(c,0, s);
   SET_STRING_ELT(c,1, s);
   SET_STRING_ELT(c,2, s);
   return c;
}
works only for the first half; but I can't manage to fill the string
and return it back.
Note that I had to use CHARACTER_DATA instead of CHARACTER,
which Thomas Lumley suggested last Saturday.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] auto.key = TRUE in xyplot() draws max. 7 groups

[Edzer J. Pebesma]

In the example below:

library(lattice)
n = 100
a = rnorm(n)
b = rnorm(n)
c = sample(c(1:7), n, repl=TRUE)
xyplot(a ~ b, groups = c, auto.key = TRUE)
c = sample(c(1:8), n, repl=TRUE)
xyplot(a ~ b, groups = c, auto.key = TRUE)
a key is drawn for the first, but not for the
second graph. I suppose it has to do with the
number of groups. If not a bug, maybe
this should go into the documentation?
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] Z aware interpolation - clarification

[Edzer J. Pebesma]

Monica,

gstat has interpolation facilities for real 3D data: IDW, as
well as variogram modelling  kriging prediction and simulation.
Try:

library(gstat)
demo(gstat3D)
gstat exists as an R package, as well as a stand alone binary.
The package is at CRAN, the rest at www.gstat.org
best regards,
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] levelplot behaviour when at cuts the z range

[Edzer J. Pebesma]

Consider the following examples:

library(lattice)
x = c(1,1,2,2)
y = c(1,2,1,2)
z = 1:4
levelplot(z~x+y,at=c(.5, 1.5, 2.5, 3.5, 4.5)) # correct
levelplot(z~x+y,at=c(.5, 1.5, 2.5))   # ?
The second plot is clearly incorrect. However, I don't know
what correct behaviour is: ignore everything above 2.5 and
issue a warning? Reject at values that do not cover the
data range, and issue an error message? I think I prefer
the second, as it does not hide extremes.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] levelplot behaviour for panel with constants

[Edzer J. Pebesma]

In the example:

x = rep(c(0,0,1,1),4)
y = rep(c(0,1,0,1),4)
z = c(1,0,1,0,0,0,1,1,0,1,0,0,1,1,1,1)
f = as.factor(c(rep(a,4),rep(b,4),rep(c,4),rep(d,4)))
levelplot(z~x+y|f,data.frame(x=x,y=y,z=z,f=f))
I noted that the last (d) plot remains empty. I guess the
reason for this is that the values are constant (1), but I consider
it more consistent if they would get the colour of 1, and would
be left blank in case they were NA's.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

Re: [R] krige in gstat() package

[Edzer J. Pebesma]

HI,
  I wonder does anyone have experience with doing sequential gaussian
simulation with krige() function in gstat?
I find it VERY slow compared to use krige() to achieve kriging function
itself..  I wonder why, is that because it has to model the variogram, and
do the kriging separately for each point to be simulated?
It does not model variograms on the way. It is slower than kriging
because it uses the sequential simulation algorithm: for each node
visited, a value is simulated, and this value is added to the
conditioning data. For this reason you _have to_ limit the search
neighbourhood (using the nmax or maxdist arguments) if you're
simulating on more than say a few hundred nodes. Taking a very
small nmax may yield fast simulations, at the expense of a good
representation of the target spatial correlation structure.
so it would be N times slower to achieve the simulation than the kriging
if the number of points to be estimated is N??
No. Experiment with the nmax argument.
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

[R] Gstat: multivariable geostatistics for S (R and S-Plus)

[Edzer J. Pebesma]

The majority of the functionality present in the gstat stand-alone
program (http://www.gstat.org/) is now available as a package/library for
the S language (R, S-Plus), again called gstat. The package provides
multivariable geostatistical modelling, prediction and simulation, as
well as several visualisation functions.  Gstat was started 10 years
ago and was released under the GPL in 1996; the original stand-alone
program is closely linked to several GIS systems. Gstat was not initially
written for teaching purposes, but for research purposes, emphasizing
flexibility, scalability and portability. It can deal with a large number
of practical issues in geostatistics, including change of support (block
kriging), simple/ordinary/universal (co)kriging, fast local neighbourhood
selection, flexible trend modelling, variables with different sampling
configurations, and efficient simulation of large spatially correlated
random fields, indicator kriging and simulation, and (directional)
variogram and cross variogram modelling. The S formula/models interface
is used to define multivariable geostatistical models.

The source and windows package for R are available from CRAN. The page on
http://www.gstat.org/s.html has links to R and S-Plus (6.x) source code,
as well as examples, graphs, and a longer list of features. (A binary
Win32 S-Plus library is planned later this year.)

On http://www.ci.tuwien.ac.at/Conferences/DSC-2003/ you can find a
draft paper further describing the package.

Any feedback is appreciated.
--
Edzer

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help