R : Copyright 2005, The R Foundation for Statistical Computing
Version 2.1.1  (2005-06-20), ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for a HTML browser interface to help.
Type 'q()' to quit R.

> ### * <HEADER>
> ###
> attach(NULL, name = "CheckExEnv")
> assign(".CheckExEnv", as.environment(2), pos = length(search())) # base
> ## add some hooks to label plot pages for base and grid graphics
> setHook("plot.new", ".newplot.hook")
> setHook("persp", ".newplot.hook")
> setHook("grid.newpage", ".gridplot.hook")
> 
> assign("cleanEx",
+        function(env = .GlobalEnv) {
+ 	   rm(list = ls(envir = env, all.names = TRUE), envir = env)
+            RNGkind("default", "default")
+ 	   set.seed(1)
+    	   options(warn = 1)
+ 	   delayedAssign("T", stop("T used instead of TRUE"),
+ 		  assign.env = .CheckExEnv)
+ 	   delayedAssign("F", stop("F used instead of FALSE"),
+ 		  assign.env = .CheckExEnv)
+ 	   sch <- search()
+ 	   newitems <- sch[! sch %in% .oldSearch]
+ 	   for(item in rev(newitems))
+                eval(substitute(detach(item), list(item=item)))
+ 	   missitems <- .oldSearch[! .oldSearch %in% sch]
+ 	   if(length(missitems))
+ 	       warning("items ", paste(missitems, collapse=", "),
+ 		       " have been removed from the search path")
+        },
+        env = .CheckExEnv)
> assign("..nameEx", "__{must remake R-ex/*.R}__", env = .CheckExEnv) # for now
> assign("ptime", proc.time(), env = .CheckExEnv)
> grDevices::postscript("gstat-Examples.ps")
> assign("par.postscript", graphics::par(no.readonly = TRUE), env = .CheckExEnv)
> options(contrasts = c(unordered = "contr.treatment", ordered = "contr.poly"))
> options(warn = 1)    
> library('gstat')
> 
> assign(".oldSearch", search(), env = .CheckExEnv)
> assign(".oldNS", loadedNamespaces(), env = .CheckExEnv)
> cleanEx(); ..nameEx <- "bpy.colors"
> 
> ### * bpy.colors
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: bpy.colors
> ### Title: blue-pink-yellow color scheme that prints well as grey tone
> ### Aliases: bpy.colors
> ### Keywords: color
> 
> ### ** Examples
> 
> bpy.colors(10)
 [1] "#000033" "#000099" "#0000FF" "#5000FF" "#9F0FF0" "#EF42BD" "#FF758A"
 [8] "#FFA857" "#FFDB24" "#FFFF60"
> p <- expand.grid(x=1:30,y=1:30)
> p$z <- p$x + p$y
> image(p, col = bpy.colors(100))
> # require(lattice)
> # trellis.par.set("regions", list(col=bpy.colors())) # make default
> 
> 
> 
> cleanEx(); ..nameEx <- "bubble"
> 
> ### * bubble
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: bubble
> ### Title: Create a bubble plot of spatial data
> ### Aliases: bubble
> ### Keywords: dplot
> 
> ### ** Examples
> 
> data(meuse)
> bubble(meuse, max = 2.5, main = "cadmium concentrations (ppm)", 
+     key.entries = c(.5,1,2,4,8,16))
> bubble(meuse, "x", "y", "zinc", main = "zinc concentrations (ppm)",
+     key.entries =  100 * 2^(0:4))
> 
> 
> 
> cleanEx(); ..nameEx <- "fit.lmc"
> 
> ### * fit.lmc
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: fit.lmc
> ### Title: Fit a Linear Model of Coregionalization to a Multivariable
> ###   Sample Variogram
> ### Aliases: fit.lmc
> ### Keywords: models
> 
> ### ** Examples
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "fit.variogram"
> 
> ### * fit.variogram
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: fit.variogram
> ### Title: Fit a Variogram Model to a Sample Variogram
> ### Aliases: fit.variogram
> ### Keywords: models
> 
> ### ** Examples
> 
> data(meuse)
> vgm1 <- variogram(log(zinc)~1, ~x+y, meuse)
> fit.variogram(vgm1, vgm(1,"Sph",300,1))
  model     psill    range
1   Nug 0.0506555   0.0000
2   Sph 0.5906009 896.9702
> 
> 
> 
> cleanEx(); ..nameEx <- "fit.variogram.reml"
> 
> ### * fit.variogram.reml
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: fit.variogram.reml
> ### Title: REML Fit Direct Variogram Partial Sills to Data
> ### Aliases: fit.variogram.reml
> ### Keywords: models
> 
> ### ** Examples
> 
> data(meuse)
> fit.variogram.reml(log(zinc)~1, ~x+y, meuse, model = vgm(1, "Sph", 900,1))
  model      psill range
1   Nug 0.02547708     0
2   Sph 0.61088826   900
> 
> 
> 
> cleanEx(); ..nameEx <- "fulmar"
> 
> ### * fulmar
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: fulmar
> ### Title: Fulmaris glacialis data
> ### Aliases: fulmar
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(fulmar)
> summary(fulmar)
      year            x                y               depth     
 Min.   :1998   Min.   :476210   Min.   :5694947   Min.   : 1.0  
 1st Qu.:1998   1st Qu.:535522   1st Qu.:5806777   1st Qu.:13.0  
 Median :1999   Median :568618   Median :5896021   Median :25.0  
 Mean   :1999   Mean   :576982   Mean   :5889112   Mean   :23.6  
 3rd Qu.:1999   3rd Qu.:604579   3rd Qu.:5946550   3rd Qu.:32.0  
 Max.   :1999   Max.   :739042   Max.   :6150942   Max.   :54.0  
     coast             fulmar      
 Min.   :  1.024   Min.   : 0.000  
 1st Qu.:  4.857   1st Qu.: 0.000  
 Median : 38.696   Median : 0.000  
 Mean   : 56.642   Mean   : 1.005  
 3rd Qu.: 89.929   3rd Qu.: 0.000  
 Max.   :263.205   Max.   :46.487  
> 
> 
> 
> cleanEx(); ..nameEx <- "get.contr"
> 
> ### * get.contr
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: get.contr
> ### Title: Calculate contrasts from multivariable predictions
> ### Aliases: get.contr
> ### Keywords: models
> 
> ### ** Examples
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "gstat"
> 
> ### * gstat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: gstat
> ### Title: Create gstat objects, or subset it
> ### Aliases: gstat print.gstat [.gstat
> ### Keywords: models
> 
> ### ** Examples
> 
> data(meuse)
> # let's do some manual fitting of two direct variograms and a cross variogram
> g <- gstat(id = "ln.zinc", formula = log(zinc)~1, locations = ~x+y, 
+         data = meuse)
> g <- gstat(g, id = "ln.lead", formula = log(lead)~1, locations = ~x+y, 
+         data = meuse)
> # examine variograms and cross variogram:
> plot(variogram(g))
> # enter direct variograms:
> g <- gstat(g, id = "ln.zinc", model = vgm(.55, "Sph", 900, .05))
> g <- gstat(g, id = "ln.lead", model = vgm(.55, "Sph", 900, .05))
> # enter cross variogram:
> g <- gstat(g, id = c("ln.zinc", "ln.lead"), model = vgm(.47, "Sph", 900, .03))
> # examine fit:
> plot(variogram(g), model = g$model, main = "models fitted by eye")
> # see also demo(cokriging) for a more efficient approach
> g["ln.zinc"]
data:
ln.zinc : formula = log(zinc)`~`1 ; locations = ~x + y ; data dim = 155 x 14
variograms:
           model psill range
ln.zinc[1]   Nug  0.05     0
ln.zinc[2]   Sph  0.55   900
> g["ln.lead"]
data:
ln.lead : formula = log(lead)`~`1 ; locations = ~x + y ; data dim = 155 x 14
variograms:
           model psill range
ln.lead[1]   Nug  0.05     0
ln.lead[2]   Sph  0.55   900
> g[c("ln.zinc", "ln.lead")]
data:
ln.zinc : formula = log(zinc)`~`1 ; locations = ~x + y ; data dim = 155 x 14
ln.lead : formula = log(lead)`~`1 ; locations = ~x + y ; data dim = 155 x 14
variograms:
                   model psill range
