## Assumed distributions of signal and noise

$$X_n \sim N(0,1) =f_n(x | \theta=0)= \frac{1}{\sqrt{2\pi}}e^{-\frac{x^2}{2}}$$

$$X_s \sim N(d',1) =f_s(x | \theta=d')= \frac{1}{\sqrt{2\pi}\sigma}e^{-\frac{(x-d')^2}{2 \sigma^2}}$$

### Example

library(tidyverse)
dprime <- 1
sigma <- 2
x <-  seq(-5, 5, .1)

dist_dat <- tibble(x) %>%
mutate(f_s = dnorm(x, dprime, sigma),
f_n = dnorm(x, 0, 1))

p <- dist_dat %>%
ggplot() +
geom_line(aes(x = x, y = f_s, lty = "signal")) +
geom_line(aes(x = x, y = f_n, lty = "noise")) +
labs(y = "f")
p

## Likelihood

$\log{\Lambda}(x) = -\frac{1}{2\sigma^2} \left( (1-\sigma^2) x^2 - 2d'x+d'^2+2\sigma^2\log(\sigma) \right)$

### Example

ll_fun <- function(x, p) - 1 / (2 * p[2]^2) * ( (1 - p[2]^2) * x^2 - 2 * p[1] * x + p[1]^2 +2*p[2]^2*log(p[2]) )

ll <- tibble(x) %>%
mutate(ll = ll_fun(x, c(dprime, sigma)))

ll %>%
ggplot() +
geom_line(aes(x = x, y = ll)) 

It can be seen that the log likelihood is non-monotonic on x. For large positive and negative values of x, it is more likely that the state of the world was signal, rather than noise.

## Decision rule

For $$\log{\Lambda}(x) > \lambda$$ choose $$a_S$$

results in

$-\frac{1}{2\sigma^2} \left( (1-\sigma^2) x^2 - 2d'x+d'^2+2\sigma^2\log(\sigma) + 2\sigma^2 \lambda \right) > 0$ $(1-\sigma^2) x^2 - 2d'x+d'^2+2\sigma^2 (\log(\sigma) + \lambda) < 0$ When $$\sigma > 1$$, the decision rule is

$\text{If } x < c_1 \, \text{ or } \, x > c_2 \text{ then }a=a_S$ $\text{If } c_1 < x < c_2 \text{ then }a=a_N$

When $$\sigma < 1$$, the decision rule is

$\text{If } x < c_1 \, \text{ or } \, x > c_2 \text{ then }a=a_N$

$\text{If } c_1 < x < c_2 \text{ then }a=a_S$

where $$c_1$$ and $$c_2$$ are given by

$c_1 = \frac{d'+\sqrt{d'^2-(1-\sigma^2)(d'^2+2\sigma^2(log(\sigma) + \lambda))}}{1-\sigma^2}$ $c_2 = \frac{d'-\sqrt{d'^2-(1-\sigma^2)(d'^2+2\sigma^2(log(\sigma) + \lambda))}}{1-\sigma^2}$ We can use the sensory space to decide, but there are 2 relevant criteria. The observer model is completely specified by $$d'$$, $$\sigma$$ and $$\lambda$$. For the model to be identifiable, we need that the observer performs the task using at least two criteria.

Notice that for having real solutions $$\lambda$$ should be larger that the minimun of the log likelihood.

### Example

For the proposed example, let’s calculate $$c_1$$ and $$c_2$$ in the sensory space for different likelihood criteria including the ideal observer ($$\lambda=0$$)

calculate_c1_c2 <- function(.dprime, .sigma, .lambda) {
sqr_root <- sqrt(.dprime^2 - (1-.sigma^2)*(.dprime^2 + 2*.sigma^2*(log(.sigma) + .lambda)))

c1 <- (.dprime + sqr_root) / (1 - .sigma^2)
c2 <- (.dprime - sqr_root) / (1 - .sigma^2)

tibble(c1, c2)
}

c1_c2 <- tibble(lambda = c(-0.8595, -0.8, -0.5 , 0, 1, 3)) %>%
group_by(lambda) %>%
summarise(calculate_c1_c2(dprime, sigma, lambda))
## summarise() ungrouping output (override with .groups argument)
p +
geom_vline(data = c1_c2, aes(xintercept = c1, color = factor(lambda))) +
geom_vline(data = c1_c2, aes(xintercept = c2, color = factor(lambda)))

For the ideal obsercer, $$c_1$$ and $$c_2$$ corresponds to the intersections of the two distributions.

You can compute $$c_1$$ and $$c_2$$ also using the function UVGSDTcrit from the book Modelling Psychophysical data in R

UVGSDTcrit <- function(dp, sig, logbeta){
TwoSigSq <- 2 * sig^2
minLam <- optimize(function(X, dp, sig)
-((1 - sig^2) * X^2 - 2 * dp * X +
dp^2 + TwoSigSq * log(sig))/TwoSigSq, c(-10, 10),
dp = dp, sig = sig)\$objective
if (logbeta < minLam) warning("complex roots")
cf <- -c(dp^2 + TwoSigSq * log(sig) + logbeta *
TwoSigSq,
-2 * dp, 1 - sig^2)/TwoSigSq
proot <- polyroot(cf)
sort(Re(proot))
}

For example, for the optimal criterion

UVGSDTcrit(dprime, sigma, 0)
## [1] -1.847545  1.180878

## Probabilities

Under the 0-1 loss model, the observer model can be reparametrized in terms of the probabilities. Let’s consider the situation $$\sigma >1$$

$$p_{FA}=P(a_S|H_N)=\int_{x \in R} f_n(x) \,dx = 1 - \left(\Phi(c_2)- \Phi(c_1)\right)$$

$$p_{H}=P(a_S|H_N)=\int_{x \in R} f_s(x) \,dx = 1 - \left( \Phi\left(\frac{c_2-d'}{\sigma}\right) - \Phi \left(\frac{c_1-d'}{\sigma}\right) \right)$$

## ROC

Let’s plot the ROC curve for the previous example

c1_c2_roc <- c1_c2 %>%
mutate(p_fa = 1 - (pnorm(c2) - pnorm(c1)),
p_h = 1 - (pnorm(c2, dprime, sigma) - pnorm(c1, dprime, sigma)))

p_c1_c2_roc <- c1_c2_roc %>%
ggplot() +
geom_line(aes(x = p_fa, y = p_h)) +
geom_point(aes(x = p_fa, y = p_h, color = factor(lambda))) +
geom_segment(color = "grey", aes(x = 0, y = 0, xend = 1, yend = 1)) +
xlim(0, 1) +
ylim(0, 1)

p_c1_c2_roc

As expected the criterion $$\lambda = 0$$ maximizes correct responses

c1_c2_roc %>%
mutate(p_m = 1 - p_h,
p_cr = 1 - p_fa,
p_correct = .5 *(p_h + p_cr))
## # A tibble: 6 x 8
##   lambda     c1      c2    p_fa   p_h     p_m   p_cr p_correct
## *  <dbl>  <dbl>   <dbl>   <dbl> <dbl>   <dbl>  <dbl>     <dbl>
## 1 -0.860 -0.362 -0.304  0.978   0.991 0.00924 0.0218     0.506
## 2 -0.8   -0.733  0.0660 0.706   0.873 0.127   0.294      0.584
## 3 -0.5   -1.31   0.646  0.354   0.694 0.306   0.646      0.670
## 4  0     -1.85   1.18   0.151   0.541 0.459   0.849      0.695
## 5  1     -2.56   1.89   0.0344  0.365 0.635   0.966      0.665
## 6  3     -3.54   2.87   0.00222 0.186 0.814   0.998      0.592

### ROC using a sensory criterion

Let’s suppose that the observer, instead of using a likelihood criterion uses a criterion in the sensory evidence

$$p_{FA} = 1 - \Phi(c)$$

$$p_{H} = 1 - \Phi\left(\frac{c-d'}{\sigma}\right)$$

c1_c2_roc_sensory_crit <- tibble(crit = seq(-10, 10, .1)) %>%
mutate(p_fa = 1 - pnorm(crit),
p_h = 1 - pnorm(crit, dprime, sigma))

p_c1_c2_roc  +
geom_line(data = c1_c2_roc_sensory_crit, aes(x = p_fa, y = p_h, lty = "sensory criterion")) +
scale_linetype_manual(values = 2)

The model predicts that for large proportion of false alarms, this prorportion is larger than proportion of hits.

Notice that if we use a more moderate value for $$\sigma$$ the ROC curves for the likelihood and sensory criterion are very similar

c1_c2_roc_smaller_sigma <- tibble(lambda = c(-0.8, -0.8, -0.5 , 0, 1, 3)) %>%
group_by(lambda) %>%
summarise(calculate_c1_c2(dprime, 1.5, lambda)) %>% # sigma = 1.5
mutate(p_fa = 1 - (pnorm(c2) - pnorm(c1)),
p_h = 1 - (pnorm(c2, dprime, sigma) - pnorm(c1, dprime, sigma)))
## summarise() regrouping output by 'lambda' (override with .groups argument)
c1_c2_roc_sensory_crit_smaller_sigma <- tibble(crit = seq(-10, 10, .1)) %>%
mutate(p_fa = 1 - pnorm(crit),
p_h = 1 - pnorm(crit, dprime, 1.5)) # sigma = 1.5

c1_c2_roc_smaller_sigma %>%
ggplot() +
geom_line(aes(x = p_fa, y = p_h)) +
geom_point(aes(x = p_fa, y = p_h, color = factor(lambda))) +
geom_segment(color = "grey", aes(x = 0, y = 0, xend = 1, yend = 1)) +
geom_line(data = c1_c2_roc_sensory_crit_smaller_sigma, aes(x = p_fa, y = p_h, lty = "sensory criterion")) +
scale_linetype_manual(values = 2) +
xlim(0, 1) +
ylim(0, 1)  

### Slope in ROC and likelihood at criterion

The slope in ROC is the likelihood at criterion independently of the decision rule used

slope_loglik_crit <- tibble(lambda = seq(-0.8595, 3, 0.005)) %>%
group_by(lambda) %>%
summarise(calculate_c1_c2(dprime, sigma, lambda)) %>%
mutate(p_fa = 1 - (pnorm(c2) - pnorm(c1)),
p_h = 1 - (pnorm(c2, dprime, sigma) - pnorm(c1, dprime, sigma)),
lik = exp(ll_fun(c1, c(dprime, sigma))),
slope = ( lag(p_h) - p_h) / (lag(p_fa) - p_fa)) %>%
select(p_fa, p_h, lik, slope) %>%
mutate(cond = "loglik crit")
## summarise() ungrouping output (override with .groups argument)
slope_sensory_crit <- tibble(crit = seq(-3, 3, .05)) %>%
mutate(p_fa = 1 - pnorm(crit),
p_h = 1 - pnorm(crit, dprime, sigma),
lik = exp(ll_fun(crit, c(dprime, sigma))),
slope = ( lag(p_h) - p_h) / (lag(p_fa) - p_fa)) %>%
select(p_fa, p_h, lik, slope) %>%
mutate(cond = "sensory crit")

slope_loglik_crit %>%
bind_rows(slope_sensory_crit) %>%
ggplot()+
geom_segment(color = "grey", aes(x = 0, y = 0, xend = 1, yend = 1)) +
geom_line(aes(x = p_fa, y =p_h, color = cond, size = "ROC")) +
geom_line(aes(x = p_fa, y = lik, color = cond, lty = "likelihood at criterion", size = "SLOPES")) +
geom_line(aes(x = p_fa, y = slope, color = cond, lty = "slope", size = "SLOPES")) +
scale_size_manual(values = c(2, .5)) +
coord_cartesian(xlim = c(0, 1), y = c(0, 1))
## Warning: Removed 2 row(s) containing missing values (geom_path).

### zROC

The zROC curve is a straight line only for the decision rule using the sensory criterion.

$$p_{FA} = 1 - \Phi(c) = \Phi(-c)$$

$$p_{H} = 1 - \Phi(\frac{c-d'}{\sigma}) = \Phi(\frac{d'-c}{\sigma})$$

then

$$\Phi^{-1}(p_{FA}) = -c$$

$$\Phi^{-1}(p_{H}) = \frac{d'-c}{\sigma}$$

That is

$$\Phi^{-1}(p_{H}) = \frac{1}{\sigma}d' + \frac{1}{\sigma}\Phi^{-1}(p_{FA})$$

The slope is $$1 / \sigma$$ and the intercept is $$d' / \sigma$$.

library(cowplot)

plot_grid(
slope_loglik_crit %>%
bind_rows(slope_sensory_crit) %>%
ggplot()+
geom_segment(color = "grey", aes(x = 0, y = 0, xend = 1, yend = 1)) +
geom_line(aes(x = p_fa, y =p_h, color = cond)) +
coord_equal() +
theme(legend.position = "top")
,
slope_loglik_crit %>%
bind_rows(slope_sensory_crit) %>%
mutate(z_fa = qnorm(p_fa),
z_h = qnorm(p_h)) %>%
ggplot()+
geom_line(aes(x = z_fa, y =z_h, color = cond)) +
coord_equal() +
theme(legend.position = "top")
)