Bayesian Analysis for Psychometric Functions

Yasuharu Okamoto, 2017.07

2AFC tasks with rating by two, three, or four categories

That three- or four-category rating method is preferable to the two-category method in constant stimuli was shown by simulation of Okamoto (2017).

An observer is usually asked to compare two stimuli, the standard stimulus and a comparison stimulus, and to respond by choosing one of the two response categories, e.g. “The left stimulus is brighter than the right one” and “The right stimulus is brighter than the left one”. However, in some cases, the observer feels difficulty to decide. Klein (2001) proposes to use more than two response categories. This website provides Stan scripts for 2AFC task with more than two rating categories, in addition to a Stan script for usual 2AFC task with binary judgment. First, discuss the usual 2AFC with binary response categories, then the one with three rating categories, and the one with four rating categories.

Two Alternative Forced Choice Task with Two Response Categories

A comparison stimulus is compared with the standard stimulus , and judgment “Which is stronger?” is required. The probability that is judged to be stronger than is denoted by . In a two-alternative procedure, an observer usually chooses one of two judgments, ( is stronger than ) and ( is weaker than ), so the probability is given by

can be considered to be a function of and is called a psychometric function (PF). A psychometric function is well approximated by a cumulative normal distribution (phi-gamma hypothesis; see Luce and Galanter, 1963). Hence, set

where is the cumulative standard normal distribution.

The point of subjective equality PSE is given by

Just noticeable difference JND is given by

where discrimination probability 0.75 is chosen for JND. satisfies the equation

Denote the comparison stimulus and its judgment in trial by and , where

Then the probability of the data is given by

where

Based on probability (1.2), the posterior distributions are given.

A Stan script for the above model is this.

================ pf2cat.stan ===============================

data {

int N;

real St[N];

int<lower = 1, upper = 2> Res[N];

}

parameters {

real mu;

real<lower = 0.0> sigma;

}