ln.zinc[1]           Nug  0.05     0
ln.zinc[2]           Sph  0.55   900
ln.lead[1]           Nug  0.05     0
ln.lead[2]           Sph  0.55   900
ln.zinc.ln.lead[1]   Nug  0.03     0
ln.zinc.ln.lead[2]   Sph  0.47   900
> g[1]
data:
ln.zinc : formula = log(zinc)`~`1 ; locations = ~x + y ; data dim = 155 x 14
variograms:
           model psill range
ln.zinc[1]   Nug  0.05     0
ln.zinc[2]   Sph  0.55   900
> g[2]
data:
ln.lead : formula = log(lead)`~`1 ; locations = ~x + y ; data dim = 155 x 14
variograms:
           model psill range
ln.lead[1]   Nug  0.05     0
ln.lead[2]   Sph  0.55   900
> 
> # Inverse distance interpolation with inverse distance power set to .5:
> # (kriging variants need a variogram model to be specified)
> data(meuse)
> data(meuse.grid)
> meuse.gstat <- gstat(id = "zinc", formula = zinc ~ 1, locations = ~ x + y,
+         data = meuse, nmax = 7, set = list(idp = .5))
> meuse.gstat
data:
zinc : formula = zinc`~`1 ; locations = ~x + y ; data dim = 155 x 14 nmax = 7
set idp = 0.5; 
> z <- predict(meuse.gstat, meuse.grid)
[inverse distance weighted interpolation]
> library(lattice) # for levelplot
> levelplot(zinc.pred~x+y, z, aspect = mapasp(z))
> # see demo(cokriging) and demo(examples) for further examples, 
> # and the manuals for predict.gstat and image
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "image"
> 
> ### * image
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: image
> ### Title: Image Gridded Coordinates in Data Frame
> ### Aliases: image.data.frame image xyz2img
> ### Keywords: dplot
> 
> ### ** Examples
> 
> data(meuse)
> data(meuse.grid)
> g <- gstat(formula=log(zinc)~1,locations=~x+y,data=meuse,model=vgm(1,"Exp",300))
> x <- predict(g, meuse.grid)
[using ordinary kriging]
> image(x, 4, main="kriging variance and data points")
> points(meuse$x, meuse$y, pch = "+")
> # non-square cell test:
> image(x[((x$y - 20) %% 80) == 0,], main = "40 x 80 cells")
> image(x[((x$x - 20) %% 80) == 0,], main = "80 x 40 cells")
> # the following works for square cells only:
> oldpin <- par("pin")
> ratio <- length(unique(x$x))/length(unique(x$y))
> par(pin = c(oldpin[2]*ratio,oldpin[2]))
> image(x, main="Exactly square cells, using par(pin)")
> par(pin = oldpin)
> library(lattice)
> levelplot(var1.var~x+y, x, aspect = mapasp(x), main = "kriging variance")
> 
> 
> 
> graphics::par(get("par.postscript", env = .CheckExEnv))
> cleanEx(); ..nameEx <- "jura"
> 
> ### * jura
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: jura
> ### Title: Jura data set
> ### Aliases: jura prediction.dat validation.dat transect.dat
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(jura)
> summary(prediction.dat)
      Xloc            Yloc          Landuse           Rock      
 Min.   :0.626   Min.   :0.580   Min.   :1.000   Min.   :1.000  
 1st Qu.:2.282   1st Qu.:1.487   1st Qu.:2.000   1st Qu.:2.000  
 Median :3.043   Median :2.581   Median :3.000   Median :2.000  
 Mean   :2.980   Mean   :2.665   Mean   :2.548   Mean   :2.699  
 3rd Qu.:3.665   3rd Qu.:3.752   3rd Qu.:3.000   3rd Qu.:3.000  
 Max.   :4.920   Max.   :5.690   Max.   :4.000   Max.   :5.000  
       Cd               Co               Cr              Cu        
 Min.   :0.1350   Min.   : 1.552   Min.   : 8.72   Min.   :  3.96  
 1st Qu.:0.6375   1st Qu.: 6.520   1st Qu.:27.44   1st Qu.: 11.02  
 Median :1.0700   Median : 9.760   Median :34.84   Median : 17.60  
 Mean   :1.3091   Mean   : 9.303   Mean   :35.07   Mean   : 23.73  
 3rd Qu.:1.7150   3rd Qu.:11.980   3rd Qu.:42.22   3rd Qu.: 27.82  
 Max.   :5.1290   Max.   :17.720   Max.   :67.60   Max.   :166.40  
       Ni              Pb               Zn        
 Min.   : 4.20   Min.   : 18.96   Min.   : 25.20  
 1st Qu.:13.80   1st Qu.: 36.52   1st Qu.: 55.00  
 Median :20.56   Median : 46.40   Median : 73.56  
 Mean   :19.73   Mean   : 53.92   Mean   : 75.08  
 3rd Qu.:25.42   3rd Qu.: 60.40   3rd Qu.: 89.92  
 Max.   :53.20   Max.   :229.56   Max.   :219.32  
> summary(validation.dat)
      Xloc            Yloc          Landuse          Rock            Cd        
 Min.   :0.491   Min.   :0.524   Min.   :1.00   Min.   :1.00   Min.   :0.3250  
 1st Qu.:2.207   1st Qu.:1.593   1st Qu.:2.00   1st Qu.:2.00   1st Qu.:0.6765  
 Median :3.001   Median :2.389   Median :3.00   Median :2.00   Median :1.1865  
 Mean   :2.921   Mean   :2.546   Mean   :2.41   Mean   :2.36   Mean   :1.2343  
 3rd Qu.:3.716   3rd Qu.:3.339   3rd Qu.:3.00   3rd Qu.:3.00   3rd Qu.:1.6350  
 Max.   :4.745   Max.   :5.285   Max.   :4.00   Max.   :5.00   Max.   :3.7800  
       Co               Cr              Cu                Ni       
 Min.   : 1.652   Min.   : 3.32   Min.   :  3.552   Min.   : 1.98  
 1st Qu.: 7.950   1st Qu.:28.44   1st Qu.:  9.150   1st Qu.:15.28  
 Median :10.060   Median :34.54   Median : 16.140   Median :21.28  
 Mean   : 9.793   Mean   :34.88   Mean   : 23.218   Mean   :20.76  
 3rd Qu.:12.490   3rd Qu.:40.59   3rd Qu.: 23.190   3rd Qu.:25.36  
 Max.   :20.600   Max.   :70.00   Max.   :154.600   Max.   :43.68  
       Pb               Zn        
 Min.   : 18.68   Min.   : 25.00  
 1st Qu.: 35.31   1st Qu.: 53.19  
 Median : 47.00   Median : 73.92  
 Mean   : 56.48   Mean   : 77.96  
 3rd Qu.: 60.10   3rd Qu.: 90.40  
 Max.   :300.00   Max.   :259.84  
> summary(transect.dat)
       X           Rock.type        Block.Ni           Cd        
 Min.   :1.000   Min.   :1.000   Min.   : 6.61   Min.   : 0.135  
 1st Qu.:2.312   1st Qu.:1.000   1st Qu.:17.94   1st Qu.: 0.655  
 Median :3.625   Median :2.000   Median :20.24   Median : 1.317  
 Mean   :3.625   Mean   :2.047   Mean   :19.99   Mean   : 1.486  
 3rd Qu.:4.938   3rd Qu.:2.000   3rd Qu.:23.57   3rd Qu.: 1.961  
 Max.   :6.250   Max.   :4.000   Max.   :37.05   Max.   : 3.925  
                                                 NA's   :96.000  
       Ni       
 Min.   : 4.20  
 1st Qu.:13.31  
 Median :20.52  
 Mean   :19.62  
 3rd Qu.:24.98  
 Max.   :43.68  
 NA's   :90.00  
