R/lensModel.r

lensModel <- function(inSet, exSet, cueSet) {
  lens.data <- data.frame(inSet, exSet, cueSet)
  inLM <- apply(inSet, 2, function(x) lm(as.formula(paste("x ~ ", paste(names(cueSet), collapse="+"))), data=lens.data, singular.ok=FALSE))
  exLM <- apply(exSet, 2, function(x) lm(as.formula(paste("x ~ ", paste(names(cueSet), collapse="+"))), data=lens.data, singular.ok=FALSE))
  coef.in <- matrix(unlist(lapply(inLM, coef)), ncol=ncol(inSet), byrow=F)
  coef.ex <- matrix(unlist(lapply(exLM, coef)), ncol=ncol(exSet), byrow=F)
  Rv <- unlist(lapply(inLM, function(x) sqrt(summary(x)$r.square)))
  Ru <- unlist(lapply(exLM, function(x) sqrt(summary(x)$r.square)))
  fit.in <- apply(t(coef.in)[,-1], 1, function(x) rowSums(matrix(x, nrow=nrow(cueSet), ncol=ncol(cueSet), byrow=T) * as.matrix(cueSet))) + matrix(rep(coef.in[1,], each=nrow(cueSet)), nrow=nrow(cueSet), ncol=ncol(coef.in), byrow=F)
  fit.ex <- apply(t(coef.ex)[,-1], 1, function(x) rowSums(matrix(x, nrow=nrow(cueSet), ncol=ncol(cueSet), byrow=T) * as.matrix(cueSet))) + matrix(rep(coef.ex[1,], each=nrow(cueSet)), nrow=nrow(cueSet), ncol=ncol(coef.ex), byrow=F)
  res.in <- inSet - fit.in
  res.ex <- exSet - fit.ex
  p.in <- matrix(unlist(lapply(inLM, function(x) coef(summary(x))[,4])), ncol=ncol(inSet), byrow=F)
  p.ex <- matrix(unlist(lapply(exLM, function(x) coef(summary(x))[,4])), ncol=ncol(inSet), byrow=F)
  coef.cor <- diag(cor(coef.in[-1,], coef.ex[-1,]))
  A <- diag(cor(inSet, exSet, use="pair"))
  G <- diag(cor(fit.in, fit.ex, use="pair"))
  C <- diag(cor(res.in, res.ex, use="pair"))
  mB.in <- colMeans(coef.in[-1,], na.rm=T)
  aaB.in <- colMeans(abs(coef.in[-1,]), na.rm=T)
  mB.ex <- colMeans(coef.ex[-1,], na.rm=T)
  aaB.ex <- colMeans(abs(coef.ex[-1,]), na.rm=T)
  stats <- data.frame(rbind(Rv, Ru, coef.cor, A, G, C, mB.in, aaB.in, mB.ex, aaB.ex), row.names=c("Validity Saturation", "Utilization Saturation", "Coefficient Correlation", "Achievement", "Linear Knowledge", "Unmodeled Knowledge", "Avg. Cue Validity", "Avg. Absolute Cue Validity", "Avg. Cue Utilization", "Avg. Absolute Cue Utilization"))
  colnames(stats) <- colnames(inSet)
  colnames(coef.in) <- colnames(inSet)
  colnames(coef.ex) <- colnames(exSet)
  rownames(coef.in) <- c("Intercept", colnames(cueSet))
  rownames(coef.ex) <- c("Intercept", colnames(cueSet))
  colnames(p.in) <- colnames(inSet)
  colnames(p.ex) <- colnames(exSet)
  rownames(p.in) <- c("Intercept", colnames(cueSet))
  rownames(p.ex) <- c("Intercept", colnames(cueSet))
  out <- list("Lens Stats"=stats, "Cue Validities"=coef.in, "Cue Utilizations"=coef.ex, "Validity p"=p.in, "Utilization p"=p.ex)
  class(out) <- c("lensMod")
  out
}

Try the multicon package in your browser

Any scripts or data that you put into this service are public.

multicon documentation built on May 2, 2019, 3:18 a.m.