In kleinschmidt/beliefupdatr: Belief updating for phonetic adaptation

This package includes a Stan model for inferring listeners prior beliefs, based on adaptation behavior.

library(tidyverse)
library(magrittr)
library(beliefupdatr)

## devtools::install_github('kleinschmidt/phonetic-sup-unsup')
library(supunsup)

library(rstan)

set.seed(1L)

Load and format data

d <- supunsup::supunsup_clean %>%
  filter(supCond == 'unsupervised') %>%
  select(vot, trueCat, respCat, trial, bvotCond, subject)

Fit the model

Getting ready

You can list the stan models included with

beliefupdatr::list_models()

Stan models are pre-compiled when the package is installed. If you want to access them directly, the compiled stanmodel objects are stored in the list:

mod <- beliefupdatr:::stanmodels[['conj_id_lapsing_fit']]

These include the original Stan code, the translated C++ code, and the compiled objects.

A few tips for speeding things up: if you have multiple cores and enough memory, you can run the chains in parallel by setting

options(mc.cores = parallel::detectCores())

Then Stan will use however many cores you have, and can run up to that many chains in parallel. Note, though, that if you only have four cores this can tie up your laptop pretty well, so I typically use one less than the number of cores I really have (e.g., parallel::detectCores()-1).

Actually fitting the model

This package provides the function infer_prior_beliefs that wraps the process of converting data from a data frame into stan format, and running the stan model sampler. It requires training data set (from which the sufficient stats of each exposure category are determined), a test data set (behavioral responses), and the names of the columns for the cue value, true category, subjects' response, and condition variable (group).

fit <- infer_prior_beliefs(df = d,
                           cue = 'vot',
                           category = 'trueCat',
                           response = 'respCat',
                           condition = 'bvotCond',
                                        # use first subject in each condition
                                        # to determine exposure statistics:
                           ranefs = 'subject',
                           chains = 4, iter = 500)

If all goes well, that shouldn't take too long. Once you've fit the model, you can eyeball the inferred prior values with

print(fit, digits=1, pars=c('kappa_0', 'nu_0', 'mu_0', 'sigma_0', 'lapse_rate'))

The first thing to look for here is the Rhat statistic. This is a measure of how well your chains have mixed, which is a way of approximating how successful your sampler was. The idea is that if your chains are sampling from the actual posterior, then they'll mix together well. Thus, if they don't mix well, you know that something is wrong. Of course, just because they mix well doesn't necessarily mean everything is hunky-dory, but it's a very important check.

The next thing to look at is the parameter values themselves. Are the of the right orders of magnitude? Are the values actually reasonable? Can you interpret them?

Inferring prior with incremental updating

The basic model above assumes that all training data has been encountered before any of the test data, and it ignores any incrementality in how listeners' behavior changes as they update their beliefs. Can we capture belief updating behavior under the assumption that listeners incrementally update their beliefs, and what to the inferred prior beliefs look like?

We have to format the data in a slightly different way here. First, instead of representing the training data in terms of raw observations, we're representing it in terms of sufficient statistics: count, mean, and standard deviation. The Stan code for the model above converts the raw observations into this format, but we gain more flexibility by doing it ourselves outside of Stan in this case.

Second, we're making a simplifying assumption that of the sufficient statistics, only the count changes over blocks. This is reasonable when modeling aggregate data because each subject gets a random trial order, so on average the statistics of each category will mirror the true statistics within each block. Specifically, we model each blocks' test observations based on the number of tokens from each category a subject has seen, on average, up to halfway through that block. This is a compromise between assuming that no updating happens throughout the block (as far as the test trials are concerned), and that the responses reflect the updated beliefs at the end of the block.

This is all handled by the same wrapper function as before. If you pass a value for the optional n_blocks argument, the incremental model will be used, assuming that many blocks:

fit_inc <-
  infer_prior_beliefs(df = d,
                      cue = 'vot',
                      category = 'trueCat',
                      response = 'respCat',
                      condition = 'bvotCond',
                                        # use first subject in each condition
                                        # to determine exposure statistics:
                      ranefs = 'subject',
                      n_blocks = 10,    # discretize trials into 10 blocks
                      iter = 500, chains = 4, init = 0) %T>%
  print(digits=1, pars=c('kappa_0', 'nu_0', 'mu_0', 'sigma_0', 'lapse_rate'))