> 
> 
> 
> cleanEx(); ..nameEx <- "krige"
> 
> ### * krige
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: krige
> ### Title: Simple, Ordinary or Universal, global or local, Point or Block
> ###   Kriging, or simulation.
> ### Aliases: krige
> ### Keywords: models
> 
> ### ** Examples
> 
> data(meuse)
> data(meuse.grid)
> m <- vgm(.59, "Sph", 874, .04)
> # ordinary kriging:
> x <- krige(log(zinc)~1, ~x+y, model = m, data = meuse, newd = meuse.grid)
[using ordinary kriging]
> library(lattice)
> levelplot(var1.pred~x+y, x, aspect = mapasp(x),
+         main = "ordinary kriging predictions")
> levelplot(var1.var~x+y, x, aspect = mapasp(x),
+         main = "ordinary kriging variance")
> # simple kriging:
> x <- krige(log(zinc)~1, ~x+y, model = m, data = meuse, newdata = meuse.grid, 
+         beta=5.9)
[using simple kriging]
> # residual variogram:
> m <- vgm(.4, "Sph", 954, .06)
> # universal block kriging:
> x <- krige(log(zinc)~x+y, ~x+y, model = m, data = meuse, newdata = 
+         meuse.grid, block = c(40,40))
[using universal kriging]
> levelplot(var1.pred~x+y, x, aspect = mapasp(x),
+         main = "universal kriging predictions")
> # add grid:
> levelplot(var1.var~x+y, x, aspect = mapasp(x),
+                 panel = function(...) {
+                         panel.levelplot(...)
+                         panel.abline(h = 0:3*1000 + 330000, v= 0:2*1000 + 179000, col = "grey")
+                 },
+         main = "universal kriging variance")
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "krige.cv"
> 
> ### * krige.cv
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: krige.cv
> ### Title: (co)kriging cross validation, n-fold or leave-one-out
> ### Aliases: krige.cv gstat.cv
> ### Keywords: models
> 
> ### ** Examples
> 
> data(meuse)
> m <- vgm(.59, "Sph", 874, .04)
> # five-fold cross validation:
> x <- krige.cv(log(zinc)~1, ~x+y, model = m, data = meuse, nmax = 40, nfold=5)
[using ordinary kriging]
[using ordinary kriging]
[using ordinary kriging]
[using ordinary kriging]
[using ordinary kriging]
> bubble(x, z = "residual", main = "log(zinc): 5-fold CV residuals")
> 
> 
> 
> cleanEx(); ..nameEx <- "makegrid"
> 
> ### * makegrid
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: makegrid
> ### Title: make regular grid with square cells and round numbers
> ### Aliases: makegrid
> ### Keywords: dplot
> 
> ### ** Examples
> 
> data(meuse)
> grd <- makegrid(meuse$x, meuse$y, 1000)
> diff(grd$x)
  [1]   110   110   110   110   110   110   110   110   110   110   110   110
 [13]   110   110   110   110   110   110   110   110   110   110   110   110
 [25]   110   110 -2860   110   110   110   110   110   110   110   110   110
 [37]   110   110   110   110   110   110   110   110   110   110   110   110
 [49]   110   110   110   110   110 -2860   110   110   110   110   110   110
 [61]   110   110   110   110   110   110   110   110   110   110   110   110
 [73]   110   110   110   110   110   110   110   110 -2860   110   110   110
 [85]   110   110   110   110   110   110   110   110   110   110   110   110
 [97]   110   110   110   110   110   110   110   110   110   110   110 -2860
[109]   110   110   110   110   110   110   110   110   110   110   110   110
[121]   110   110   110   110   110   110   110   110   110   110   110   110
[133]   110   110 -2860   110   110   110   110   110   110   110   110   110
[145]   110   110   110   110   110   110   110   110   110   110   110   110
[157]   110   110   110   110   110 -2860   110   110   110   110   110   110
[169]   110   110   110   110   110   110   110   110   110   110   110   110
[181]   110   110   110   110   110   110   110   110
> diff(grd$y)
  [1]   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
 [19]   0   0   0   0   0   0   0   0 110   0   0   0   0   0   0   0   0   0
 [37]   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0 110
 [55]   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
 [73]   0   0   0   0   0   0   0   0 110   0   0   0   0   0   0   0   0   0
 [91]   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0 110
[109]   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
[127]   0   0   0   0   0   0   0   0 110   0   0   0   0   0   0   0   0   0
[145]   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0 110
[163]   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
[181]   0   0   0   0   0   0   0   0
> summary(grd)
       x                y         
 Min.   :178570   Min.   :329670  
 1st Qu.:179230   1st Qu.:329780  
 Median :180000   Median :330000  
 Mean   :180000   Mean   :330000  
 3rd Qu.:180770   3rd Qu.:330220  
 Max.   :181430   Max.   :330330  
> grd <- makegrid(meuse$x, meuse$y, cell.size = 40)
> diff(grd$x)
   [1]    40    40    40    40    40    40    40    40    40    40    40    40
  [13]    40    40    40    40    40    40    40    40    40    40    40    40
  [25]    40    40    40    40    40    40    40    40    40    40    40    40
  [37]    40    40    40    40    40    40    40    40    40    40    40    40
  [49]    40    40    40    40    40    40    40    40    40    40    40    40
  [61]    40    40    40    40    40    40    40    40    40    40 -2800    40
  [73]    40    40    40    40    40    40    40    40    40    40    40    40
  [85]    40    40    40    40    40    40    40    40    40    40    40    40
  [97]    40    40    40    40    40    40    40    40    40    40    40    40
 [109]    40    40    40    40    40    40    40    40    40    40    40    40
 [121]    40    40    40    40    40    40    40    40    40    40    40    40
 [133]    40    40    40    40    40    40    40    40    40 -2800    40    40
 [145]    40    40    40    40    40    40    40    40    40    40    40    40
 [157]    40    40    40    40    40    40    40    40    40    40    40    40
 [169]    40    40    40    40    40    40    40    40    40    40    40    40
 [181]    40    40    40    40    40    40    40    40    40    40    40    40
 [193]    40    40    40    40    40    40    40    40    40    40    40    40
 [205]    40    40    40    40    40    40    40    40 -2800    40    40    40
 [217]    40    40    40    40    40    40    40    40    40    40    40    40
 [229]    40    40    40    40    40    40    40    40    40    40    40    40
 [241]    40    40    40    40    40    40    40    40    40    40    40    40
 [253]    40    40    40    40    40    40    40    40    40    40    40    40
 [265]    40    40    40    40    40    40    40    40    40    40    40    40
 [277]    40    40    40    40    40    40    40 -2800    40    40    40    40
 [289]    40    40    40    40    40    40    40    40    40    40    40    40
 [301]    40    40    40    40    40    40    40    40    40    40    40    40
 [313]    40    40    40    40    40    40    40    40    40    40    40    40
 [325]    40    40    40    40    40    40    40    40    40    40    40    40
 [337]    40    40    40    40    40    40    40    40    40    40    40    40
 [349]    40    40    40    40    40    40 -2800    40    40    40    40    40
 [361]    40    40    40    40    40    40    40    40    40    40    40    40
 [373]    40    40    40    40    40    40    40    40    40    40    40    40
 [385]    40    40    40    40    40    40    40    40    40    40    40    40
 [397]    40    40    40    40    40    40    40    40    40    40    40    40
 [409]    40    40    40    40    40    40    40    40    40    40    40    40
 [421]    40    40    40    40    40 -2800    40    40    40    40    40    40
 [433]    40    40    40    40    40    40    40    40    40    40    40    40
 [445]    40    40    40    40    40    40    40    40    40    40    40    40
 [457]    40    40    40    40    40    40    40    40    40    40    40    40
 [469]    40    40    40    40    40    40    40    40    40    40    40    40
 [481]    40    40    40    40    40    40    40    40    40    40    40    40
 [493]    40    40    40    40 -2800    40    40    40    40    40    40    40
 [505]    40    40    40    40    40    40    40    40    40    40    40    40
 [517]    40    40    40    40    40    40    40    40    40    40    40    40
 [529]    40    40    40    40    40    40    40    40    40    40    40    40
 [541]    40    40    40    40    40    40    40    40    40    40    40    40
 [553]    40    40    40    40    40    40    40    40    40    40    40    40
 [565]    40    40    40 -2800    40    40    40    40    40    40    40    40
 [577]    40    40    40    40    40    40    40    40    40    40    40    40
 [589]    40    40    40    40    40    40    40    40    40    40    40    40
 [601]    40    40    40    40    40    40    40    40    40    40    40    40
 [613]    40    40    40    40    40    40    40    40    40    40    40    40
 [625]    40    40    40    40    40    40    40    40    40    40    40    40
 [637]    40    40 -2800    40    40    40    40    40    40    40    40    40
 [649]    40    40    40    40    40    40    40    40    40    40    40    40
 [661]    40    40    40    40    40    40    40    40    40    40    40    40
 [673]    40    40    40    40    40    40    40    40    40    40    40    40
 [685]    40    40    40    40    40    40    40    40    40    40    40    40
 [697]    40    40    40    40    40    40    40    40    40    40    40    40
 [709]    40 -2800    40    40    40    40    40    40    40    40    40    40
 [721]    40    40    40    40    40    40    40    40    40    40    40    40
 [733]    40    40    40    40    40    40    40    40    40    40    40    40
 [745]    40    40    40    40    40    40    40    40    40    40    40    40
 [757]    40    40    40    40    40    40    40    40    40    40    40    40
 [769]    40    40    40    40    40    40    40    40    40    40    40    40
 [781] -2800    40    40    40    40    40    40    40    40    40    40    40
 [793]    40    40    40    40    40    40    40    40    40    40    40    40
 [805]    40    40    40    40    40    40    40    40    40    40    40    40
 [817]    40    40    40    40    40    40    40    40    40    40    40    40
 [829]    40    40    40    40    40    40    40    40    40    40    40    40
 [841]    40    40    40    40    40    40    40    40    40    40    40 -2800
 [853]    40    40    40    40    40    40    40    40    40    40    40    40
 [865]    40    40    40    40    40    40    40    40    40    40    40    40
 [877]    40    40    40    40    40    40    40    40    40    40    40    40
 [889]    40    40    40    40    40    40    40    40    40    40    40    40
 [901]    40    40    40    40    40    40    40    40    40    40    40    40
 [913]    40    40    40    40    40    40    40    40    40    40 -2800    40
 [925]    40    40    40    40    40    40    40    40    40    40    40    40
 [937]    40    40    40    40    40    40    40    40    40    40    40    40
 [949]    40    40    40    40    40    40    40    40    40    40    40    40
 [961]    40    40    40    40    40    40    40    40    40    40    40    40
 [973]    40    40    40    40    40    40    40    40    40    40    40    40
 [985]    40    40    40    40    40    40    40    40    40 -2800    40    40
 [997]    40    40    40    40    40    40    40    40    40    40    40    40
