In ncss-tech/soilReports: R Markdown Reports and Convenience Functions for Soil Survey

    .report.name <- 'shiny-pedon-summary'
    .report.version <- '1.0'
    .report.description <- 'Interactively subset and summarize NASIS pedon data from one or more map units'

knitr::opts_chunk$set(echo = FALSE, message=FALSE, warning = FALSE, fig.align='center', fig.retina=2, dev='png', antialias='cleartype')
#if(demo_mode) {
  # uncomment for stand-alone (no shiny)

  # source("config.R")
  # source("util.R")
  # f <- loafergopher
  # input <- data.frame(1)
  # input$reportname <- "Loafergopher"
  # f$musym <- rep("<missing>", length(loafergopher))
#} else {
  loaded <<- FALSE
  f <- peds

  message(sprintf("generating report for n=%s profiles; peiids:",paste0(profile_id(f), collapse=","), length(f)))
#}

# compute summaries
#s <- summarize.component(f)

# determine max number of profiles:
max.comp.profiles <- s$n

Component Report

r format(Sys.time(), "%Y-%m-%d")

r input$reportname

ranges are (r p.low.rv.high) percentiles

Taxon Names and Pedon Types

Check to make sure that pedons used within this report have been correctly assigned to this component. If not, please fix in NASIS.

wzxhzdk:2

MUSYM Summary

wzxhzdk:3

Hillslope Position Summary

wzxhzdk:4

Geomorphic Component Summaries

wzxhzdk:5

wzxhzdk:6

Geomorphons Summary

this.data <- categorical.prop.table(f$gis_geomorphons)
this.align <- rep('c', times=ncol(this.data) + 1)

print(xtable(this.data, align=this.align), type='html', include.rownames=FALSE, table.placement="H", caption.placement="top", html.table.attributes='cellpadding="3" cellspacing="3"')

Drainage Class Summary

wzxhzdk:8

Surface Shape Summary

wzxhzdk:9

Ecosite Summary

wzxhzdk:10

Generalized Horizon Classification

These tables describe the mapping between field-described horizonation (top row) and generalized horizonation (first column). Numbers describe the number of times any given field-described horizon has been allocated to a generalized horizon. If present, values in the "NA" row should be further investigated.

wzxhzdk:11

cols <- c(rev(brewer.pal(11, 'Spectral')))
col.palette <- colorRampPalette(cols)

this.data <- t(s$ct)
this.data[this.data == 0] <- NA

levelplot(this.data, col.regions=col.palette, colorkey=list(tick.number=15), 
          xlab = 'Original Horizon Designation', ylab='GHL', 
          main = 'GHL Assignment Evaluation', 
          scales = list(alternating=3), 
          panel = function(x, y, z, ...) {
            panel.levelplot(x, y, z, ...)
            idx <- which(!is.na(z))
            panel.text(x[idx], y[idx], z[idx], font=2)
            panel.abline(h=seq(from=0.5, to=length(y), by=1), col=grey(0.45))
            panel.abline(v=seq(from=0.5, to=length(x), by=1), col=grey(0.45))
          })

GHL assignment as a network graph.

this.data <- t(s$ct)

# convert contingency table -> adj. matrix
m <- genhzTableToAdjMat(this.data)

# plot using a function from the sharpshootR package
par(mar=c(1,1,1,1))
plotSoilRelationGraph(m, graph.mode = 'directed', edge.arrow.size=0.5, vertex.label.family='sans')

# clay box-whisker plot, grouped by genhz, over-printed with original hz names
# subset data
h.i <- horizons(f)
h.i.sub <- subset(h.i, subset=!is.na(clay), drop=TRUE)
# hack: reset factor levels, to accomodate filtered O horizons
h.i.sub$genhz <- factor(h.i.sub$genhz)

# plotting style
tps <- list(box.umbrella=list(col=grey(0.4)), 
                    box.rectangle=list(col=grey(0.4)), 
                        box.dot=list(col=grey(0.4), cex=0.75), 
                        plot.symbol=list(col=grey(0.4), cex=0.5)
)
# plot
print(bwplot(genhz ~ clay, data=h.i.sub, main=f, par.settings=tps) + layer(panel.text(x=h.i.sub$clay, y=jitter(as.numeric(h.i.sub$genhz), factor=1.5), label=h.i.sub$hzname, cex=0.75, font=2, col='RoyalBlue')))

Maximum-Likelihood Horizonation

The figure below describes the most likely horizonation, based on the collection of pedons associated with this component. This is only an estimate, expert knowledge should be used to adjust these values as needed. When pedon numbers are low or horizonation is not consistent, overlap can occur. Values in square brackets are related to Brier Scores, smaller values suggest more consistent horizonation within the collection.

trellis.par.set(list(superpose.line=list(lwd=2)))
print(s$ml.hz.plot)