transformed parameters {

vector[2] p[N];

for (i in 1:N){

p[i][2] = Phi((St[i] - mu) / sigma);

p[i][1] = 1.0 - p[i][2];

}

model {

for (i in 1:N)

Res[i] ~ categorical(p[i]);

}

==========================================================

The above Stan script is used by the following Python script. The scripts files and a sample data file are archived into the file PF2Cat.zip, which can be down loaded by clicking the file name PF2Cat.zip.

================= main.py ================================

from subprg import *

import pystan

from pystan import StanModel

import pickle

X, n, meanSt, sgmSt, fout, nm_fout = InputData()

N = n

St = []

Res = []

for v in X:

St.append(v[0])

Res.append(v[1])

print('\nSt = \n', St)

print('\nRes = \n', Res)

Data = {'N': N, 'St': St, 'Res': Res}

minSt = St[0]

maxSt = St[0]

for i in range(N):

if minSt > St[i]:

minst = St[i]

if maxSt < St[i]:

maxSt = St[i]

def f_init():

return dict(mu = (minSt + maxSt) / 2.0, sigma = (maxSt - minSt) / 2.0)

fit = pystan.stan(file = 'pf2cat.stan', data = Data, seed = 999, chains =4,

iter = 5000, warmup = 2000,

init = f_init, pars = ['mu', 'sigma'], n_jobs = 1)

print(fit)

mrg_mu = fit['mu']

mrg_sgm = fit['sigma']

mrg_mu.sort()

mrg_sgm.sort()

NSimu = len(mrg_mu)

print('\nNSimu = ', NSimu)

DispResults( NSimu, X, n, mrg_mu, mrg_sgm, fout)

fout.close()

print("\nOutput file {} was saved.".format(nm_fout))

=======================================================

Run the program main.py, then the name of the input data file is asked (Figure 1.1).

Figure 1.1

The format of the data file is shown in Figure 1.2.

Figure 1.2

The data are written after the line with slash / at the head. Each line of the data includes a set of dataID, (denoted by “St”), and in each trial. DataID is simply any string and used to identify the individual data, but is not used in analysis. After the last data, a line with slash / at the head is put to indicate the end of the data.

After the input data file name is set and the Enter key is pressed down, then an output text file name is asked. Set a text file name, then calculation starts.

When calculation by Stan ends, a window, in which the PF is drawn, is displayed (Figure 1.3).

Figure 1.3

Small green circles represent data points s. To avoid overlapping of circles, small random fluctuations are added to get the drawn points s, where

is a small random fluctuation.

After closing the window in Figure 1.3, the next window, in which a histogram of samples for PSE () is drawn, will be displayed (Figure 1.4).

Figure 1.4

After the window in Figure 1.4 is closed, the window, in which a histogram for JND is drawn, is displayed (Figure 1.5).

Figure 1.5

Close the window in Figure 1.5, then the program ends (Figure 1.6).

Figure 1.6

After the program ends, the output file can be opened by an editor. Figure 1.7 shows the contents of the output file.

・

Figure 1.7

At the end of the file, medians, quartile points, and 95% CIs are written.

Two Alternative Forced Choice Task with Three Rating Categories

A task with three rating categories allows “don’t know” or “indifferent” judgment as an additional one to the binary decision judgment (Böckenholt, 2001; Kaernbach, 2001). Hence, the three categories of judgments are as follows

: which means that is judged to be weaker than

: which means that is judged to equal

: which means that is judged to be stronger than

Combine the two categories and into one category, which is denoted by .

Set the following two PFs

Then we have

and

JND is given by

Assuming that C1 and C2 are at an equal distance from PSE, we have

Denote the comparison stimulus and its judgment in trial by and , where

Then the probability of the data is given by

where

A Stan script for model (2.2) is this.

================= pf3cat.stan ============================

data {

int N;

real St[N];

int<lower = 1, upper = 3> Res[N];

}

parameters {

ordered[2] C;

real<lower = 0.0> sigma;

}

transformed parameters {

vector[3] p[N];

real pse;

pse = (C[1] + C[2]) / 2.0;

for (i in 1:N){

p[i][3] = Phi((St[i] - C[2]) / sigma);

p[i][2] = Phi((St[i] - C[1]) / sigma) - Phi((St[i] - C[2]) / sigma);

p[i][1] = 1.0 -Phi((St[i] - C[1]) / sigma);

}

model {

for (i in 1:N)

Res[i] ~ categorical(p[i]);

}

==================================================

The above Stan script was used by the following Python script. The script files were archived into the file PF3Cat.zip, which can be downloaded by clicking the file name PF3Cat.zip.

================= main.py ========================

import matplotlib.pyplot as plt

import numpy as np

import math

from RNGen import *

from subprg import *

import pystan

X, n, meanSt, sgmSt, fout, nm_fout = InputData()

N = n

St = []

Res = []

for v in X:

St.append(v[0])

Res.append(v[1])

print('\nSt = ', St)

print('\nRes = ', Res)

Data = {'N': N, 'St': St, 'Res': Res}

minSt = St[0]

maxSt = St[0]

for v in St:

if minSt > v:

minSt = v

if maxSt < v:

maxSt = v

def i_func():

return dict(C = [(3.0 * minSt + 2.0 * maxSt) / 5.0,

(2.0 * minSt + 3.0 * maxSt) / 5.0],

sigma = (maxSt - minSt) / 2.0)

fit = pystan.stan(file = 'pf3cat.stan', data = Data, seed = 999,

iter = 5000, warmup = 2000, pars = ['C', 'pse', 'sigma'],

init = i_func, n_jobs = 1)

print(fit)

mrg_sgm = fit['sigma']

mrg_C = fit['C']

mrg_C1 = []

mrg_C2 = []

for v in mrg_C:

mrg_C1.append(v[0])

mrg_C2.append(v[1])

mrg_pse = fit['pse']

mrg_C1.sort()

mrg_C2.sort()

mrg_pse.sort()

mrg_sgm.sort()

NSimu = len(mrg_sgm)

n = N

DispResults( NSimu, X, n, mrg_C1, mrg_C2, mrg_pse, mrg_sgm, fout )

fout.close()

print("\nOutput file {} was saved.".format(nm_fout))

====================================================

Run this Python script, then the name of the input data file is asked (Figure 2.1).

Figure 2.1

The format of the data file is shown in Figure 2.2.

・

Figure 2.2

The data are written after the line with slash / at the head. Each line of the data includes a set of dataID, S_comp (denoted by “St”), and Res in each trial. DataID is simply any string and used to identify the set, but is not used in analysis. After the last data, a line with slash / at the head is put to indicate the end of the data.

After the input data file name is set and the Enter key is pressed down, then an output text file name is asked. Set a text file name, then calculation starts.

When calculation by Stan ends, a window, in which the three PFs are drawn, is displayed (Figure 2.3).

Figure 2.3

Two PFs denoted by “PF1” and “PF2” are corresponding functions to C1 and C2. The PF denoted by “pse” is the one, which is estimated by the value PSE given by Assumption (2.1).

Small green circles represent data points s. To avoid overlapping of circles, small random fluctuations are added to get the drawn points s, where

is a small random fluctuation.

After closing the window in Figure 2.3, the next window, in which frequency polygons of samples for C1, pse and C2 are drawn, will be displayed (Figure 2.4).

Figure 2.4

After the window in Figure 2.4 is closed, the window, in which a histogram for JND is drawn, is displayed (Figure 2.5).

Figure 2.5

Close the window in Figure 2.5, then the program ends.

After the program ends, the output file can be opened by an editor. Figure 2.6 shows the contents of the output file.

・

Figure 2.6

At the end of the file, medians, quartile points, and 95% CIs for C1, pse, C2 and JND are written.

Two Alternative Forced Choice Task with Four Rating Categories

A 2AFC task with four rating categories allows four response categories:

: which means that is judged to be weaker than

: which means that is probably weaker than

: which means that is probably stronger than

: which means that is judged to be stronger than

Combine the three response categories , , and into one category, which is denoted by .

Combine the two categories and into one category, which is denoted by , and the two categories and into one category, which is denoted by .

Then, set the following three PFs

We have the following probabilities

PSE is given by

and JND is given by

Denote the comparison stimulus and its judgment in trial by and , where

Then the probability of the data is given by

where

A Stan script for model (3.1) is this.

=================== pf4cat.stan ==========================

data {

int N;

real St[N];

int<lower = 1, upper = 4> Res[N];

}

parameters {

ordered[3] C;

real<lower = 0.0> sigma;

}

transformed parameters {

vector[4] p[N];

for (i in 1:N){

p[i][4] = Phi((St[i] - C[3]) / sigma);

p[i][3] = Phi((St[i] - C[2]) / sigma) - Phi((St[i] - C[3]) / sigma);

p[i][2] = Phi((St[i] - C[1]) / sigma) - Phi((St[i] - C[2]) / sigma);

p[i][1] = 1.0 -Phi((St[i] - C[1]) / sigma);

}

model {

for (i in 1:N)

Res[i] ~ categorical(p[i]);

}

=====================================================

A Python script using the above Stan script is as follows. The script files were archived into the file PF4Cat.zip, which can be down loaded by clicking the file name PF4Cat.zip.

=====================================================

from subprgs4Cat import *

import pystan

X, n, meanSt, sgmSt, fout, nm_fout = InputData()

N = n

St = []

Res = []

for v in X:

St.append(v[0])

Res.append(v[1])

print('\nSt = \n', St)

print('\nRes = \n', Res)

Data = {'N': N, 'St': St, 'Res': Res}

minSt = St[0]

maxSt = St[0]

for v in St:

if minSt > v:

minSt = v

if maxSt < v:

maxSt = v

def i_func():

return dict(C = [(3.0 * minSt + 2.0 * maxSt) / 5.0,

(minSt + maxSt) / 2.0,

(2.0 * minSt + 3.0 * maxSt) / 5.0],

sigma = (maxSt - minSt) / 2.0)

fit = pystan.stan(file = 'pf4cat.stan', data = Data, init = i_func,

iter = 5000, warmup = 2000, pars = ['C', 'sigma'], n_jobs = 1)

print(fit)

mrg_sgm = fit['sigma']

mrg_C = fit['C']

mrg_C1 = []

mrg_Mu = []

mrg_C2 = []

for v in mrg_C:

mrg_C1.append(v[0])

mrg_Mu.append(v[1])

mrg_C2.append(v[2])

mrg_sgm.sort()

mrg_C1.sort()

mrg_Mu.sort()

mrg_C2.sort()

NSimu = len(mrg_sgm)

n = N

DispResults( NSimu, X, n, mrg_C1, mrg_Mu, mrg_C2, mrg_sgm, fout )

fout.close()

print("\nOutput file {} was saved.".format(nm_fout))

====================================================

Run the above Python script, then the name of the input data file is asked (Figure 3.1).

Figure 3.1

The format of the data file is shown in Figure 3.2.

・

Figure 3.2

After the input data file name is set and the Enter key is pressed down, then an output text file name is asked. Set a text file name, then calculation starts.

When calculation ends, a window, in which the three PFs are drawn, is displayed (Figure 3.3).

Figure 3.3

Small green circles represent data points s. To avoid overlapping of circles, small random fluctuations are added to get the drawn points s, where

is a small random fluctuation.

After closing the window in Figure 3.3, the next window, in which frequency polygons of samples for C1, PSE, and C2 are drawn, will be displayed (Figure 3.4).

Figure 3.4

After the window in Figure 3.4 is closed, the window, in which a histogram for JND is drawn, is displayed (Figure 3.5).

Figure 3.5

Close the window in Figure 3.5, then the program ends.

After the program ends, the output file can be opened by an editor. Figure 3.6 shows the contents of the output file.

Figure 3.6

At the end of the file, medians, quartile points, and 95% CIs for C1, PSE, C2 and JND are written.

Reference

Böckenholt,U. (2001). Threshold and intransitivities in pairwise judgments: A multilevel analysis. Journal of Educational and Behavioral Sciences, 26, 269-282.

Kaernbach, C. (2001). Adaptive threshold estimation with unforced-choice tasks. Perception & Psychophysics, 63, 1377-1388.

Klein, S. A. (2001). Measuring, estimating, and understanding the psychometric function: A commentary. Perception & Psychophysics, 63, 1421-1455.

Luce, R. D. & Galanter, E. (1963). Discrimination. In P. Suppes, J. L. Zinnes, R. R. Bush, E. Galanter, R.D. Luce, W. J. McGill, A. Newell, and H. A. Simon (Eds.), Handbook of Mathematical Psychology, Vol. I (pp. 191-243). New York: John Wiley and Sons, Inc.

Back