[1009]    40    40    40    40    40    40    40    40    40    40    40    40
[1021]    40    40    40    40    40    40    40    40    40    40    40    40
[1033]    40    40    40    40    40    40    40    40    40    40    40    40
[1045]    40    40    40    40    40    40    40    40    40    40    40    40
[1057]    40    40    40    40    40    40    40    40 -2800    40    40    40
[1069]    40    40    40    40    40    40    40    40    40    40    40    40
[1081]    40    40    40    40    40    40    40    40    40    40    40    40
[1093]    40    40    40    40    40    40    40    40    40    40    40    40
[1105]    40    40    40    40    40    40    40    40    40    40    40    40
[1117]    40    40    40    40    40    40    40    40    40    40    40    40
[1129]    40    40    40    40    40    40    40 -2800    40    40    40    40
[1141]    40    40    40    40    40    40    40    40    40    40    40    40
[1153]    40    40    40    40    40    40    40    40    40    40    40    40
[1165]    40    40    40    40    40    40    40    40    40    40    40    40
[1177]    40    40    40    40    40    40    40    40    40    40    40    40
[1189]    40    40    40    40    40    40    40    40    40    40    40    40
[1201]    40    40    40    40    40    40
> diff(grd$y)
   [1]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
  [25]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
  [49]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0
  [73]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
  [97]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [121]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0
 [145]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [169]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [193]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0
 [217]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [241]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [265]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0
 [289]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [313]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [337]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0
 [361]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [385]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [409]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0
 [433]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [457]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [481]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0
 [505]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [529]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [553]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0
 [577]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [601]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [625]  0  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0
 [649]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [673]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [697]  0  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0  0
 [721]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [745]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [769]  0  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0  0  0
 [793]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [817]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [841]  0  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0  0  0  0
 [865]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [889]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [913]  0  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0  0  0  0  0
 [937]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [961]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
 [985]  0  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1009]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1033]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1057]  0  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1081]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1105]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1129]  0  0  0  0  0  0  0 40  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1153]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1177]  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
[1201]  0  0  0  0  0  0
> summary(grd)
       x                y         
 Min.   :178600   Min.   :329680  
 1st Qu.:179280   1st Qu.:329840  
 Median :180000   Median :330000  
 Mean   :180000   Mean   :330000  
 3rd Qu.:180720   3rd Qu.:330160  
 Max.   :181400   Max.   :330320  
> 
> 
> 
> cleanEx(); ..nameEx <- "meuse"
> 
> ### * meuse
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: meuse
> ### Title: Meuse river data set
> ### Aliases: meuse
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(meuse)
> summary(meuse)
       x                y             cadmium           copper      
 Min.   :178605   Min.   :329714   Min.   : 0.200   Min.   : 14.00  
 1st Qu.:179371   1st Qu.:330762   1st Qu.: 0.800   1st Qu.: 23.00  
 Median :179991   Median :331633   Median : 2.100   Median : 31.00  
 Mean   :180005   Mean   :331635   Mean   : 3.246   Mean   : 40.32  
 3rd Qu.:180630   3rd Qu.:332463   3rd Qu.: 3.850   3rd Qu.: 49.50  
 Max.   :181390   Max.   :333611   Max.   :18.100   Max.   :128.00  
                                                                    
      lead            zinc             elev             dist        
 Min.   : 37.0   Min.   : 113.0   Min.   : 5.180   Min.   :0.00000  
 1st Qu.: 72.5   1st Qu.: 198.0   1st Qu.: 7.546   1st Qu.:0.07569  
 Median :123.0   Median : 326.0   Median : 8.180   Median :0.21184  
 Mean   :153.4   Mean   : 469.7   Mean   : 8.165   Mean   :0.24002  
 3rd Qu.:207.0   3rd Qu.: 674.5   3rd Qu.: 8.955   3rd Qu.:0.36407  
 Max.   :654.0   Max.   :1839.0   Max.   :10.520   Max.   :0.88039  
                                                                    
       om         ffreq  soil   lime       landuse       dist.m      
 Min.   : 1.000   1:84   1:97   0:111   W      :50   Min.   :  10.0  
 1st Qu.: 5.300   2:48   2:46   1: 44   Ah     :39   1st Qu.:  80.0  
 Median : 6.900   3:23   3:12           Am     :22   Median : 270.0  
 Mean   : 7.478                         Fw     :10   Mean   : 290.3  
 3rd Qu.: 9.000                         Ab     : 8   3rd Qu.: 450.0  
 Max.   :17.000                         (Other):25   Max.   :1000.0  
 NA's   : 2.000                         NA's   : 1                   