Slice-Wise Plot

    if(length(peds) & length(input$thematic_field)) {
      if(!input$thematic_field %in% 
         c("moist_soil_color","dry_soil_color")) {

        s <- slab(peds, fm = as.formula(sprintf( " ~ %s", 
                                                input$thematic_field)))
        # max depth is <1% contributing fraction
        max.idx <- which(s$contributing_fraction <= 0.01)
        if(!length(max.idx))
          max.idx <- nrow(s)
        s.sub <- s[1:max.idx[1],]
        s.sub <- s.sub[complete.cases(s.sub),]
        #print((s.sub))

        a <- xyplot(top ~ p.q50, data=s.sub, ylab='Depth',
                    xlab=paste0(input$thematic_field,
                         '\nmedian bounded by 5th and 95th percentiles'), 
                    lower=s$p.q5, upper=s$p.q95, 
                    ylim=c(max(s.sub$bottom, na.rm=TRUE),-5),
                    panel=panel.depth_function, 
                    prepanel=prepanel.depth_function,
                    cf=s$contributing_fraction,
                    layout=c(1,1), scales=list(x=list(alternating=1)))

        b <- slab(peds[1,], fm = as.formula(paste0(" ~ ",input$thematic_field)))

        if(length(input$modal_pedon) & !is.na(input$modal_pedon)) {
          modalped <- horizons(peds[which(input$modal_pedon == paste(site(peds)$peiid,
                                                                     site(peds)$upedonid,
                                                                     sep=":")),])

          #TODO: select modal pedon from RV copedon
          modalped$middepth <- modalped$hzdept + (modalped$hzdepb - modalped$hzdept)/2

          foo <- 1:(length(modalped$middepth)*3)
          modalpedplot <- data.frame(foo)

          modalpedplot$y <- 1:length(foo)
          modalpedplot$y[which((foo %% 3) == 1)] <- modalped$hzdept
          modalpedplot$y[which((foo %% 3) == 2)] <- modalped$middepth
          modalpedplot$y[which((foo %% 3) == 0)] <- modalped$hzdepb

          modalpedplot$x <- rep(modalped[[input$thematic_field]],each=3)
          b <- xyplot(y ~ x, data=modalpedplot, type="l",
                      col="RED", lwd=3, ylim=c(250,-5), layout=c(1,1),
                      scales=list(x=list(alternating=1)),
                      par.settings = list(superpose.line = list(lwd=3)))
          (a + as.layer(b))
        } else {
          a
        }
      } else {
        print("Color quantiles by depth coming soon")
      }
    } else return(-1)

Profile Plot

# this resets the default image width according to the number of profiles
opts_chunk$set(fig.width=max.comp.profiles * 1.25)
opts_chunk$set(fig.height=4)

wzxhzdk:18

Texture Class Summary Tables

These tables describe the frequency of textural classes, summarized by component, map unit and generalized horizon. Values within parenthesis are the fraction of horizons associated with each texture class.

wzxhzdk:19

Morphologic Summary Tables

These table describe low-rv-high values for morphologic properties, summarized by component. The low values are the r p.low.rv.high[1] percentile, RV values are the r p.low.rv.high[2] percentile, and the high values are the r p.low.rv.high[3] percentile.

wzxhzdk:20

Aggregate Color Summary, dry

par(mar=c(4.5, 2, 0, 0))
aggregateColorPlot(aggregateColor(f, groups = 'genhz', col = 'dry_soil_color'), label.font = 2, label.cex = 0.95, print.n.hz = TRUE)

Aggregate Color Summary, moist

par(mar=c(4.5, 2, 0, 0))
aggregateColorPlot(aggregateColor(f, groups = 'genhz', col = 'moist_soil_color'), label.font = 2, label.cex = 0.95, print.n.hz = TRUE)

Morphologic Summary by Map Unit

Whiskers extend from the 5th to 95th percentiles, the body represents the 25th through 75th percentiles, and the dot is the 50th percentile.

print(s$pmg)

Surface Fragment Summary Tables

These table describe low-rv-high values for surface rock fragments, summarized by component and map unit. The low values are the r p.low.rv.high[1] percentile, RV values are the r p.low.rv.high[2] percentile, and the high values are the r p.low.rv.high[3] percentile.

s$sf

Diagnostic feature summary

The low values are the r p.low.rv.high[1] percentile, RV values are the r p.low.rv.high[2] percentile, and the high values are the r p.low.rv.high[3] percentile.

s$dt

diagnosticPropertyPlot2(f, v=c('lithic.contact', 'paralithic.contact', 'argillic.horizon', 'cambic.horizon', 'ochric.epipedon', 'mollic.epipedon', 'very.shallow', 'shallow', 'mod.deep', 'deep', 'very.deep'), k=3)

Pedon GIS Summary

The low values are the r p.low.rv.high[1] percentile, RV values are the r p.low.rv.high[2] percentile, and the high values are the r p.low.rv.high[3] percentile. These values were sampled from raster data sources, at each pedon location. Arrows on the circular histogram of field-measured aspect values are related to percentiles and "mean resultant length", on a circular basis. Grey arrows are the r p.low.rv.high[1] and r p.low.rv.high[3] percentiles and the red arrow is the r p.low.rv.high[2] percentile. Longer arrows suggest an aspect-affected pattern or aspect-biased sampling site selection.

s$pg

# this resets the default image width according to the number of profiles
opts_chunk$set(fig.width=4.5)

par(mar=c(0,0,0,0))
aspect.plot(f$aspect, q=p.low.rv.high, plot.title=input$reportname, pch=21, bg='RoyalBlue', col='black', arrow.col=c('grey', 'red', 'grey'))

# try(unlink('this.component.Rda'))

---

This document is based on aqp version r utils::packageDescription("aqp", field="Version") and soilDB version r utils::packageDescription("soilDB", field="Version").

ncss-tech/soilReports documentation built on March 6, 2025, 1:15 a.m.