soil              package:RandomFields              R Documentation

_S_o_i_l _d_a_t_a _o_f _N_o_r_t_h _B_a_v_a_r_i_a, _G_e_r_m_a_n_y

_D_e_s_c_r_i_p_t_i_o_n:

     Soil physical and chemical data collected on a field in the
     Weissenstaedter Becken, Germany

_U_s_a_g_e:

     data(soil)

_F_o_r_m_a_t:

     This data frame contains the following columns:

     _x._c_o_o_r_d x coordinates given in cm

     _y._c_o_o_r_d y coordinates given in cm

     _n_r number of the samples, which were taken in this order

     _m_o_i_s_t_u_r_e moisture content [Kg/Kg * 100%]

     _N_O_3._N nitrate nitrogen [mg/Kg]

     _T_o_t_a_l._N total nitrogen [mg/Kg]

     _N_H_4._N ammonium nitrogen [mg/Kg]

     _D_O_C dissolved organic carbon [mg/Kg]

     _N_2_0_N nitrous oxide [mg/Kg dried substance]

_D_e_t_a_i_l_s:

     For technical reasons some of the data were obtained as
     differences of two measurements (which are not available anymore).
     Therefore, some of the data have negative values.

_S_o_u_r_c_e:

     The data were collected by Wolfgang Falk, Soil Physics Group,
     <URL: http://www.geo.uni-bayreuth.de/bodenphysik/Welcome.html>,
     University of Bayreuth, Germany.

_R_e_f_e_r_e_n_c_e_s:

     Falk, W. (2000) _Kleinskalige raeumliche Variabilitaet von Lachgas
     und bodenchemischen Parameters [Small Scale Spatial Variability of
     Nitrous Oxide and Pedo-Chemical Parameters]_, Master thesis,
     University of Bayreuth, Germany.

_E_x_a_m_p_l_e_s:

     ################################################################
     ##                                                            ##
     ##         a geostatistical analysis that demonstrates        ##
     ##         features of the package `RandomFields'             ##
     ##                                                            ##
     ##        *** THIS EXAMPLE IS NOT RUN AUTOMATICALLY ***       ##
     ##        ***      SINCE IT IS TIME CONSUMING       ***       ##
     ##                                                            ##
     ################################################################

     ## Not run: 

     data(soil)
     names(soil)
     pts <- soil[,c(1,2)]
     d <- soil$moisture

     ## define some graphical parameters first
     close.screen(close.screen())
     maxbin <- max(dist(pts)) / 2
     (zlim <- range(d))
     cn <- 100
     colour <- rainbow(cn)
     par(bg="white",cex=1, cex.lab=1.4, cex.axis=1.4, mar=c(4.3,4.3,0.8,0.8))
     lu.x <- min(soil$x)
     lu.y <- max(soil$y)
     y <- x <- seq(min(soil$x), max(soil$x), l=121) 

     ## ... and a certain appearance of the legend
     my.legend <- function(lu.x, lu.y, zlim, col, cex=1) {
       ## uses already the legend code of R-1.3.0
       cn <- length(col)
       filler <- vector("character", length=(cn-3)/2)
       legend(lu.x, lu.y, y.i=0.03, x.i=0.1, 
              legend=c(format(zlim[2], dig=2), filler,
              format(mean(zlim), dig=2), filler,
              format(zlim[1], dig=2)),
              lty=1, col=col[length(col):1],cex=cex)
     }

     ## plot the data first
     plot(pts, col=colour[1+99*((d-min(d))/diff(zlim))], pch=16,
          xlab="x [cm]", ylab="y [cm]")
     my.legend(lu.x, lu.y, zlim=zlim, col=colour, cex=1.4)

     ## empirical variogram
     ev <- EmpiricalVariogram(pts, data=d, grid=FALSE,
                              bin=c(-1,seq(0,maxbin,l=30)))

     ## show all models,
     by.eye <- ShowModels(x=x, y=y, emp=ev, col=colour, zlim=zlim,
                          Mean=mean(d))

     ## fit parameters of the whittlematern model by MLE
     mle <- fitvario(x=pts, data=d, model="whittle", par=rep(NA,5))

     ## plot the fitted model and the empirical variogram
     plot(ev$c, ev$emp.var, ylim=c(0,11), ylab="variogram", xlab="lag")
     gx <- seq(0.001, max(ev$c), l=100)
     if(!is.null(by.eye)) lines(gx, Variogram(gx, model=by.eye)) 
     lines(gx, Variogram(gx, model=mle$mle), col=2)
     lines(gx, Variogram(gx, model=mle$lsq), col=3)
     lines(gx, Variogram(gx, model=mle$nlsq), col=4)
     lines(gx, Variogram(gx, model=mle$slsq), col=5)
     legend(120, 4, c("empirical", "by eye", "MLE", "lsq", "sqrt(n) lsq",
                    "sd^-0.5 lsq"),
            lty=c(-1, rep(1, 5)), pch=c(1, rep(-1, 5)),
            col=c(1, 1, 2, 3, 4, 5), cex=1.4)

     ## map of expected values
     k <- Kriging("O", x=x, y=y, grid=TRUE, model=mle$mle, given=pts, data=d)
     image(x, y, k, col=colour, zlim=zlim, xlab="x [cm]", ylab="y [cm]")
     par(bg="white")
     my.legend(lu.x, lu.y, zlim=zlim, col=colour, cex=1.4)

     ## what is the probability that at no point of the
     ## grid given by x and y the moisture is greater than 24 percent?
     cs <- CondSimu("O", x=x, y=y, grid=TRUE, model=mle$mle, given=pts,
                    data=d, n=10) # better n=100 or n=1000
     image(x, y, cs[,,1], col=colour, zlim=zlim, xlab="x [cm]",
           ylab="y [cm]")
     my.legend(lu.x, lu.y, zlim=zlim, col=colour, cex=1.4)
     mean(apply(cs<=24, 3, all)) ## about 40 percent ...

     ## End(Not run)

