Setup

# stable packages from CRAN
install.packages('aqp')
install.packages('soilDB')
install.packages('sharpshootR')
install.packages('FactoMineR')
install.packages('cowplot')

# latest versions from GitHub
remotes::install_github("ncss-tech/aqp", dependencies=FALSE, upgrade_dependencies=FALSE)
remotes::install_github("ncss-tech/soilDB", dependencies=FALSE, upgrade_dependencies=FALSE)
remotes::install_github("ncss-tech/sharpshootR", dependencies=FALSE, upgrade_dependencies=FALSE)

Hillslope Position / Soil Morphology

library(aqp)
library(soilDB)
library(sharpshootR)

# soil of interest
# https://casoilresource.lawr.ucdavis.edu/gmap/?loc=44.10219,-122.38589,z14
s <- 'peavine'
# s <- 'clarion'

# all siblings
sib <- siblings(s)

# all soils of interest
soils <- c(s, sib$sib$sibling)

# get morphology + extended summaries
osds <- fetchOSD(soils, extended = TRUE)

There may be missing series from either / both SPC and geomorphic summaries. TBC.

# there may be records missing from SPC / hill slope position
nm <- intersect(profile_id(osd$SPC), osd$hillpos$series)

# keep only those series that exist in both
sub <- subset(osd$SPC, profile_id(osd$SPC) %in% nm)

## inverse problem: extra records in hill slope summaries
# subset hillpos
hillpos.sub <- subset(osd$hillpos, subset = series %in% profile_id(sub))

## TODO: use subset objects below
res <- vizHillslopePosition(osds$hillpos)

# re-order hillslope proportions according to clustering
hp <- osds$hillpos[res$order, ]
nm <- names(hp[, 2:6])
series Toeslope Footslope Backslope Shoulder Summit n shannon_entropy
REMOTE 0.00 0.00 0.53 0.13 0.33 30 1.40
DIGGER 0.00 0.00 0.59 0.19 0.22 64 1.38
FORMADER 0.00 0.00 0.43 0.36 0.21 14 1.53
BOHANNON 0.00 0.00 0.46 0.37 0.18 136 1.49
SHIVIGNY 0.00 0.00 0.59 0.34 0.06 32 1.23
MEDA 0.00 0.25 0.25 0.25 0.25 4 2.00
PEAVINE 0.13 0.13 0.31 0.23 0.21 39 2.24
SLICKROCK 0.03 0.25 0.33 0.15 0.24 72 2.08
PREACHER 0.10 0.21 0.40 0.06 0.23 113 2.07
HONEYGROVE 0.19 0.29 0.26 0.00 0.26 80 1.98
HEMCROSS 0.18 0.38 0.21 0.00 0.23 39 1.93
PANTHER 0.18 0.82 0.00 0.00 0.00 146 0.68
MINNIECE 0.57 0.43 0.00 0.00 0.00 7 0.99
print(res$fig)

par(mar=c(0,0,1,1))
SoilTaxonomyDendrogram(osds$SPC, scaling.factor = 0.015, width = 0.28, name.style = 'center-center')

# colors
par(mar = c(1, 1, 1, 1))
hp.cols <- RColorBrewer::brewer.pal(n = 5, name = 'Set1')[c(2, 3, 4, 5, 1)]
soilPalette(hp.cols, lab = nm)

Composite Figures

Version 1: Lines

par(mar = c(0.5, 0, 0, 2))
layout(matrix(c(1,2)), widths = c(1,1), heights = c(2,1))
plotProfileDendrogram(osds$SPC, res$clust, dend.y.scale = 3, scaling.factor = 0.012, y.offset = 0.2, width = 0.32, name.style = 'center-center', cex.names = 0.7, shrink = TRUE, cex.id = 0.55)

matplot(y = hp[, 2:6], type = 'b', lty = 1, pch = 16, axes = FALSE, col = hp.cols, xlab = '', ylab = '', xlim = c(0.5, length(osds$SPC) + 1))
# grid(nx = 0, ny = NULL)
axis(side = 4, line = -1, las = 1, cex.axis = 0.7)
# axis(side = 2, line = -3, las = 1, cex.axis = 0.7)
legend('topleft', legend = rev(nm), col = rev(hp.cols), pch = 16, bty = 'n', cex = 0.8, pt.cex = 2, horiz = TRUE, inset = c(0.01, 0.01))
mtext('Probability', side = 2, line = -2, font = 2)

Version 2: Stacked Bars

par(mar = c(0.5, 0, 0, 2))
layout(matrix(c(1,2)), widths = c(1,1), heights = c(2,1))
plotProfileDendrogram(osds$SPC, res$clust, dend.y.scale = 3, scaling.factor = 0.012, y.offset = 0.2, width = 0.32, name.style = 'center-center', cex.names = 0.7, shrink = TRUE, cex.id = 0.55)

sp <- c(1.5, rep(1, times = length(osds$SPC) - 1))
barplot(height = t(as.matrix(hp[, 2:6])), beside = FALSE, width = 0.5, space = sp, col = hp.cols,  axes = FALSE, xlab = '', ylab = '', xlim = c(0.5, length(osds$SPC) + 1), ylim = c(0, 1.2))

legend(x = 0.75, y = 1.2, legend = rev(nm), col = rev(hp.cols), pch = 15, bty = 'n', cex = 0.8, pt.cex = 1.25, horiz = TRUE)
mtext('Probability', side = 2, line = -2, font = 2)

3D Surface Shape

res <- vizSurfaceShape(osds$shape_across, title = 'Surface Shape (Across)')

# re-order proportions according to clustering
ss <- osds$shape_across[res$order, ]
nm <- names(ss[, 2:6])
series Concave Linear Convex Complex Undulating n shannon_entropy
SHIVIGNY 0.00 0.00 1.00 0 0 32 0.00
FORMADER 0.12 0.12 0.76 0 0 17 1.02
DIGGER