R/psychometric_tools.R

In calenwalshe/detectability_maps: What the Package Does (One Line, Title Case)

#' Maximimum likelihood objective function
#'
#' @param df
#'
#' @return
#'
#' @examples
f.search.mle <- function(df) {
#' Title
#'
#' @param p_kSp
#' @param p_k0
#' @param p_Beta
#' @param gamma
#'
#' @return
#' @export
#'
#' @examples
  f <- function(p_kSp, p_k0, p_Beta, gamma = 0) {
    df <- lapply(df, FUN = function(x) x %>%
                   group_by(tAmp, SUBJECT, type, experiment_name, condition, X, Y, eccentricity, bCon) %>%
                   summarize(N_trial = n(), .groups = "keep"))
    # log likelihood
    df$HIT$nll   <-  pnorm(1/2  * f.dprime.uni(df$HIT$tAmp, df$HIT$X, df$HIT$Y, p_Beta,p_kSp,p_k0) - gamma, log.p = T)
    df$MISS$nll  <-  pnorm(-1/2 * f.dprime.uni(df$MISS$tAmp, df$MISS$X, df$MISS$Y, p_Beta,p_kSp,p_k0) - gamma, log.p = T)
    df$FA$nll    <-  pnorm(-1/2 * f.dprime.uni(df$FA$tAmp, df$FA$X, df$FA$Y, p_Beta,p_kSp,p_k0) + gamma, log.p = T)
    df$CR$nll    <-  pnorm(1/2  * f.dprime.uni(df$CR$tAmp, df$CR$X, df$CR$Y, p_Beta,p_kSp,p_k0) + gamma, log.p = T)

    df.all <- do.call(rbind, list(df$HIT, df$MISS, df$CR, df$FA))

    df.nll.norm <- df.all %>%
      ungroup() %>%
      group_by(experiment_name, eccentricity) %>%
      dplyr::mutate(nll = N_trial * nll)

    nll <- -sum(df.nll.norm$nll)

    if(!is.nan(nll)){
      return(nll)
    } else {
      return(10^5)
    }
  }

}

#' Uni function in d' units.
#'
#' @param c_Amp
#' @param c_x
#' @param c_y
#' @param p_Beta
#' @param p_kSp
#' @param p_k0
#'
#' @return
#' @export
#'
#' @examples
f.dprime.uni <- function(c_Amp, c_x, c_y, p_Beta, p_kSp, p_k0) {
#' Title
#'
#' @param c_Amp
#' @param c_x
#' @param c_y
#' @param p_Beta
#' @param p_kSp
#' @param p_k0
#'
#' @return
#' @export
#'
#' @examples
  f <- function(c_Amp, c_x, c_y, p_Beta, p_kSp, p_k0) {
    fov <- f.foveation(c_x, c_y, p_kSp, p_k0)

    uncertainty   <- f.uncertainty(c_x, c_y)

    dprime.no_uni <- (c_Amp * fov)^p_Beta # c_Amp is normalized, so it is d'
    dprime.uni    <- log((exp(dprime.no_uni) + uncertainty) / (1 + uncertainty))
  }

  f.1 <- Vectorize(f)
  return(f.1(c_Amp, c_x, c_y, p_Beta, p_kSp, p_k0))
}

#' d' function without the effect of uncertainty.
#'
#' @param c_Amp
#' @param c_x
#' @param c_y
#' @param p_Beta
#' @param p_kSp
#' @param p_k0
#'
#' @return
#' @export
#'
#' @examples
f.dprime.no.uni <- function(c_Amp, c_x, c_y, p_Beta, p_kSp, p_k0) {
#' Title
#'
#' @param c_Amp
#' @param c_x
#' @param c_y
#' @param p_Beta
#' @param p_kSp
#' @param p_k0
#'
#' @return
#' @export
#'
#' @examples
  f <- function(c_Amp, c_x, c_y, p_Beta, p_kSp, p_k0) {
    fov <- f.foveation(c_x, c_y, p_kSp, p_k0)

    uncertainty <- 0
    dprime.no_uni <- (c_Amp * fov)^p_Beta # c_Amp is normalized, so it is d'
    dprime.uni    <- log((exp(dprime.no_uni) + uncertainty) / (1 + uncertainty))
  }
  f.1 <- Vectorize(f)
  return(f.1(c_Amp, c_x, c_y, p_Beta, p_kSp, p_k0))
}

#' For an x,y coordinate we get the uncertainty constant at that location. From a precomputed matrix.
#'
#' @param x
#' @param y
#' @param umap6cpdmatrix
#'
#' @return
#' @export
#'
#' @examples
f.uncertainty <- function(x,y, umap6cpdmatrix=get("umap6cpdmatrix", globalenv())) {

  xOffset <- 1025
  yOffset <- 1025

  x <- x  + xOffset # the purpose is to index the uncertainty matrix correctly. in this script x = 0 means center of the screen/fovea
  y <- -y + yOffset

  uVal <- umap6cpdmatrix[cbind(y,x)]
}

#' For an x,y coordinate produce the ganglion cell spacing. As per Drasdo et al.
#'
#' @param x
#' @param y
#'
#' @return
#' @export
#'
#' @examples
f.spacing <- function(x, y) {
  x <- x / 120 # supplied in pixels but the watson function takes deg.
  y <- y / 120
  #f <- searchR::watson_spacing
  f <- searchR::sp_xy
  f.1 <- Vectorize(f)
  #return(f.1(x, y, binoc = 1, sp_scal = 0.9436633, f0 = 0))
  return(searchR::sp_xy(x, y))
}

#' For a given x,y and fitted psychometric constants we can interpolate the effect of foveation on d'
#'
#' @param p_kSp
#' @param p_k0
#' @param x
#' @param y
#'
#' @return
#'
#' @examples
f.foveation <- function(x, y, p_kSp, p_k0) {

  f <- function(x, y, p_kSp, p_k0) {
    space_targ <- f.spacing(x, y) * 120
    space_fov  <- f.spacing(0, 0) * 120 # spacing in the fovea, units of pixels 120 ppd

    foveation  <- p_k0 * (1 + p_kSp * space_fov) / (1 + p_kSp * space_targ)
  }

  f.1 <- Vectorize(f)
  return(f.1(x, y, p_kSp, p_k0))
}