# Get p-value about contrast hypothesis for rectangular matrix

I am trying to create a contrast in R for a rectangular matrix.

pos_c1 = c(1/2,1/2,-1/6,-1/6,-1/6,-1/6,-1/6,-1/6)
posmat.temp = rbind(constant = 1/8,pos_c1)


in order to solve the matrix, I have read somewhere to use the function 'pseudoinverse'

posmat = pseudoinverse(posmat.temp)
posmat = posmat[, -1]
posmat


Then I have associated my matrix with the data:

contrasts(datasheet.complete$LastPosition) = posmat  but when I make my mixed-effect analysis with lmer(), I do not get the result for the contrast, but for each condition separately. To give an idea, with a square-matrix that calls the function solve(), in the output I get a p_value that tells me how much my contrast hypothesis is confirmed by the data analysis. cl_c1 = c(0,0,-1,1) cl_c2 = c(0,-1,0,1) cl_c3 = c(-1,0,0,1) closuremat.temp = rbind(constant = 1/4,cl_c1,cl_c2,cl_c3) closuremat = solve(closuremat.temp) closuremat = closuremat[, -1] closuremat contrasts(datasheet.complete$Closure) = closuremat

model = lmer(Score~Closure*ExpertiseType*LastPosition+(1|Participant)+(1|Item), data = datasheet.complete, REML = TRUE)
summary(model)

-> Closurecl_c1   ....   ....   p_value


I would expect therefore in the output a parameter called LastPositionpos_c1, but instead I get:

LastPosition1   ....   .....   p_value
Lastposition2   ....   .....   p_value
Lastposition3   ....   .....   p_value
....


Do you know what is the right function to use to solve a matrix which is not squared, and to obtain the result for the contrast hypothesis? Or can you suggest me how to move?

Ok, if I write the following:

conditions <- data.frame(Pos=factor(c(2,4,1,2,5,6,7,2,2,2,5,6,3,3,3,8,5,3,4,2,1,4,3,3,2,6,1,8,3,7,5,7,8,3,6,6,1,6,3,7)))
pos_c1 = c(1/2,1/2,-1/6,-1/6,-1/6,-1/6,-1/6,-1/6)
posmat.temp = cbind(pos_c1)
posmat = pseudoinverse(posmat.temp)
posmat = posmat[, -1, drop = FALSE]
posmat


I get:

     [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]
[1,] 0.75 -0.25 -0.25 -0.25 -0.25 -0.25 -0.25


Looking at the other squared contrasts that work fine, it looks to me that this contrast should present 8 rows (conditions) and 1 column (the contrast hypothesis). Am I correct?

Following this fact, I get an error, which makes sense:

