3

I am trying to reproduce the SPSS output for significance a linear trend among means when equal variances are not assumed. I have gratefully used code from http://www-personal.umich.edu/~gonzo/coursenotes/file3.pdf to create a function for calculating separate variances, which based on my searching I understand as the “equal variances not assumed” output in SPSS.

My problem/goal: I am only assessing polynomial orthogonal trends (mostly linear). I want to adapt the code creating the function so that the contrast argument can take pre-made contrast matrices rather than manually specifying the coefficients each time (room for typos!).

… I have tried those exact commands but receive Error in contrast %*% means : non-conformable arguments . I have played around with the code but I can’t get it to work.

Code for creating the function from the notes:

sepvarcontrast <- function(dv, group, contrast) {
  means <- c(by(dv, group, mean))
  vars <- c(by(dv, group, var))
  ns <- c(by(dv, group, length))
  ihat <- contrast %*% means
  t.denominator <- sqrt(contrast^2 %*% (vars/ns))
  t.welch <- ihat/ t.denominator
  num.contrast <- ifelse(is.null(dim(contrast)),1,dim(contrast)[1])
  df.welch <- rep(0, num.contrast)
  if (is.null(dim(contrast))) contrast <- t(as.matrix(contrast))
  for (i in 1:num.contrast) {
    num <- (contrast[i,]^2 %*% (vars))^2
    den <- sum((contrast[i,]^2 * vars)^2 / (ns-1))
    df.welch[i] <- num/den
    }
  p.welch <- 2*(1- pt(abs(t.welch), df.welch))
  result <- list(ihat = ihat, se.ihat = t.denominator, t.welch = t.welch,
                 df.welch = df.welch, p.welch = p.welch)
  return(result)
}

I would like to be able to use the function like this:

# Create a polynomial contrast matrix for 5 groups, then save 

contr.mat5 <- contr.poly(5)


# Calculate separate variance

sepvarcontrast(dv, group, contrast = contr.mat5)

I have tried those exact commands to see if they would work but receive Error in contrast %*% means : non-conformable arguments.

All suggestions are appreciated! I am still learning how to create a reprex...

Cassandra
  • 137
  • 1
  • 9
  • 2
    This is hard to test without example data. But it looks like the function is set up to take multiple contrasts in a matrix if each row is a contrast. `contr.poly()` (and other contrast functions) give each contrast in a column. So if you just transpose with `contr.mat5 <- t(contr.poly(5))` I think this may work. – Marius Jun 19 '19 at 00:23
  • Thank you for the suggestion! Yes, transposing the matrix worked exactly as I wanted. – Cassandra Jun 19 '19 at 17:54
  • `tapply` can possibly replace `c` + `by`. Can check with [reproducible example](https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example/5965451). – Parfait Jun 19 '19 at 18:26

0 Answers0