soilDB Demo: Processing SSURGO Attribute Data with SDA_query()

Submitted by dylan on Thu, 2012-04-26 23:18.

Mapping near Paloma, CAMapping near Paloma, CA This image has nothing to do with the following content.

A quick example of how to use the USDA-NRCS soil data access query facility (SDA), via the soilDB package for R. The following code describes how to get component-level soils data for Yolo County (survey area CA113) from SDA and compute representative sub-order level classification for each map unit. This example requires an understanding of SQL, US Soil Taxonomy and the SSURGO database.

library(soilDB)
library(plyr)

# format query
q <- "SELECT
component.mukey, cokey, comppct_r, compname, taxclname, taxorder, taxsuborder, taxgrtgroup, taxsubgrp
FROM legend
INNER JOIN mapunit ON mapunit.lkey = legend.lkey
LEFT OUTER JOIN component ON component.mukey = mapunit.mukey
WHERE legend.areasymbol = 'CA113'"


# run query, process results, and return as data.frame object
res <- SDA_query(q)


# function for computing total component percentages by strata (suborder)
f.sum <- function(i) {
   n <- nrow(i)
   s <- sum(i$comppct_r)
   return(data.frame(pct=s, n=n))
}

# function for picking the largest suborder from within each map unit (mukey)
f.largest <- function(i) {
   i.sorted <- i[order(i$pct, decreasing=TRUE), ]
   top.suborder <- i.sorted$taxsuborder[1]
   top.suborder.pct <- i.sorted$pct[1]
   return(data.frame(suborder=top.suborder, pct=top.suborder.pct))
}

# tabulate percentage of suborder-level taxa within each map unit (mukey)
comp.suborder.sums <- ddply(res, .(mukey, taxsuborder), f.sum, .progress='text')

# keep the largest suborder, and its associated total percentage
comp.suborder <- ddply(comp.suborder.sums, .(mukey), f.largest, .progress='text')

# tabulate occurance of suborders... not that this does not take into acocunt map unit polygon area
table(comp.suborder$suborder)

( categories: )