Exploring Geomorphic Summaries
D.E. Beaudette
2023-02-26
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 belowres <- 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 |