getSurprisal: Get surprisal

In soundgen: Sound Synthesis and Acoustic Analysis

Get surprisal

Description

Tracks the unpredictability of spectro-temporal changes in a sound over time, returning continuous contours of Shannon surprisal ($info), Bayesian surprise ($kl for Kullback-Leibler divergence), and autocorrelation-based surprisal ($surprisal). This is an attempt to track auditory salience over time - that is, to identify parts of a sound that are likely to involuntarily attract the listeners' attention.

Usage

getSurprisal(
  x,
  samplingRate = NULL,
  scale = NULL,
  from = NULL,
  to = NULL,
  winSurp = 2000,
  input = c("audSpec", "env", "melspec", "spectrogram", "pspec")[1],
  takeLog = TRUE,
  audSpec_pars = list(nFilters = 8, step = 15, minFreq = 60),
  spec_pars = list(windowLength = 20, step = 20),
  env_pars = list(windowLength = 40, step = 20),
  melfcc_pars = list(windowLength = 20, step = 20, maxfreq = NULL, nbands = NULL),
  method = c("acf", "np")[1],
  sameLagAllFreqs = FALSE,
  weightByAmpl = TRUE,
  weightByPrecision = TRUE,
  onlyPeakAutocor = TRUE,
  rescale = FALSE,
  summaryFun = "mean",
  reportEvery = NULL,
  cores = 1,
  plot = TRUE,
  savePlots = NULL,
  osc = c("none", "linear", "dB")[2],
  heights = c(3, 1),
  ylim = NULL,
  contrast = 0.2,
  brightness = 0,
  maxPoints = c(1e+05, 5e+05),
  padWithSilence = TRUE,
  colorTheme = c("bw", "seewave", "heat.colors", "...")[1],
  col = NULL,
  extraContour = NULL,
  xlab = NULL,
  ylab = NULL,
  xaxp = NULL,
  mar = c(5.1, 4.1, 4.1, 2),
  main = NULL,
  grid = NULL,
  width = 900,
  height = 500,
  units = "px",
  res = NA,
  ...
)

Arguments