> 
> 
> 
> cleanEx(); ..nameEx <- "meuse.all"
> 
> ### * meuse.all
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: meuse.all
> ### Title: Meuse river data set - original, full data set
> ### Aliases: meuse.all
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(meuse.all)
> summary(meuse.all)
     sample             x                y             cadmium      
 Min.   :  1.00   Min.   :178605   Min.   :329714   Min.   : 0.000  
 1st Qu.: 41.75   1st Qu.:179358   1st Qu.:330771   1st Qu.: 0.800  
 Median : 82.50   Median :179945   Median :331558   Median : 1.900  
 Mean   : 82.50   Mean   :179989   Mean   :331614   Mean   : 3.109  
 3rd Qu.:123.25   3rd Qu.:180626   3rd Qu.:332410   3rd Qu.: 3.725  
 Max.   :164.00   Max.   :181390   Max.   :333611   Max.   :18.100  
                                                                    
     copper            lead             zinc             elev       
 Min.   : 14.00   Min.   : 27.00   Min.   : 107.0   Min.   : 0.000  
 1st Qu.: 23.00   1st Qu.: 68.75   1st Qu.: 191.8   1st Qu.: 7.390  
 Median : 29.50   Median :116.00   Median : 307.5   Median : 8.124  
 Mean   : 39.42   Mean   :148.55   Mean   : 464.6   Mean   : 7.775  
 3rd Qu.: 48.00   3rd Qu.:201.75   3rd Qu.: 662.5   3rd Qu.: 8.915  
 Max.   :128.00   Max.   :654.00   Max.   :1839.0   Max.   :10.520  
                                                                    
     dist.m             om             ffreq            soil      
 Min.   :  10.0   Min.   : 1.000   Min.   :1.000   Min.   :1.000  
 1st Qu.:  80.0   1st Qu.: 5.000   1st Qu.:1.000   1st Qu.:1.000  
 Median : 270.0   Median : 6.550   Median :1.000   Median :1.000  
 Mean   : 294.2   Mean   : 7.291   Mean   :1.604   Mean   :1.463  
 3rd Qu.: 450.0   3rd Qu.: 8.950   3rd Qu.:2.000   3rd Qu.:2.000  
 Max.   :1000.0   Max.   :17.000   Max.   :3.000   Max.   :3.000  
                  NA's   : 2.000                                  
      lime           landuse     in.pit        in.meuse155      in.BMcD       
 Min.   :0.0000   W      :54   Mode :logical   Mode :logical   Mode :logical  
 1st Qu.:0.0000   Ah     :42   FALSE:156       FALSE:9         FALSE:66       
 Median :0.0000   Am     :22   TRUE :8         TRUE :155       TRUE :98       
 Mean   :0.2988   Fw     :10                                                  
 3rd Qu.:1.0000   Ab     : 8                                                  
 Max.   :1.0000   (Other):27                                                  
                  NA's   : 1                                                  
> 
> 
> 
> cleanEx(); ..nameEx <- "meuse.alt"
> 
> ### * meuse.alt
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: meuse.alt
> ### Title: Meuse river altitude data set
> ### Aliases: meuse.alt
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(meuse.alt)
> library(lattice)
> xyplot(y~x, meuse.alt, aspect = "iso")
> 
> 
> 
> cleanEx(); ..nameEx <- "meuse.grid"
> 
> ### * meuse.grid
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: meuse.grid
> ### Title: Prediction Grid for Meuse Data Set
> ### Aliases: meuse.grid
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(meuse.grid)
> library(lattice)
> xyplot(y~x, meuse.grid, asp=mapasp(meuse.grid), pch="+")
> 
> 
> 
> cleanEx(); ..nameEx <- "ncp.grid"
> 
> ### * ncp.grid
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: ncp.grid
> ### Title: Grid for the NCP, the Dutch part of the North Sea
> ### Aliases: ncp.grid
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(ncp.grid)
> summary(ncp.grid)
       x                y               depth           coast       
 Min.   :466500   Min.   :5699000   Min.   : 1.00   Min.   :  1.00  
 1st Qu.:531500   1st Qu.:5859000   1st Qu.:25.00   1st Qu.: 38.00  
 Median :566500   Median :5954000   Median :31.00   Median : 76.00  
 Mean   :572625   Mean   :5940147   Mean   :31.79   Mean   : 89.72  
 3rd Qu.:606500   3rd Qu.:6024000   3rd Qu.:40.00   3rd Qu.:136.00  
 Max.   :736500   Max.   :6129000   Max.   :59.00   Max.   :248.00  
      area      
 Min.   : 1.00  
 1st Qu.: 1.00  
 Median : 1.00  
 Mean   : 1.87  
 3rd Qu.: 2.00  
 Max.   :19.00  
> 
> 
> 
> cleanEx(); ..nameEx <- "ossfim"
> 
> ### * ossfim
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: ossfim
> ### Title: Kriging standard errors as function of grid spacing and block
> ###   size
> ### Aliases: ossfim
> ### Keywords: models
> 
> ### ** Examples
> 
> x <- ossfim(1:15,1:15, model = vgm(1,"Exp",15))
> library(lattice)
> levelplot(kriging.se~spacing+block.size, x, 
+   main = "Ossfim results, variogram 1 Exp(15)")
> # if you wonder about the decrease in the upper left corner of the graph,
> # try the above with nmax set to 100, or perhaps 200.
> 
> 
> 
> cleanEx(); ..nameEx <- "oxford"
> 
> ### * oxford
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: oxford
> ### Title: Oxford soil samples
> ### Aliases: oxford
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(oxford)
> summary(oxford)
    PROFILE           XCOORD        YCOORD          ELEV       PROFCLASS
 Min.   :  1.00   Min.   :100   Min.   : 100   Min.   :540.0   Cr:19    
 1st Qu.: 32.25   1st Qu.:200   1st Qu.: 600   1st Qu.:558.0   Ct:36    
 Median : 63.50   Median :350   Median :1100   Median :573.0   Ia:71    
 Mean   : 63.50   Mean   :350   Mean   :1100   Mean   :573.6            
 3rd Qu.: 94.75   3rd Qu.:500   3rd Qu.:1600   3rd Qu.:584.5            
 Max.   :126.00   Max.   :600   Max.   :2100   Max.   :632.0            
 MAPCLASS      VAL1            CHR1           LIME1            VAL2     
 Cr:31    Min.   :2.000   Min.   :1.000   Min.   :0.000   Min.   :4.00  
 Ct:36    1st Qu.:3.000   1st Qu.:2.000   1st Qu.:1.000   1st Qu.:4.00  
 Ia:59    Median :4.000   Median :2.000   Median :4.000   Median :8.00  
          Mean   :3.508   Mean   :2.468   Mean   :2.643   Mean   :6.23  
          3rd Qu.:4.000   3rd Qu.:3.000   3rd Qu.:4.000   3rd Qu.:8.00  
          Max.   :4.000   Max.   :4.000   Max.   :4.000   Max.   :8.00  
      CHR2       LIME2          DEPTHCM         DEP2LIME         PCLAY1     
 Min.   :2   Min.   :0.000   Min.   :10.00   Min.   :20.00   Min.   :10.00  
 1st Qu.:2   1st Qu.:4.000   1st Qu.:25.00   1st Qu.:20.00   1st Qu.:20.00  
 Median :2   Median :5.000   Median :36.00   Median :20.00   Median :24.50  
 Mean   :3   Mean   :3.889   Mean   :46.25   Mean   :30.32   Mean   :24.44  
 3rd Qu.:4   3rd Qu.:5.000   3rd Qu.:64.75   3rd Qu.:40.00   3rd Qu.:28.00  
 Max.   :6   Max.   :5.000   Max.   :91.00   Max.   :90.00   Max.   :37.00  
     PCLAY2           MG1              OM1              CEC1      
 Min.   :10.00   Min.   : 19.00   Min.   : 2.600   Min.   : 7.00  
 1st Qu.:10.00   1st Qu.: 44.00   1st Qu.: 4.100   1st Qu.:12.00  
 Median :10.00   Median : 72.00   Median : 5.350   Median :15.00  
 Mean   :14.76   Mean   : 93.53   Mean   : 5.995   Mean   :18.88  
 3rd Qu.:20.00   3rd Qu.:123.25   3rd Qu.: 7.175   3rd Qu.:25.25  
 Max.   :40.00   Max.   :308.00   Max.   :13.100   Max.   :43.00  
      PH1            PHOS1             POT1      
 Min.   :4.200   Min.   : 1.700   Min.   : 83.0  
 1st Qu.:7.200   1st Qu.: 6.200   1st Qu.:127.0  
 Median :7.500   Median : 8.500   Median :164.0  
 Mean   :7.152   Mean   : 8.752   Mean   :181.7  
 3rd Qu.:7.600   3rd Qu.:10.500   3rd Qu.:194.8  
 Max.   :7.700   Max.   :25.000   Max.   :847.0  
