2

For a particular tool I need a model matrix which allows me to build desired contrasts and which is a full rank, i.e. no columns are linear combinations of other columns.

The experimental design involves four factors:

  • type (T) with two levels, A and B
  • group (G) with two levels, C and T
  • time point (TP) with three levels
  • subject id (SID)

Both A and B samples were taken from each subject. Subjects either belong to group C or T (control vs treatment). Multiple samples were taken from each subject at different time points.

The comparisons I desire to make are always within a type (no comparisons between types). I want to test interaction between time points and groups, so for example (T.TP1-C.TP1)-(T.TP0-C.TP0).

The only problem is that the model matrix must be full rank, and I don't know how to achieve it.

Here is a mock example:

mock <- data.frame(
  ID=paste0("ID", 1:16), 
  type=rep(c("A", "B"), each=8), 
  treatment=rep(c("C", "T"), each=4), 
  tp=c("T1", "T2"), 
  PID=rep(paste0("P.", 1:8), each=2))

which gives

     ID type treatment tp PID
1   ID1    A         C T1 P.1
2   ID2    A         C T2 P.1
3   ID3    A         C T1 P.2
4   ID4    A         C T2 P.2
5   ID5    A         T T1 P.3
6   ID6    A         T T2 P.3
7   ID7    A         T T1 P.4
8   ID8    A         T T2 P.4
9   ID9    B         C T1 P.5
10 ID10    B         C T2 P.5
11 ID11    B         C T1 P.6
12 ID12    B         C T2 P.6
13 ID13    B         T T1 P.7
14 ID14    B         T T2 P.7
15 ID15    B         T T1 P.8
16 ID16    B         T T2 P.8

Normally, without the repeated measures, I would do something like

mock$ttt <- with(mock, paste(type, treatment, tp, sep="_"))
mm <- model.matrix(~ 0 + ttt, mock)

...and then define contrasts (B_T_T2-B_C_T2)-(B_T_T1-B_C_T1) to test an interaction between time points and treatment within type B.

However, I am at loss how to do it with the repeated measures. I tried the following:

mock$type_pid <- paste0(mock$type, "_", mock$PID)
mm <- model.matrix(~ 0 + type_pid + type:tp:treatment, mock)

I get a matrix which is not fully ranked, however I have the coefficients I need for my contrasts. How can I get a fully ranked matrix with the necessary coefficients?

Please note that I am not trying to fit a mixed random model (despite the repeated measures), because my particular setup does not allow for this.

kjetil b halvorsen
  • 63,378
  • 26
  • 142
  • 467
January
  • 6,999
  • 1
  • 32
  • 55
  • You can always guarantee full rank by appending a nonzero multiple of the identity matrix (or any other full-rank matrix of its dimensions) to the model matrix. See https://stats.stackexchange.com/a/164546/919 for details. – whuber May 23 '20 at 15:03

1 Answers1

2

The problem is that it's not possible to tell (purely syntactically) from the formula that the matrix is not of full rank, so model.matrix() can't know this.

I think you either need to write a function that knows more about the structure of your experiment, or you need linear algebra: one possibility is to do

xqr<-qr(X)

qr.X(xqr, ncol=xqr$rank)
Thomas Lumley
  • 21,784
  • 1
  • 22
  • 73