x	path to a folder, one or more wav or mp3 files c('file1.wav', 'file2.mp3'), Wave object, numeric vector, or a list of Wave objects or numeric vectors
samplingRate	sampling rate of x (only needed if x is a numeric vector)
scale	maximum possible amplitude of input used for normalization of input vector (only needed if x is a numeric vector)
from, to	if NULL (default), analyzes the whole sound, otherwise from...to (s)
winSurp	surprisal analysis window, ms (Inf = from sound onset)
input	audSpec = auditory spectrogram (audSpectrogram, speed with default settings ~= 1x), spectrogram = STFT spectrogram with (spectrogram, speed ~= 0.25x), pspec = STFT power spectrogram with (melfcc, speed ~= 0.2x), melspec = STFT mel-spectrogram with (melfcc, speed ~= 0.45x), env = analytic envelope (getRMS, speed ~= 27x) Any custom spectrogram-like matrix of features (time in columns labeled in s, features in rows) is also accepted (see examples)
takeLog	if TRUE, the input is log-transformed prior to calculating surprisal. Negative values are treated as in audSpectrogram - note that the chosen dynamic range affects this normalization (the default is 80 dB). If input = audSpec or input = spectrogram, there can be other internal preprocessing like modifying contrast based on audSpec_pars or spec_pars
audSpec_pars, spec_pars, melfcc_pars, env_pars	a list of parameters passed to audSpectrogram (if input = 'audSpec'), spectrogram (if input = 'spectrogram'), melfcc (if input = 'melspec' or 'pspec'), or getRMS (if input = 'env')
method	(for $surprisal only, has no effect on $info and $kl) acf = change in maximum autocorrelation after adding the final point; np = nonlinear prediction (see nonlinPred - works but is VERY slow); none = do not calculate $surprisal to save time and only return $info and $kl
sameLagAllFreqs	(only for method = 'acf') if TRUE, the bestLag is calculated by averaging the ACFs of all channels, and the same bestLag is used to calculate the surprisal in each frequency channel (we expect the same "rhythm" for all frequencies); if FALSE, the bestLag is calculated separately for each frequency channel (we can track different "rhythms" at different frequencies)
weightByAmpl	if TRUE, ACFs and surprisal are weighted by max amplitude per frequency channel
weightByPrecision	if TRUE, surprisal is weighted by the current autocorrelation, so deviations from a previous pattern are more surprising if this pattern is strong
onlyPeakAutocor	if TRUE, only peaks of ACFs are considered (so bestLag can never be 1, and the first change after a string of static values results in surprisal = NA)
rescale	if TRUE, surprisal is normalized from (-Inf, Inf) to [-1, 1]
summaryFun	functions used to summarize each acoustic characteristic, eg "c('mean', 'sd')"; user-defined functions are fine (see examples); NAs are omitted automatically for mean/median/sd/min/max/range/sum, otherwise take care of NAs yourself
reportEvery	when processing multiple inputs, report estimated time left every ... iterations (NULL = default, NA = don't report)
cores	number of cores for parallel processing
plot	if TRUE, plots the auditory spectrogram and the suprisalLoudness contour
savePlots	full path to the folder in which to save the plots (NULL = don't save, ” = same folder as audio)
osc	"none" = no oscillogram; "linear" = on the original scale; "dB" = in decibels
heights	a vector of length two specifying the relative height of the spectrogram and the oscillogram (including time axes labels)
ylim	frequency range to plot, kHz (defaults to 0 to Nyquist frequency). NB: still in kHz, even if yScale = bark, mel, or ERB
contrast	controls the sharpness or contrast of the image: <0 = decrease contrast, 0 = no change, >0 increase contrast. Recommended range approximately (-1, 1). The spectrogram is raised to the power of exp(3 * contrast)
brightness	makes the image lighter or darker: <0 = darker, 0 = no change, >0 = lighter, range (-1, 1). The color palette is preserved, so "brightness" works by capping an increasing proportion of image at the lightest or darkest color. To lighten or darken the palette, just change the colors instead
maxPoints	the maximum number of "pixels" in the oscillogram (if any) and spectrogram; good for quickly plotting long audio files; defaults to c(1e5, 5e5); does not affect reassigned spectrograms
padWithSilence	if TRUE, pads the sound with just enough silence to resolve the edges properly (only the original region is plotted, so the apparent duration doesn't change)
colorTheme	black and white ('bw'), as in seewave package ('seewave'), matlab-type palette ('matlab'), or any palette from palette such as 'heat.colors', 'cm.colors', etc
col	actual colors, eg rev(rainbow(100)) - see ?hcl.colors for colors in base R (overrides colorTheme)
extraContour	a vector of arbitrary length scaled in Hz (regardless of yScale, but nonlinear yScale also warps the contour) that will be plotted over the spectrogram (eg pitch contour); can also be a list with extra graphical parameters such as lwd, col, etc. (see examples)
xlab, ylab, main, mar, xaxp	graphical parameters for plotting
grid	if numeric, adds n = grid dotted lines per kHz
width, height, units, res	graphical parameters for saving plots passed to png
...	other graphical parameters

Details

Algorithm: the sound is transformed into some spectrogram-like representation (e.g., an auditory spectrogram, a mel-warped STFT spectrogram, etc.) or an RMS amplitude envelope. Using just the envelope is very fast, but then we discard all spectral information. For each frequency channel, a sliding window is analyzed to compare the actually observed final value with its expected value. There are many ways to extrapolate / predict time series and thus perform this comparison. The resulting per-channel surprisal contours are aggregated by taking their mean - optionally, weighted by the maximum amplitude of each frequency channel across the analysis window. Because increases in loudness are known to be important predictors of auditory salience, loudness per frame is also returned, as well as the product of its positive changes and surprisal.

Value

Returns a list with surprisal statistics per frame ($detailed) and per file ($summary). Calculated measures:

loudness

subjective loudness in sone, as per getLoudness

surprisal

surprisal calculated as the change in autocorrelation or as the nonlinear prediction error

surprisalLoudness

the product of surprisal and the first derivative of loudness with respect to time, treating negative values of dLoudness as zero

info

Shannon information calculated as -log(p), where p = density of Gaussian distribution at the next observation

infoW

windowed Shannon information: same as "info", but after applying a half-Gaussian taper that prioritizes more recent observations

kl

Bayesian log-surprisal: Kullback–Leibler divergence between the Gaussian distributions per frequency channel before vs. after observing the next datapoint

klW

windowed Bayesian log-surprisal

Under $detailed, the function also returns several "_mat" objects that give the same statistics with a separate value for each time-frequency bin, as well as the best lag used to calculate autocorrelation (see examples).

Examples

# A quick example
s = soundgen(nSyl = 2, sylLen = 50, pauseLen = 25, addSilence = 15)
surp = getSurprisal(s, samplingRate = 16000)
surp

## Not run: 
# A couple of more meaningful examples

## Example 1: a temporal deviant
s0 = soundgen(nSyl = 8, sylLen = 150,
              pauseLen = c(rep(200, 7), 450), pitch = c(200, 150),
              temperature = .05, plot = FALSE)
sound = c(rep(0, 4000),
          addVectors(rnorm(16000 * 3.5, 0, .02), s0, insertionPoint = 4000),
          rep(0, 4000))
spectrogram(sound, 16000, yScale = 'ERB')

# long window  (Inf = from the beginning)
surp = getSurprisal(sound, 16000, winSurp = Inf)
# Which frequency-time bins are surprising?
filled.contour(x = as.numeric(colnames(surp$detailed$surprisal_mat)) / 1000,
               y = as.numeric(rownames(surp$detailed$surprisal_mat)),
               z = t(surp$detailed$surprisal_mat),
               xlab = 'Time, s',
               ylab = 'Frequency, kHz')
hist(surp$detailed$bestLag, xlab = 'Period, s')
abline(v = .35, lty = 3, lwd = 3, col = 'blue')  # true period = 350 ms
filled.contour(x = as.numeric(colnames(surp$detailed$bestLag)) / 1000,
               y = as.numeric(rownames(surp$detailed$bestLag)),
               z = t(surp$detailed$bestLag),
               xlab = 'Time, s',
               ylab = 'Frequency, kHz')

# just use the amplitude envelope instead of an auditory spectrogram
surp = getSurprisal(sound, 16000, winSurp = Inf, input = 'env')

# increase spectral and temporal resolution (very slow)
surp = getSurprisal(sound, 16000, winSurp = 2000,
  audSpec_pars = list(nFilters = 50, step = 10,
  yScale = 'bark', bandwidth = 1/4))

# weight by increase in loudness
spectrogram(sound, 16000, extraContour = surp$detailed$surprisalLoudness /
  max(surp$detailed$surprisalLoudness, na.rm = TRUE) * 8000)

par(mfrow = c(3, 1))
plot(surp$detailed$surprisal, type = 'l', xlab = '',
  ylab = '', main = 'surprisal')
abline(h = 0, lty = 2)
plot(surp$detailed$dLoudness, type = 'l', xlab = '',
  ylab = '', main = 'd-loudness')
abline(h = 0, lty = 2)
plot(surp$detailed$surprisalLoudness, type = 'l', xlab = '',
  ylab = '', main = 'surprisal * d-loudness')
par(mfrow = c(1, 1))

# short window = amnesia (every event is equally surprising)
getSurprisal(sound, 16000, winSurp = 250)

# add bells and whistles
surp = getSurprisal(sound, samplingRate = 16000,
  yScale = 'mel',
  osc = 'dB',  # plot oscillogram in dB
  heights = c(2, 1),  # spectro/osc height ratio
  brightness = -.1,  # reduce brightness
  # colorTheme = 'heat.colors',  # pick color theme...
  col = rev(hcl.colors(30, palette = 'Viridis')),  # ...or specify the colors
  cex.lab = .75, cex.axis = .75,  # text size and other base graphics pars
  ylim = c(0, 5),  # always in kHz
  main = 'Audiogram with surprisal contour', # title
  extraContour = list(col = 'blue', lty = 2, lwd = 2)
  # + axis labels, etc
)

## Example 2: a spectral deviant
s1 = soundgen(
  nSyl = 11, sylLen = 150, invalidArgAction = 'ignore',
  formants = NULL, lipRad = 0,  # so all syls have the same envelope
  pauseLen = 90, pitch = c(1000, 750), rolloff = -20,
  pitchGlobal = c(rep(0, 5), 18, rep(0, 5)),
  temperature = .01, pitchCeiling = 7000,
  plot = TRUE, windowLength = 35)
surp = getSurprisal(s1, 16000, winSurp = 1500)
filled.contour(x = as.numeric(colnames(surp$detailed$surprisal_mat)) / 1000,
               y = as.numeric(rownames(surp$detailed$surprisal_mat)),
               z = t(surp$detailed$surprisal_mat),
               xlab = 'Time, s',
               ylab = 'Frequency, kHz')
# deviant surprising both at 1 kHz (expected tone omitted) and at the new freq
surp = getSurprisal(s1, 16000, winSurp = 1500,
  input = 'env')  # doesn't work - need spectral info

s2 = soundgen(
  nSyl = 11, sylLen = 150, invalidArgAction = 'ignore',
  formants = NULL, lipRad = 0,  # so all syls have the same envelope
  pauseLen = 90, pitch = c(200, 150),  rolloff = -20,
  pitchGlobal = c(rep(18, 5), 0, rep(18, 5)),
  temperature = .01, plot = TRUE, windowLength = 35, yScale = 'ERB')
surp = getSurprisal(s2, 16000, winSurp = 1500)

## Example 3: different rhythms in different frequency bins
s6_1 = soundgen(nSyl = 23, sylLen = 100, pauseLen = 50, pitch = 1200,
  rolloffExact = 1, invalidArgAction = 'ignore', plot = TRUE)
s6_2 = soundgen(nSyl = 10, sylLen = 250, pauseLen = 100, pitch = 400,
  rolloffExact = 1, invalidArgAction = 'ignore', plot = TRUE)
s6_3 = soundgen(nSyl = 5, sylLen = 400, pauseLen = 200, pitch = 3400,
  rolloffExact = 1, invalidArgAction = 'ignore', plot = TRUE)
s6 = addVectors(s6_1, s6_2)
s6 = addVectors(s6, s6_3)

surp = getSurprisal(s6, 16000, winSurp = Inf, sameLagAllFreqs = TRUE,
  audSpec_pars = list(nFilters = 32))
surp = getSurprisal(s6, 16000, winSurp = Inf, sameLagAllFreqs = FALSE,
  audSpec_pars = list(nFilters = 32))  # learns all 3 rhythms
filled.contour(x = as.numeric(colnames(surp$detailed$surprisal_mat)) / 1000,
               y = as.numeric(rownames(surp$detailed$surprisal_mat)),
               z = t(surp$detailed$surprisal_mat),
               xlab = 'Time, s',
               ylab = 'Frequency, kHz')

## Example 4: different time scales
s8 = soundgen(nSyl = 4, sylLen = 75, pauseLen = 50)
s8 = rep(c(s8, rep(0, 2000)), 8)
getSurprisal(s8, 16000, input = 'env', winSurp = Inf)
# ACF picks up first the fast rhythm, then after a few cycles switches to
# the slow rhythm

# Custom input: produce a nice spectrogram first, then feed it into ssm()
sp = spectrogram(s0, 16000, windowLength = 10, step = 10, contrast = .3,
  output = 'processed')  # return the modified spectrogram
colnames(sp) = as.numeric(colnames(sp)) / 1000  # convert ms to s
getSurprisal(s0, 16000, input = sp, takeLog = FALSE)

# Custom input: use acoustic features returned by analyze()
an = analyze(s0, 16000, windowLength = 20)
input_an = t(an$detailed[, 4:ncol(an$detailed)]) # or select pitch, HNR, ...
input_an = t(apply(input_an, 1, scale))  # z-transform all variables
input_an[is.na(input_an)] = 0  # get rid of NAs
colnames(input_an) = an$detailed$time / 1000  # time stamps in s
rownames(input_an) = 1:nrow(input_an)
image(t(input_an))  # not a spectrogram, just a feature matrix
getSurprisal(s0, 16000, input = input_an, takeLog = FALSE)

# analyze all sounds in a folder
surp = getSurprisal('~/Downloads/temp/', savePlots = '~/Downloads/temp/surp')
surp$summary

## End(Not run)

soundgen documentation built on Feb. 24, 2026, 5:08 p.m.