> 
> 
> 
> cleanEx(); ..nameEx <- "pcb"
> 
> ### * pcb
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: pcb
> ### Title: PCB138 measurements in sediment at the NCP, the Dutch part of
> ###   the North Sea
> ### Aliases: pcb
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(pcb)
> library(lattice)
> xyplot(y~x|as.factor(yf), pcb, aspect = "iso")
> # demo(pcb)
> 
> 
> 
> cleanEx(); ..nameEx <- "plot.gstatVariogram"
> 
> ### * plot.gstatVariogram
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: plot.gstatVariogram
> ### Title: Plot a Sample Variogram
> ### Aliases: plot.gstatVariogram plot.variogramMap
> ### Keywords: dplot
> 
> ### ** Examples
> 
> data(meuse)
> vgm1 <- variogram(log(zinc)~1, ~x+y, meuse)
> plot(vgm1)
> model.1 <- fit.variogram(vgm1,vgm(1,"Sph",300,1))
> plot(vgm1, model=model.1)
> plot(vgm1, plot.numbers = TRUE, pch = "+")
> vgm2 <- variogram(log(zinc)~1, ~x+y, meuse, alpha=c(0,45,90,135))
> plot(vgm2)
> # the following demonstrates plotting of directional models:
> model.2 <- vgm(.59,"Sph",926,.06,anis=c(0,0.3))
> plot(vgm2, model=model.2)
>  
> g = gstat(id="zinc < 200", form=I(zinc<200)~1,loc=~x+y,data=meuse)
> g = gstat(g, id="zinc < 400", form=I(zinc<400)~1,loc=~x+y,data=meuse)
> g = gstat(g, id="zinc < 800", form=I(zinc<800)~1,loc=~x+y,data=meuse)
> # calculate multivariable, directional variogram:
> v = variogram(g, alpha=c(0,45,90,135))
> plot(v, group.id = FALSE, auto.key = TRUE) # id and id pairs panels
> plot(v, group.id = TRUE, auto.key = TRUE)  # direction panels
> 
> if (require(sp)) {
+   plot(variogram(g, cutoff=1000, width=100, map=TRUE),
+         main = "(cross) semivariance maps")
+   plot(variogram(g, cutoff=1000, width=100, map=TRUE), np=TRUE,
+         main = "number of point pairs")
+ }
Loading required package: sp

Attaching package: 'sp'


	The following object(s) are masked from package:gstat :

	 bpy.colors bubble makegrid map.to.lev mapasp point.in.polygon select.spatial zerodist 

> 
> 
> 
> 
> cleanEx(); ..nameEx <- "plot.pointPairs"
> 
> ### * plot.pointPairs
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: plot.pointPairs
> ### Title: Plot a point pairs, identified from a variogram cloud
> ### Aliases: plot.pointPairs
> ### Keywords: dplot
> 
> ### ** Examples
> 
> ### The following requires interaction, and is therefore outcommented
> #data(meuse)
> #vgm1 <- variogram(log(zinc)~1, ~x+y, meuse, cloud = TRUE)
> #pp <- plot(vgm1, id = TRUE)
> ### Identify the point pairs
> #plot(pp, data = meuse) # meuse has x and y as coordinates
> 
> 
> 
> cleanEx(); ..nameEx <- "plot.variogramCloud"
> 
> ### * plot.variogramCloud
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: plot.variogramCloud
> ### Title: Plot and Identify Data Pairs on Sample Variogram Cloud
> ### Aliases: plot.variogramCloud
> ### Keywords: dplot
> 
> ### ** Examples
> 
> data(meuse)
> plot(variogram(log(zinc)~1, loc=~x+y, data=meuse, cloud=TRUE))
> ## commands that require interaction:
> # x <- variogram(log(zinc)~1, loc=~x+y, data=meuse, cloud=TRUE)
> # plot(plot(x, idendify = TRUE), meuse)
> # plot(plot(x, digitize = TRUE), meuse)
> 
> 
> 
> cleanEx(); ..nameEx <- "point.in.polygon"
> 
> ### * point.in.polygon
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: point.in.polygon
> ### Title: do point(s) fall in a given polygon?
> ### Aliases: point.in.polygon
> ### Keywords: models
> 
> ### ** Examples
> 
> # open polygon:
> point.in.polygon(1:10,1:10,c(3,4,4,3),c(3,3,4,4))
 [1] FALSE FALSE  TRUE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE
> # closed polygon:
> point.in.polygon(1:10,1:10,c(3,4,4,3,3),c(3,3,4,4,3))
 [1] FALSE FALSE  TRUE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE
> 
> 
> 
> cleanEx(); ..nameEx <- "predict.gstat"
> 
> ### * predict.gstat
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: predict.gstat
> ### Title: Multivariable Geostatistical Prediction and Simulation
> ### Aliases: predict.gstat
> ### Keywords: models
> 
> ### ** Examples
> 
> # generate 5 conditional simulations
> data(meuse)
> data(meuse.grid)
> v <- variogram(log(zinc)~1,~x+y, meuse)
> m <- fit.variogram(v, vgm(1, "Sph", 300, 1))
> plot(v, model = m)
> set.seed(131)
> sim <- krige(formula = log(zinc)~1, locations = ~x+y, model = m, data 
+         = meuse, newdata = meuse.grid, nmax = 15, beta = 5.9, nsim = 5)
[using conditional Gaussian simulation]
> # show all 5 simulation, using map.to.lev to rearrange sim:
> if (require(sp) == FALSE) {
+ library(lattice)
+ levelplot(z~x+y|name, map.to.lev(sim, z=c(3:7)), aspect = mapasp(sim))
+ 
+ # calculate generalised least squares residuals w.r.t. constant trend:
+ g <- gstat(id = "log.zinc", formula = log(zinc)~1, locations = ~x+y, 
+      model = m, data = meuse)
+ blue0 <- predict(g, newdata = meuse, BLUE = TRUE)
+ blue0$blue.res <- log(meuse$zinc) - blue0$log.zinc.pred
+ bubble(blue0, zcol = "blue.res", main = "GLS residuals w.r.t. constant")
+ 
+ # calculate generalised least squares residuals w.r.t. linear trend:
+ m <- fit.variogram(variogram(log(zinc)~sqrt(dist.m),~x+y, meuse), 
+         vgm(1, "Sph", 300, 1))
+ g <- gstat(id = "log.zinc", formula = log(zinc)~sqrt(dist.m), 
+         locations = ~x+y, model = m, data = meuse)
+ blue1 <- predict(g, newdata = meuse, BLUE = TRUE)
+ blue1$blue.res <- log(meuse$zinc) - blue1$log.zinc.pred
+ bubble(blue1, zcol = "blue.res", 
+         main = "GLS residuals w.r.t. linear trend")
+ 
+ # unconditional simulation on a 100 x 100 grid
+ xy <- expand.grid(1:100, 1:100)
+ names(xy) <- c("x","y")
+ g.dummy <- gstat(formula = z~1, locations = ~x+y, dummy = TRUE, beta = 0,
+         model = vgm(1,"Exp",15), nmax = 20)
+ yy <- predict(g.dummy, newdata = xy, nsim = 4)
+ # show one realisation:
+ levelplot(sim1~x+y, yy, aspect = mapasp(yy))
+ # show all four:
+ levelplot(z~x+y|name, map.to.lev(yy, z=c(3:6)), aspect = mapasp(yy))
+ }
Loading required package: sp

Attaching package: 'sp'


	The following object(s) are masked from package:gstat :

	 bpy.colors bubble makegrid map.to.lev mapasp point.in.polygon select.spatial zerodist 

