SoilGrids-class {GSIF}R Documentation

A class for SoilGrids — soil property and/or class maps

Description

A class containing predictions and prediction error (or multiple realizations) of some of the target global soil property at six standard depths. Standard depths used are based on the GlobalSoilMap.net specifications: sd1 = 2.5 cm (0–5), sd2 = 10 cm (5–15), sd3 = 22.5 cm (15–30), sd4 = 45 cm (30–60), sd5 = 80 cm (60–100), sd6 = 150 cm (100–200).

Slots

varname:

object of class "character"; abbreviated variable name registered in the Global Soil Data registry

TimeSpan:

object of class "list"; contains begin and end of the sampling period of class "POSIXct"

sd1:

object of class "SpatialPixelsDataFrame"; predictions and variances, or number of realizations of the target variable at depth 2.5 cm (0–5)

sd2:

object of class "SpatialPixelsDataFrame"; predictions and variances, or number of realizations of the target variable at depth 10 cm (5–15)

sd3:

object of class "SpatialPixelsDataFrame"; predictions and variances, or number of realizations of the target variable at depth 22.5 cm (15–30)

sd4:

object of class "SpatialPixelsDataFrame"; predictions and variances, or number of realizations of the target variable at depth 45 cm (30–60)

sd5:

object of class "SpatialPixelsDataFrame"; predictions and variances, or number of realizations of the target variable at depth 80 cm (60–100)

sd6:

object of class "SpatialPixelsDataFrame"; predictions and variances, or number of realizations of the target variable at depth 150 cm (100–200)

Gridded data submitted to sd* slots of the "SoilGrids" class must satisfy all of the following requirements (class validity):

Methods

summary

signature(x = "SoilGrids"): generates summary statistics for the object

Author(s)

Tomislav Hengl and Robert A. MacMillan

References

See Also

GlobalSoilMap-class, SpatialComponents-class, geosamples-class

Examples

# load soil samples from the plotKML package: 
library(plotKML)
library(aqp)
library(plyr)
library(splines)
library(rgdal)
library(raster)

data(eberg)
## subset data to 10%:
eberg <- eberg[runif(nrow(eberg)) < .1,]
## sites table:
s.lst <- c("ID", "soiltype", "TAXGRSC", "X", "Y")
h.lst <- c("UHDICM","LHDICM","SNDMHT","SLTMHT","CLYMHT")
sites <- eberg[,s.lst]
## get horizons table:
horizons <- getHorizons(eberg, idcol="ID", sel=h.lst)
## create object of type "SoilProfileCollection"
eberg.spc <- join(horizons, sites, type='inner')
depths(eberg.spc) <- ID ~ UHDICM + LHDICM
site(eberg.spc) <- as.formula(paste("~", paste(s.lst[-1], collapse="+"), sep=""))
coordinates(eberg.spc) <- ~X+Y
proj4string(eberg.spc) <- CRS("+init=epsg:31467")
## convert to logits:
eberg.spc@horizons$SNDMHT.t <- log((eberg.spc@horizons$SNDMHT/100)/
    (1-eberg.spc@horizons$SNDMHT/100))
## convert to geosamples:
eberg.geo <- as.geosamples(eberg.spc)
## load gridded data:
data(eberg_grid)
gridded(eberg_grid) <- ~x+y
proj4string(eberg_grid) <- CRS("+init=epsg:31467")
## derive spc's:
formulaString <- ~ PRMGEO6+DEMSRT6+TWISRT6+TIRAST6
eberg_spc <- spc(eberg_grid, formulaString)
## build a 3D "gstatModel": 
glm.formulaString = as.formula(paste("SNDMHT.t ~ ", 
     paste(names(eberg_spc@predicted), collapse="+"), "+ ns(altitude, df=4)"))
## Not run: 
SNDMHT.m <- fit.gstatModel(observations=eberg.geo, glm.formulaString, 
     covariates=eberg_spc@predicted)
summary(SNDMHT.m@regModel)
SNDMHT.m@vgmModel
## prepare new locations (6 standard depths): 
new3D <- sp3D(eberg_spc@predicted)
## Make predictions at six depths:
sd.l <- lapply(new3D, FUN=function(x){predict(SNDMHT.m, predictionLocations=x, nfold=0)})
## back-transformation function:
invlogit = function(x){exp(x)/(1+exp(x))*100}
## for the back-transformation for the mean value see Diggle and Ribeiro, 2007, p. 148:
invlogit.m = function(x, v){((1+exp(-x))^(-1)-.5*v*exp(-x)*(1-exp(-x))*(1+exp(-x))^(-3) )*100}
## back-transform values from logits:
for(j in 1:length(sd.l)){ 
    sd.l[[j]]@predicted$M <- round(invlogit.m(sd.l[[j]]@predicted$SNDMHT.t,
       sd.l[[j]]@predicted$var1.var))
    sd.l[[j]]@predicted$L <- round(invlogit(sd.l[[j]]@predicted$SNDMHT.t
     - 1.645*sqrt(sd.l[[j]]@predicted$var1.var)))
    sd.l[[j]]@predicted$U <- round(invlogit(sd.l[[j]]@predicted$SNDMHT.t
     + 1.645*sqrt(sd.l[[j]]@predicted$var1.var))) 
}
str(sd.l[[1]]@predicted@data)

## reproject to WGS84 system (100 m resolution):
p = get("cellsize", envir = GSIF.opts)[1]
s = get("stdepths", envir = GSIF.opts)
sd.ll <- sapply(1:length(sd.l), FUN=function(x){ 
     make.3Dgrid(sd.l[[x]]@predicted[c("L","M","U")],
     pixsize=p, stdepths=s[x])})
## save to a "SoilGrids" object:
SNDMHT.gsm <- SoilGrids(obj=sd.ll, varname="SNDPPT", 
             TimeSpan=list(begin="1999-02-01", end="2001-07-01"))
str(SNDMHT.gsm, max.level=2)
## visualize all maps in Google Earth:
data(R_pal)
z0 = mean(eberg_grid$DEMSRT6, na.rm=TRUE)
## export grids:
for(j in 1:length(sd.ll)){
  kml(slot(SNDMHT.gsm, paste("sd", j, sep="")), 
     folder.name = paste("eberg_sd", j, sep=""),
     file = paste("SNDMHT_sd", j, ".kml", sep=""), 
     colour = M, z.lim=c(10,85),
     raster_name = paste("SNDMHT_sd", j, ".png", sep=""), 
     altitude = z0+5000+(s[j]*2500))
}

## End(Not run)

[Package GSIF version 0.5-4 Index]