the model matrix is not full rank : This is a classic question which a biologist face without clear understanding of the model design

Saw this answer biostar. Tried to make my metadata as such but still the

"Error in checkFullRank(modelMatrix) : "


This is my coldata

dput(coldata)
structure(list(Group = structure(c(2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L,
6L, 6L, 6L, 6L), .Label = c("HSC", "Blast", "CMP", "GMP", "LSC",
"Mono"), class = "factor"), Mutation = structure(c(10L, 4L, 3L,
5L, 2L, 3L, 9L, 11L, 1L, 1L, 7L, 7L, 10L, 7L, 8L, 2L, 10L, 10L,
10L, 10L, 2L, 2L, 7L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L,
6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 2L, 2L, 5L, 2L,
9L, 1L, 1L, 8L, 6L, 6L, 6L, 6L, 6L, 6L), .Label = c("DNMT3A-R882H",
"FLT3-ITD", "IDH2-R140Q", "KRAS-G13D", "MAU2-Q98H", "No", "NPM1-L287ins(TCTG)",
"NRAS-G13D", "STAG2-R614", "TET2-E1357", "TET2-R550"), class = "factor"),
Satus = structure(c(2L, 5L, 4L, 6L, 3L, 4L, 9L, 10L, 1L,
1L, 7L, 7L, 2L, 7L, 8L, 3L, 2L, 2L, 2L, 2L, 3L, 3L, 7L, 16L,
16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L,
16L, 16L, 17L, 17L, 17L, 17L, 17L, 17L, 17L, 12L, 12L, 13L,
12L, 9L, 11L, 11L, 14L, 15L, 15L, 15L, 15L, 15L, 15L), .Label = c("Blast-DNMT3A-R882H",
"Blast-E1357", "Blast-FLT-ITD", "Blast-IDH2-R140Q", "Blast-KRAS-G13D",
"Blast-MAU2-Q98H", "Blast-NPM1-L287ins(TCTG)", "Blast-NRAS-G13D",
"Blast-STAG2-R614", "Blast-TET2-R550", "LSC-FLT-DNMT3A-R882H",
"LSC-FLT-ITD", "LSC-MAU2-Q98H", "LSC-NRAS-G13D", "Mature-Normal",
"Progenitor-Normal", "Stem-Normal"), class = "factor")), row.names = c("Blast1",
"Blast10", "Blast11", "Blast12", "Blast13", "Blast14", "Blast15",
"Blast16", "Blast17", "Blast18", "Blast19", "Blast20", "Blast2",
"Blast21", "Blast22", "Blast23", "Blast3", "Blast4", "Blast5",
"Blast6", "Blast7", "Blast8", "Blast9", "CMP1", "CMP2", "CMP3",
"CMP4", "CMP5", "CMP6", "CMP7", "CMP8", "GMP1", "GMP2", "GMP3",
"GMP4", "GMP5", "GMP6", "GMP7", "HSC1", "HSC2", "HSC3", "HSC4",
"HSC5", "HSC6", "HSC7", "LSC1", "LSC2", "LSC3", "LSC4", "LSC5",
"LSC6", "LSC7", "LSC8", "Mono1", "Mono2", "Mono3", "Mono4", "Mono5",
"Mono6"), class = "data.frame")

dds <- DESeqDataSetFromMatrix(countData=df2, colData=coldata, design= ~ Group + Mutation)


This is the design I'm trying to incorporate. What i get is my design is nested in one of the column by looking at numerous example.

What i would like to do is I'm trying to run a LRT test on this dataset which works fine if i give "Group" as my design and then i run this

dds_lrt <- DESeq(dds, test="LRT", reduced = ~ 1)


which gives me diffrenitaly expressed across multiple groups. Now there are two question

1. Is it possible to include interaction such as in my case Group+Mutation and do the lrt test?
2. Or the above is conceptually wrong to do ?
• what is the hypothesis you are trying to test? Dec 18 '20 at 13:13
• I would like to know how these specific mutation affect gene expression across the cell types with respect to reference which is my "HSC" as its a stem cell . But from @ATpoint answer it seems I need to subset the data more over what i tired above perhaps is a wrong approach if i get it as there are multiple types of mutation. It would be more useful if i do pairwise comparison for specific cell types i guess
– kcm
Dec 18 '20 at 13:16
• @StupidWolf for example I would like to see what set of genes which changes or are significant when i traverse from HSC to down the lineage tree such as CMP,GMP,Monocyte, but I have also included LSC(Leukemic stem cells) and Blast. But not sure if this is right comparison, i did normalize and took the data for clustering that is fine ,but hypothesis testing across heterogeneous group which are diseased with mutation would that be correct?
– kcm
Dec 18 '20 at 13:20
• Looking at your design with ATpoint's table, there's a major issue. some of your mutations have n=1 in some of the interactions, for example TET2-R550 or NRAS-G13D. or NPM1-L287ins(TCTG) is exclusive for blast. So it's hard for you to make these comparisons say HSC.. when it is so incomplete Dec 19 '20 at 1:04
• I think you basically need to think about how to design this Dec 19 '20 at 1:05

> table(doldat$$Group, coldat$$Mutation)

DNMT3A-R882H FLT3-ITD IDH2-R140Q KRAS-G13D MAU2-Q98H No NPM1-L287ins(TCTG) NRAS-G13D STAG2-R614 TET2-E1357 TET2-R550
HSC              0        0          0         0         0  7                  0         0          0          0         0
Blast            2        4          2         1         1  0                  4         1          1          6         1
CMP              0        0          0         0         0  8                  0         0          0          0         0
GMP              0        0          0         0         0  7                  0         0          0          0         0
LSC              2        3          0         0         1  0                  0         1          1          0         0
Mono             0        0          0         0         0  6                  0         0          0          0         0


Look at this table, you see by eye how some mutations are nested with Group. If you want to include it then you will need to subset the experiment into multiple separate runs, depending on the question you want to answer.

• Lets say if i want to see the effect of mutation+Group across HSC , LSC and Blast would the above metadata design work ?
– kcm
Dec 18 '20 at 14:57