> 
> 
> 
> 
> cleanEx(); ..nameEx <- "select.spatial"
> 
> ### * select.spatial
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: select.spatial
> ### Title: select points spatially
> ### Aliases: select.spatial
> ### Keywords: models
> 
> ### ** Examples
> 
> data(meuse)
> ## the following command requires user interaction: left mouse
> ## selects points, right mouse ends digitizing
> # select.spatial(data=meuse)
> 
> 
> 
> cleanEx(); ..nameEx <- "show.vgms"
> 
> ### * show.vgms
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: show.vgms
> ### Title: Plot Variogram Model Functions
> ### Aliases: show.vgms
> ### Keywords: dplot
> 
> ### ** Examples
> 
> show.vgms()
> show.vgms(models = c("Exp", "Mat", "Gau"), nugget = 0.1)
> # show a set of Matern models with different smoothness:
> show.vgms(kappa.range = c(.1, .2, .5, 1, 2, 5, 10), max = 10)
> # show a set of Exponential class models with different shape parameter:
> show.vgms(kappa.range = c(.05, .1, .2, .5, 1, 1.5, 1.8, 1.9, 2), models = "Exc", max = 10)
> 
> 
> 
> cleanEx(); ..nameEx <- "sic2004"
> 
> ### * sic2004
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: sic2004
> ### Title: Spatial Interpolation Comparison 2004 data set: Natural Ambient
> ###   Radioactivity
> ### Aliases: sic2004 sic.train sic.pred sic.grid sic.test sic.val
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(sic2004) 
> # FIGURE 1. Locations of the 200 monitoring stations for the 11 data sets. 
> # The values taken by the variable are known.
> plot(y~x,sic.train,pch=1,col="red", asp=1)
> 
> # FIGURE 2. Locations of the 808 remaining monitoring stations at which 
> # the values of the variable must be estimated.
> plot(y~x,sic.pred,pch="?", asp=1, cex=.8) # Figure 2
> 
> # FIGURE 3. Locations of the 1008 monitoring stations (exhaustive data sets). 
> # Red circles are used to estimate values located at the questions marks
> plot(y~x,sic.train,pch=1,col="red", asp=1)
> points(y~x, sic.pred, pch="?", cex=.8)
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "variogram"
> 
> ### * variogram
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: variogram
> ### Title: Calculate Sample or Residual Variogram or Variogram Cloud
> ### Aliases: variogram variogram.gstat variogram.formula variogram.default
> ###   print.gstatVariogram print.variogramCloud
> ### Keywords: models
> 
> ### ** Examples
> 
> data(meuse)
> # no trend:
> variogram(log(zinc)~1, loc=~x+y, meuse)
    np       dist     gamma dir.hor dir.ver   id
1   57   79.29244 0.1234479       0       0 var1
2  299  163.97367 0.2162185       0       0 var1
3  419  267.36483 0.3027859       0       0 var1
4  457  372.73542 0.4121448       0       0 var1
5  547  478.47670 0.4634128       0       0 var1
6  533  585.34058 0.5646933       0       0 var1
7  574  693.14526 0.5689683       0       0 var1
8  564  796.18365 0.6186769       0       0 var1
9  589  903.14650 0.6471479       0       0 var1
10 543 1011.29177 0.6915705       0       0 var1
11 500 1117.86235 0.7033984       0       0 var1
12 477 1221.32810 0.6038770       0       0 var1
13 452 1329.16407 0.6517158       0       0 var1
14 457 1437.25620 0.5665318       0       0 var1
15 415 1543.20248 0.5748227       0       0 var1
> # residual variogram w.r.t. a linear trend:
> variogram(log(zinc)~x+y, loc=~x+y, meuse)
    np       dist     gamma dir.hor dir.ver   id
1   57   79.29244 0.1060834       0       0 var1
2  299  163.97367 0.1829983       0       0 var1
3  419  267.36483 0.2264256       0       0 var1
4  457  372.73542 0.2847192       0       0 var1
5  547  478.47670 0.3162418       0       0 var1
6  533  585.34058 0.3571578       0       0 var1
7  574  693.14526 0.3701742       0       0 var1
8  564  796.18365 0.4201289       0       0 var1
9  589  903.14650 0.4216983       0       0 var1
10 543 1011.29177 0.4772549       0       0 var1
11 500 1117.86235 0.5075874       0       0 var1
12 477 1221.32810 0.4617632       0       0 var1
13 452 1329.16407 0.5512305       0       0 var1
14 457 1437.25620 0.4352155       0       0 var1
15 415 1543.20248 0.4556815       0       0 var1
> # directional variogram:
> variogram(log(zinc)~x+y, loc=~x+y, meuse, alpha=c(0,45,90,135))
    np       dist      gamma dir.hor dir.ver   id
1   12   84.36080 0.04114593       0       0 var1
2   76  165.59800 0.19091543       0       0 var1
3  109  270.29441 0.21867508       0       0 var1
4  134  371.27824 0.23112878       0       0 var1
5  158  478.06480 0.38337565       0       0 var1
6  154  583.35601 0.35513567       0       0 var1
7  159  692.50911 0.35709265       0       0 var1
8  158  797.52941 0.46221222       0       0 var1
9  156  901.86529 0.47081724       0       0 var1
10 156 1011.55318 0.50937290       0       0 var1
11 137 1115.24492 0.57358764       0       0 var1
12 135 1220.31674 0.43193998       0       0 var1
13 109 1328.07859 0.68882673       0       0 var1
14 120 1436.93237 0.53015452       0       0 var1
15  96 1544.68559 0.66909962       0       0 var1
16  11   82.06663 0.07619858      45       0 var1
17  91  165.75829 0.11957011      45       0 var1
18 118  266.93093 0.20557549      45       0 var1
19 136  374.24886 0.27864922      45       0 var1
20 172  479.40618 0.23932562      45       0 var1
21 177  587.53554 0.28038440      45       0 var1
22 209  693.02620 0.34028114      45       0 var1
23 226  796.37554 0.37201935      45       0 var1
24 283  905.25038 0.36146985      45       0 var1
25 264 1012.26326 0.36891951      45       0 var1
26 274 1121.20926 0.36831067      45       0 var1
27 275 1221.63704 0.33875319      45       0 var1
28 282 1330.93431 0.33848846      45       0 var1
29 297 1438.21262 0.31476883      45       0 var1
30 299 1542.75515 0.31707228      45       0 var1
31  16   78.75466 0.07583160      90       0 var1
32  70  160.01667 0.20149652      90       0 var1
33  97  267.68973 0.20686187      90       0 var1
34  98  372.02688 0.28167260      90       0 var1
35 118  479.76226 0.30366429      90       0 var1
36  98  585.85589 0.46344817      90       0 var1
37 115  691.04342 0.36401272      90       0 var1
38 100  796.22142 0.36912878      90       0 var1
39  88  901.26201 0.50261434      90       0 var1
40  72 1004.66642 0.56369456      90       0 var1
41  68 1109.43463 0.77219638      90       0 var1
42  51 1223.73294 0.79679699      90       0 var1
43  44 1322.80887 0.82262644      90       0 var1
44  30 1430.99001 0.80073011      90       0 var1
45  16 1544.27842 1.17421050      90       0 var1
46  18   74.69621 0.19452856     135       0 var1
47  62  163.83075 0.24550456     135       0 var1
48  95  264.21071 0.28119200     135       0 var1
49  89  373.39690 0.37803627     135       0 var1
50  99  475.98691 0.35772223     135       0 var1
51 104  584.05805 0.39065627     135       0 var1
52  91  697.18636 0.46947283     135       0 var1
53  80  792.93648 0.53667425     135       0 var1
54  62  899.44175 0.45817366     135       0 var1
55  51 1014.81674 0.81777411     135       0 var1
56  21 1118.55839 1.03741404     135       0 var1
57  16 1216.88607 1.75971197     135       0 var1
58  17 1323.20745 2.49557308     135       0 var1
59  10 1431.53529 1.77666963     135       0 var1
60   4 1536.74264 2.82057119     135       0 var1
> 
> # GLS residual variogram:
> v = variogram(log(zinc)~x+y,~x+y,meuse)
> v.fit = fit.variogram(v, vgm(1, "Sph", 700, 1))
> v.fit
  model     psill    range
1   Nug 0.0823414    0.000
2   Sph 0.3886635 1098.557
> set = list(gls=1)
> v
    np       dist     gamma dir.hor dir.ver   id
1   57   79.29244 0.1060834       0       0 var1
2  299  163.97367 0.1829983       0       0 var1
3  419  267.36483 0.2264256       0       0 var1
4  457  372.73542 0.2847192       0       0 var1
5  547  478.47670 0.3162418       0       0 var1
6  533  585.34058 0.3571578       0       0 var1
7  574  693.14526 0.3701742       0       0 var1
8  564  796.18365 0.4201289       0       0 var1
9  589  903.14650 0.4216983       0       0 var1
10 543 1011.29177 0.4772549       0       0 var1
11 500 1117.86235 0.5075874       0       0 var1
12 477 1221.32810 0.4617632       0       0 var1
13 452 1329.16407 0.5512305       0       0 var1
14 457 1437.25620 0.4352155       0       0 var1
15 415 1543.20248 0.4556815       0       0 var1
> g = gstat(NULL, "log-zinc", log(zinc)~x+y,~x+y, meuse, model=v.fit, set = set)
> variogram(g)
    np       dist     gamma dir.hor dir.ver       id
