Question

DESeq2 model matrix question

0

Entering edit mode

dagurim • 0

@dagurim-22721

Last seen 4.1 years ago

Hello I'm quit newbie to RNAseq and R scripts. I was trying to figure out how to manage DESeq's design model matrix. I have 18 samples to analyze and main target is to get DEGs between Con and Disease. This is how design looks:

Batch Donor Disease
1 A Con
2 A Con
3 A Con
1 B Con
2 B Con
3 B Con
1 C Con
2 C Con
3 C Con
1 D Disease
2 D Disease
3 D Disease
1 E Disease
2 E Disease
3 E Disease
1 F Disease
2 F Disease
3 F Disease

I tried ~Disease and ~Disease + Donor:Batch and results was totally different. I'm not sure about using Donor:Batch and don't know about batch effect reducing. Can anyone give me an advice, please? Thank you.

deseq2 • 1.3k views

ADD COMMENT • link updated 4.9 years ago by swbarnes2 ★ 1.4k • written 4.9 years ago by dagurim • 0

score 0 · Answer 1 · 2020-01-16

0

Entering edit mode

swbarnes2 ★ 1.4k

@swbarnes2-14086

Last seen 3 days ago

San Diego

You also need to put your whole DESeq command, so we know what contrast you called.

Did your PCA results suggest to you that you had significant batch or donor effect?

Why did you choose "Donor:Batch" that instead of doing ~Donor + Batch + Disease?

ADD COMMENT • link 4.9 years ago swbarnes2 ★ 1.4k

0

Entering edit mode

Thank you for the comment!

Here is my whole DESeq command

library("tximport")
library("readr")
library("DESeq2")
library(GenomicFeatures)
library("magrittr")

#salmon Tximport & make sample matrix
dir <- file.path("D:/")
samples <- read.table("samples.txt", header=TRUE)
samples$Disease <- factor(rep(c("Control","Disease"),each=9))
rownames(samples) <- samples$run
files <- file.path(dir,"quant", samples$run, "quant.sf")
names(files) <- samples$run

samples$Disease%<>% relevel("Control")

txdb <- makeTxDbFromGFF("gencode.v32.annotation.gtf")
k <- keys(txdb, keytype="TXNAME")
tx_map <- select(txdb, keys = k, columns="GENEID", keytype = "TXNAME")
tx2gene <- tx_map
txi <- tximport(files, type="salmon", tx2gene=tx2gene)

#DESeq2
dds <- DESeqDataSetFromTximport(txi,
                                   colData = samples,
                                   design = ~ Disease + Donor:Batch)

PCA results suggested I had both of batch & donor effects
~Donor + Batch + Disease does not work

ADD REPLY • link 4.9 years ago dagurim • 0

0

Entering edit mode

How does that design "not work"? What would be more helpful that the details of how you loaded your data is the DESeq command and the results command.

ADD REPLY • link 4.9 years ago swbarnes2 ★ 1.4k

0

Entering edit mode

dds <- DESeqDataSetFromTximport(txi,
                                colData = samples,
                                design = ~ Donor + Batch + Disease)
the design formula contains one or more numeric variables with integer values,
specifying a model with increasing fold change for higher values.
did you mean for this to be a factor? if so, first convert
this variable to a factor using the factor() function
Error in checkFullRank(modelMatrix) : 
the model matrix is not full rank, so the model cannot be fit as specified.
One or more variables or interaction terms in the design formula are linear
combinations of the others and must be removed.

Please read the vignette section 'Model matrix not full rank':

As you can see, "model matrix not full rank" message comes out.

I load my salmon quant files which was generated in linux through vignette recomended method. Make them as txi files and annotate with db loaded in gtf .

ADD REPLY • link 4.9 years ago dagurim • 0

0

Entering edit mode

The model matrix for ~Donor + Batch + Disease will not be full rank, is what the OP is saying. This is one of those situations where you have to use the trickeration outlined in the section titled Group-specific condition effects, individuals nested within groups in the DESeq2 vignette.

ADD REPLY • link 4.9 years ago James W. MacDonald 67k

0

Entering edit mode

And for the OP, you do see the part where it says

Please read the vignette section 'Model matrix not full rank':

Right?

ADD REPLY • link 4.9 years ago James W. MacDonald 67k

0

Entering edit mode

Also the part about how your batch variables are numeric rather than factor?

ADD REPLY • link 4.9 years ago James W. MacDonald 67k

0

Entering edit mode

Yes, I am re-reading the vignette. I think I missed Group-specific condition effects, individuals nested within groups part. Thank you very much.

and Batch is factor :)

ADD REPLY • link 4.9 years ago dagurim • 0

0

Entering edit mode

I just added

samples$donor.n <- factor(rep(rep(1:3,each=3),2))

and commanded

dds <- DESeqDataSetFromTximport(txi,
                                   colData = samples,
                                   design = ~ Disease + Disease:donor.n + Disease:Batch)

worked perfectly!

Thanks!

ADD REPLY • link 4.9 years ago dagurim • 0