Source-Filter

Suppose you have a confusion matrix, a matrix of values with rows corresponding to stimuli and columns to responses. Suppose further that you want to analyze the similarity between your stimuli and the bias to give one response or another by applying the similarity choice model to your confusion matrix.

In the SCM, the probability of giving response r to stimulus s is (where is the bias to give response and is the similarity between stimulus and stimulus , and is the number of stimuli):

(1)  

You might want to use this model because it has convenient, closed-form solutions for the similarity and bias parameters:

(2)  

Whence the values? You could just divide your confusion counts by the row totals, or you could, since you care about statistical rigor, first obtain maximum likelihood estimates of your confusion counts, only then normalizing by the row sums to obtain MLEs of the confusion probabilities.

And you can obtain MLEs of your confusion counts by using iterative proportional fitting. I learned about iterative proportional fitting from Rob Nosofsky in one of the best classes I took as a graduate student (or any other kind of student, come to think of it). And just today, I found myself needing to use iterative proportional fitting once again.

So I wrote an R function:

iter.prop.fit <- function(M,delta = .001){

  nr <- nrow(M)
  nc <- ncol(M)
  M.h <- matrix(1,nrow=nr,ncol=nc)
  dimnames(M.h) <- dimnames(M)
  d.t <- 1
  nz <- .001
  while(d.t >= delta){
    M.h.a <- M.h
    for(ri in 1:nr){
      for(ci in 1:nc){
        M.h[ri,ci] <- M.h[ri,ci]*sum(M[ri,])/sum(M.h[ri,]+nz)
        M.h[ri,ci] <- M.h[ri,ci]*sum(M[,ci])/sum(M.h[,ci]+nz)
        M.h[ri,ci] <- M.h[ri,ci]*(M[ri,ci]+M[ci,ri])/(M.h[ri,ci]+M.h[ci,ri]+nz)
      }
    }
    d.t <- max(abs(M.h.a-M.h),na.rm=T)
  }
  return(M.h)

}

M is your confusion matrix, delta is a criterion that determines when the algorithm stops (smaller delta means less change from one step of the algorithm to the next), and M.h (h = hat, so M.h = ) is the matrix of maximum likelihood confusion counts.

The basic idea is to take a matrix full of ones and gradually adjust it to match, as closely as possible, the counts in the observed confusion matrix.

Because we’re using iterative proportion fitting to get SCM MLEs, there are three steps for each iteration of the algorithm. First, M.h is adjusted by row (i.e., stimuli). Second, M.h is adjusted by column (i.e., responses). Third, and finally, M.h is adjusted for each unique pair of stimuli (i.e., symmetric similarity).

When the largest change from one step to the next is smaller in magnitude than delta, the algorithm stops and returns M.h. Once you have M.h, you can normalize by row to get maximum likelihood confusion probability estimates (i.e., values) which you can plug into the equations above to obtain maximum likelihood estimates of the similarity and bias parameters.

What you do with the these is, of course, up to you.

I find this match of algorithm and model interesting, primarily because it seems kind of backwards. When you want to find MLEs for the parameters of, say, GRT, you can run a Newton-Raphson algorithm that navigates parameter space to find the MLEs for a given data set.

This makes intuitive sense to me – you have data, you make an educated guess about initial values, and you have rules for moving from your initial guess to your best guess.

With iterative proportional fitting and the SCM, on the other hand, you find MLEs of the data first, and only after you have this do you calculate parameter values. You don’t make an initial guess about your parameters, and you don’t explicitly navigate parameter space in order to find your best guess. Of course, because you have closed form solutions for the SCM parameters, for each step in the algorithm, you could calculate parameter values if you wanted to, so, even though iterative proportional fitting is expressed entirely in terms of the data, the algorithm is effectively following a simple set of rules for implicitly navigating parameter space.

I’m at a loss to understand why it should be problematic to calculate false alarm rates with a Bayesian statistical model while it is acceptable to calculate statistical power with that same model.

In null-hypothesis significance testing (NHST), the false alarm rate (or, if you’re a glutton for confusing terminology, the Type I error rate) is the rate at which a true null model is incorrectly rejected. Power, on the other hand, is the probability of correctly rejecting a false null model (or 1 – the Type II error rate).

Even if we’re doing a very simple test (e.g., comparing the means of two samples), the “null” model can, in principle, be pretty much whatever you want. Of course, there are conventional choices for various tests (e.g., for comparing two means, ). In a Bayesian framework, the model is more elaborate, since you have priors, at least, and maybe also a “region of practical equivalence,” or ROPE, which is an interval rather than a point null (see, e.g., Kruschke’s BEST).

In both NHST and Bayesian frameworks, you also have a decision procedure. In NHST, for example, you would reject a null hypothesis if an observed test statistic exceeds a certain criterion. The criterion is supposed to be set, a priori, based on some determination of acceptably low false alarm rates (α). In a Bayesian framework, the decision procedure might be to reject a model (or a particular assumption about a model) if, say, the 95% highest density interval (HDI) for the posterior distribution of a parameter of interest is entirely outside the ROPE associated with that model (assumption).

Okay, so, for calculating power, whether you’re working as a Bayesian or not, you also have an alternative (to the null) model. Using any of a number of methods (e.g., analytic math, simulations), you figure out how likely you are to reject a particular null model using your decision procedure and assuming that your alternative model is true. Fairly straightforward, and pretty uncontroversial (at least if it’s prospective power).

For calculating false alarm rates, you do something very much analogous to calculating power. In NHST, you assume that your null model is true, and then you take your decision procedure and calculate how likely you are to incorrectly reject that null model. For a given false alarm rate (α), you can figure out associated criterion value(s) of your test statistic, and if your observed test statistic exceeds the criterion, you reject the null.

In a Bayesian framework, you can do more or less the same thing, I would think. You assume your (more elaborate) null(ish) model is true – priors, ROPE, and all – and you use your decision procedure to estimate the proportion of tests that would incorrectly reject the null(ish) model.

For Kruschke, the crucial bit that makes power okay is that, unlike false alarm rates, “power is not being used to make decisions about data that are actually obtained, nor is power being used to set decision criteria for data that are actually obtained.”

So, here’s what I don’t get. If you’ve done the work to carry out a Bayesian t-test, then you’ve already got everything you need to calculate a false alarm rate: a null(ish) model with priors, a ROPE, and a HDI-based decision procedure. Which is to say that for every model sufficient to the task of testing a statistical hypothesis with a particular data set, there is a false alarm rate. It may not be the false alarm rate for Bayesian t-test’s in general, but it is the false alarm rate for the model and decision procedure you’ve chosen to use.

It seems to me that this false alarm rate provides useful information. If your model has a 0.50 false alarm rate, for example, you might want to know that, figure out why it’s so high, and then change your model to lower that number. Which is to say, in part, that it is appropriate to take the false alarm rate for a given model and decision procedure into account when making decisions about obtained data.

For a given model (or parameter therein) for which you want to make a decision, the utility of knowing the relevant false alarm rate is much like the utility of knowing the relevant power. Whether you’re a Bayesian or a devotee of NHST, if you’re carrying out statistical tests, power and false alarm rates are both relevant.

Full disclosure: I took multiple statistics and cognitive modeling courses from John Krushcke, and I was the grader for his Bayesian Data Analysis class multiple times as a doctoral student at IU. When I was his grader, he was developing his BDA textbook, and I’m even acknowledged in the book (next to one Thomas Wisdom, who I don’t know, but who has a much better last name than me). So maybe this blog post is kind of like me becoming Darth Vader to John’s Obi Wan Kenobi. Or maybe that’s overstating the magnitude of this disagreement.