1   57   79.29244 0.1059824       0       0 log-zinc
2  299  163.97367 0.1826061       0       0 log-zinc
3  419  267.36483 0.2256106       0       0 log-zinc
4  457  372.73542 0.2839249       0       0 log-zinc
5  547  478.47670 0.3156089       0       0 log-zinc
6  533  585.34058 0.3566521       0       0 log-zinc
7  574  693.14526 0.3686389       0       0 log-zinc
8  564  796.18365 0.4203338       0       0 log-zinc
9  589  903.14650 0.4212183       0       0 log-zinc
10 543 1011.29177 0.4766290       0       0 log-zinc
11 500 1117.86235 0.5089493       0       0 log-zinc
12 477 1221.32810 0.4637837       0       0 log-zinc
13 452 1329.16407 0.5501710       0       0 log-zinc
14 457 1437.25620 0.4388562       0       0 log-zinc
15 415 1543.20248 0.4580369       0       0 log-zinc
> 
> 
> 
> cleanEx(); ..nameEx <- "variogram.line"
> 
> ### * variogram.line
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: variogram.line
> ### Title: Semivariance Values For a Given Variogram Model
> ### Aliases: variogram.line
> ### Keywords: models
> 
> ### ** Examples
> 
> variogram.line(vgm(5, "Exp", 10, 5), 10, 10)
       dist    gamma
1   0.00001 5.000005
2   1.11112 5.525807
3   2.22223 5.996316
4   3.33334 6.417346
5   4.44445 6.794100
6   5.55556 7.131234
7   6.66667 7.432915
8   7.77778 7.702871
9   8.88889 7.944439
10 10.00000 8.160603
> # anisotropic variogram, plotted in E-W direction:
> variogram.line(vgm(1, "Sph", 10, anis=c(0,0.5)), 10, 10)
       dist     gamma
1   0.00001 0.0000030
2   1.11112 0.3278489
3   2.22223 0.6227728
4   3.33334 0.8518530
5   4.44445 0.9821677
6   5.55556 1.0000000
7   6.66667 1.0000000
8   7.77778 1.0000000
9   8.88889 1.0000000
10 10.00000 1.0000000
> # anisotropic variogram, plotted in N-S direction:
> variogram.line(vgm(1, "Sph", 10, anis=c(0,0.5)), 10, 10, dir=c(0,1,0))
       dist     gamma
1   0.00001 0.0000015
2   1.11112 0.1659821
3   2.22223 0.3278475
4   3.33334 0.4814824
5   4.44445 0.6227716
6   5.55556 0.7475999
7   6.66667 0.8518521
8   7.77778 0.9314130
9   8.88889 0.9821674
10 10.00000 1.0000000
> 
> 
> 
> cleanEx(); ..nameEx <- "vgm"
> 
> ### * vgm
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: vgm
> ### Title: Generate, or Add to Variogram Model
> ### Aliases: vgm print.variogramModel
> ### Keywords: models
> 
> ### ** Examples
> 
> vgm()
   short                        long
1    Nug                Nug (nugget)
2    Exp           Exp (exponential)
3    Sph             Sph (spherical)
4    Gau              Gau (gaussian)
5    Exc Exclass (Exponential class)
6    Mat                Mat (Matern)
7    Cir              Cir (circular)
8    Lin                Lin (linear)
9    Bes                Bes (bessel)
10   Pen        Pen (pentaspherical)
11   Per              Per (periodic)
12   Hol                  Hol (hole)
13   Log           Log (logarithmic)
14   Pow                 Pow (power)
15   Spl                Spl (spline)
16   Err     Err (Measurement error)
17   Int             Int (Intercept)
> vgm(10, "Exp", 300)
  model psill range
1   Exp    10   300
> x <- vgm(10, "Exp", 300)
> vgm(10, "Nug", 0)
  model psill range
1   Nug    10     0
> vgm(10, "Exp", 300, 4.5)
  model psill range
1   Nug   4.5     0
2   Exp  10.0   300
> vgm(10, "Mat", 300, 4.5, kappa = 0.7)
  model psill range kappa
1   Nug   4.5     0   0.0
2   Mat  10.0   300   0.7
> vgm( 5, "Exp", 300, add.to = vgm(5, "Exp", 60, nugget = 2.5))
  model psill range
1   Nug   2.5     0
2   Exp   5.0    60
3   Exp   5.0   300
> vgm(10, "Exp", 300, anis = c(30, 0.5))
  model psill range ang1 anis1
1   Exp    10   300   30   0.5
> vgm(10, "Exp", 300, anis = c(30, 10, 0, 0.5, 0.3))
  model psill range ang1 ang2 ang3 anis1 anis2
1   Exp    10   300   30   10    0   0.5   0.3
> # Matern variogram model:
> vgm(1, "Mat", 1, kappa=.3)
  model psill range kappa
1   Mat     1     1   0.3
> x <- vgm(0.39527463, "Sph", 953.8942, nugget = 0.06105141)
> x
  model      psill    range
1   Nug 0.06105141   0.0000
2   Sph 0.39527463 953.8942
> print(x, digits = 3);
  model psill range
1   Nug 0.061     0
2   Sph 0.395   954
> # to see all components, do
> print.data.frame(x)
  model      psill    range kappa ang1 ang2 ang3 anis1 anis2
1   Nug 0.06105141   0.0000   0.0    0    0    0     1     1
2   Sph 0.39527463 953.8942   0.5    0    0    0     1     1
> 
> 
> 
> cleanEx(); ..nameEx <- "walker"
> 
> ### * walker
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: walker
> ### Title: Walker Lake sample data set
> ### Aliases: walker
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(walker)
> summary(walker)
       Id              X               Y               V         
 Min.   :  1.0   Min.   :  8.0   Min.   :  8.0   Min.   :   0.0  
 1st Qu.:118.2   1st Qu.: 51.0   1st Qu.: 80.0   1st Qu.: 184.6  
 Median :235.5   Median : 89.0   Median :139.5   Median : 424.0  
 Mean   :235.5   Mean   :111.1   Mean   :141.2   Mean   : 435.3  
 3rd Qu.:352.8   3rd Qu.:170.0   3rd Qu.:208.0   3rd Qu.: 640.9  
 Max.   :470.0   Max.   :251.0   Max.   :291.0   Max.   :1528.1  
                                                                 
       U                 T        
 Min.   :   0.00   Min.   :1.000  
 1st Qu.:  82.15   1st Qu.:2.000  
 Median : 319.30   Median :2.000  
 Mean   : 604.08   Mean   :1.904  
 3rd Qu.: 844.55   3rd Qu.:2.000  
 Max.   :5190.10   Max.   :2.000  
 NA's   : 195.00                  
> 
> 
> 
> cleanEx(); ..nameEx <- "zerodist"
> 
> ### * zerodist
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: zerodist
> ### Title: find point pairs with equal spatial coordinates
> ### Aliases: zerodist
> ### Keywords: dplot
> 
> ### ** Examples
> 
> if (require(sp) == FALSE) {
+ data(meuse)
+ # pick 10 rows
+ n <- 10
+ ran10 <- sample(nrow(meuse), size = n, replace = TRUE)
+ meusedup <- rbind(meuse, meuse[ran10, ])
+ zerodist(meusedup$x, meusedup$y)
+ zd <- zerodist(meusedup$x, meusedup$y)
+ sum(abs(zd[1:n,1] - sort(ran10))) # 0!
+ # remove the duplicate rows:
+ meusedup2 <- meusedup[-zd[,2], ]
+ # find point pairs within 500 m distance of each other:
+ zerodist(meuse$x, meuse$y, 500)
+ }
Loading required package: sp

Attaching package: 'sp'


	The following object(s) are masked from package:gstat :

	 bpy.colors bubble makegrid map.to.lev mapasp point.in.polygon select.spatial zerodist 

> 
> 
> 
> ### * <FOOTER>
> ###
> cat("Time elapsed: ", proc.time() - get("ptime", env = .CheckExEnv),"\n")
Time elapsed:  35.37 0.85 36.87 0 0 
> grDevices::dev.off()
null device 
          1 
> ###
> ### Local variables: ***
> ### mode: outline-minor ***
> ### outline-regexp: "\\(> \\)?### [*]+" ***
> ### End: ***
> quit('no')