contrasts(conditions$Pos) = posmat Error in contrasts<-(*tmp*, value = c(0.75, -0.25, -0.25, -0.25, -0.25, : wrong number of contrast matrix rows  Ok, here is an example of the issue: library(lme4) library(lmerTest) library(corpcor) database <- data.frame( Clos=factor(c(4,4,1,4,4,3,2,1,2,1,2,2,4,3,1,2,1,4,1,3,2,2,4,4,4,4,2,1,4,2,2,1,4,2,4,2,1,4,4,3)), Pos=factor(c(2,4,1,2,5,6,7,2,2,2,5,6,3,3,3,8,5,3,4,2,1,4,3,3,2,6,1,8,3,7,5,7,8,3,6,6,1,6,3,7)), RF=c(8,6,2,9,7,1,7,6,3,4,6,4,5,2,5,5,3,4,1,3,1,2,3,1,2,2,3,1,8,5,2,2,7,1,9,4,5,6,4,2), Score=c(4,3,3,5,4,3,2,4,5,2,2,3,3,4,4,4,3,2,3,3,5,4,3,4,4,2,3,4,3,4,1,2,2,2,3,4,5,3,1,2) ) clos_c1 = c(0,0,-1,1) clos_c2 = c(0,-1,0,1) clos_c3 = c(-1,0,0,1) closmat.temp = rbind(constant = 1/4,clos_c1,clos_c2,clos_c3) closmat = solve(closmat.temp) closmat = closmat[, -1] closmat pos_c1 = c(1/2,1/2,-1/6,-1/6,-1/6,-1/6,-1/6,-1/6) posmat.temp = rbind(pos_c1) posmat = pseudoinverse(posmat.temp) posmat contrasts(database$Clos) = closmat
contrasts(database$Pos) = posmat model = lmer(Score~Clos+Pos+(1|RF), data = database, REML = TRUE) summary(model)  The result of the summary(model) is Fixed effects: Estimate Std. Error df t value Pr(>|t|) (Intercept) 3.16682 0.23250 8.79795 13.621 3.26e-07 *** Closclos_c1 -0.04582 0.61284 28.45101 -0.075 0.9409 Closclos_c2 -0.25304 0.45563 28.59711 -0.555 0.5830 Closclos_c3 -0.07088 0.49369 28.72373 -0.144 0.8868 Pos1 0.99056 0.39208 26.60481 2.526 0.0178 * Pos2 0.03166 0.45087 28.87030 0.070 0.9445 Pos3 0.45739 0.56255 26.22149 0.813 0.4235 Pos4 -0.38971 0.51039 27.77335 -0.764 0.4516 Pos5 0.19800 0.48213 28.97190 0.411 0.6843 Pos6 -0.46234 0.52282 26.99040 -0.884 0.3843 Pos7 0.28973 0.56143 25.70950 0.516 0.6102  As you can see, I have the p_value for the 3 contrast hypothesis named Closclos_c1, _c2 and _c3. I would therefore expect to have the p_value also for the rectangular contrast hypothesis called pos_c1, but instead I have the p_value for 7 out of 8 conditions. What I would expect instead would be something like: Closclos_c1 Closclos_c2 Closclos_c3 Pospos_c1  • Probably posmat[, -1, drop = FALSE] will do. However, you don't need to rbind the constant when using the pseudoinverse. So posmat <- pseudoinverse(rbind(pos_c1)) or equivalently MASS::ginv(rbind(pos_c1)) should be all you need. – statmerkur May 18 '18 at 11:22 • Also, you can name the contrast via colnames(posmat) <- ".pos_c1". – statmerkur May 18 '18 at 11:32 • Thank you @statmerkur, now I have another problem. When setting the contrast for my parameter, I get the error "Error in cm[, 1L:nc] <- value : replacement has length zero" – Luca Danieli May 18 '18 at 12:31 • I can have a look into it if you provide a reproducible example. – statmerkur May 18 '18 at 12:55 • @statmerkur, you are right. I guess that there was a problem that something was "remembered", and got a false feedback. Now that I had to restart R, I get a different message, First of all, I understand that my 'posmat' (8 conditions) presents 7 columns (that I guess should be 8 rows). Secndly, the result of the contrast is "Error in contrasts<-(*tmp*, value = c(0.75, -0.25, -0.25, -0.25, -0.25, : wrong number of contrast matrix rows". I update now this problem in the main text – Luca Danieli May 19 '18 at 18:18 ## 1 Answer When computing the generalized inverse (aka pseudoinverse) of a matrix in order to specify custom contrasts in R one does not need to rbind (not to be confused with cbind) a constant to that matrix. The following code shows how the desired matrix can be computed. Note that I use MASS::ginv() here but corpcor::pseudoinverse() works equally well. conditions <- data.frame(Pos=factor(c(2,4,1,2,5,6,7,2,2,2,5,6,3,3,3,8,5,3,4,2,1,4,3,3,2,6,1,8,3,7,5,7,8,3,6,6,1,6,3,7))) pos_c1 <- c(1/2,1/2,-1/6,-1/6,-1/6,-1/6,-1/6,-1/6) posmat.temp1 <- rbind(pos_c1) posmat1 <- MASS::ginv(posmat.temp1) posmat1 [,1] [1,] 0.75 [2,] 0.75 [3,] -0.25 [4,] -0.25 [5,] -0.25 [6,] -0.25 [7,] -0.25 [8,] -0.25  If you really want to rbind the constant to the matrix and delete it afterwards, this can be achieved with the following code. drop = FALSE is important in this case as it will prevent R from coercing the matrix to a vector. posmat.temp2 <- rbind(constant = 1/8, pos_c1) posmat2 <- MASS::ginv(posmat.temp2) posmat2 <- posmat2[, -1, drop = FALSE] all.equal(posmat1, posmat2) [1] TRUE  It is easy to name the contrast with: colnames(posmat1) <- ".pos_c1"  A nice feature of R comes in handy now, it automatically adds some orthogonal contrasts to fill in the missing columns of the contrast matrix: contrasts(conditions$Pos) <- posmat1
contrasts(conditions$Pos) .pos_c1 1 0.75 -0.2886751 -0.2886751 -0.2886751 -0.2886751 -0.2886751 -0.2886751 2 0.75 0.2886751 0.2886751 0.2886751 0.2886751 0.2886751 0.2886751 3 -0.25 0.8333333 -0.1666667 -0.1666667 -0.1666667 -0.1666667 -0.1666667 4 -0.25 -0.1666667 0.8333333 -0.1666667 -0.1666667 -0.1666667 -0.1666667 5 -0.25 -0.1666667 -0.1666667 0.8333333 -0.1666667 -0.1666667 -0.1666667 6 -0.25 -0.1666667 -0.1666667 -0.1666667 0.8333333 -0.1666667 -0.1666667 7 -0.25 -0.1666667 -0.1666667 -0.1666667 -0.1666667 0.8333333 -0.1666667 8 -0.25 -0.1666667 -0.1666667 -0.1666667 -0.1666667 -0.1666667 0.8333333  Regarding your last example: You do get the p-value for the contrast pos_c1 but corpcor::pseudoinverse() and MASS::ginv(), unlike solve(), drop the column/row names. Therefore you have to manually set the desired name for your contrast. If you add colnames(posmat) <- "pos_c1" before you assign the matrix to the factor you get the same fixed effects table, but now "Pos1" has the expected name "Pospos_c1". [...] pos_c1 = c(1/2,1/2,-1/6,-1/6,-1/6,-1/6,-1/6,-1/6) posmat.temp = rbind(pos_c1) posmat = pseudoinverse(posmat.temp) colnames(posmat) <- "pos_c1" contrasts(database$Clos) = closmat
contrasts(database\$Pos) = posmat
model = lmer(Score~Clos+Pos+(1|RF), data = database, REML = TRUE)
summary(model)
[...]
Fixed effects:
Estimate Std. Error       df t value Pr(>|t|)
(Intercept)  3.16682    0.23250  8.79800  13.621 3.26e-07 ***
Closclos_c1 -0.04582    0.61284 28.45100  -0.075   0.9409
Closclos_c2 -0.25304    0.45563 28.59700  -0.555   0.5830
Closclos_c3 -0.07088    0.49369 28.72400  -0.144   0.8868
Pospos_c1    0.99056    0.39208 26.60500   2.526   0.0178 *
Pos          0.03166    0.45087 28.87000   0.070   0.9445
Pos          0.45739    0.56255 26.22100   0.813   0.4235
Pos         -0.38971    0.51039 27.77300  -0.764   0.4516
Pos          0.19800    0.48213 28.97200   0.411   0.6843
Pos         -0.46234    0.52282 26.99000  -0.884   0.3843
Pos          0.28973    0.56143 25.71000   0.516   0.6102
[...]

• that's great! I have solved one issue. I am going towards my goal. Unfortunately I am still unable to get the p_value for the contrast. Any advice? – Luca Danieli May 20 '18 at 10:16
• @Luca Danieli: Please, again, provide a reproducible example. Are you using lmerTest? – statmerkur May 20 '18 at 17:24
• ok, I have added a complete reproducible example above. Yes, I am using lmerTest – Luca Danieli May 21 '18 at 16:14
• @Luca Danieli: see my edit. Hope this clears things up. – statmerkur May 21 '18 at 17:06
• Yes! This is exactly the situation I had in the beginning, but as I did not understand that the name is fundamental to read the output, I did not give it much importance. I get indeed an error about 'duplicate rows': 'Error in row.names<-.data.frame(*tmp*, value = value) : duplicate 'row.names' are not allowed'. How's that you don't? – Luca Danieli May 21 '18 at 17:31