DESeq2 coefficient names
1
0
Entering edit mode
jkhudyakov • 0
@jkhudyakov-23010
Last seen 10 days ago
United States

I am having trouble getting specific coefficient names after running DESeq2. The coefficients I get when I run resultsNames(dds) are "Treat1", "Treat2", etc instead of "Treat_B_vs_A", "Treat_C_vs_A", etc. Has anyone else had this issue? The code that I was running and the output are below:

samples$Treat <- factor(samples$Treat, levels = c("A","B","C","D","E","F"))
samples$Treat

# A   A   A   B   B   B   C   C   C   D   D   D   E   E   E   F   F   F
# Levels: A   B   C   D   E   F

dds <- DESeqDataSetFromTximport(txi, 
                                colData = samples, 
                                design = ~ Subject + Treat)

keep <- rowSums(counts(dds) >= 10) >= 3
dds <- dds[keep,]

dds$Treat <- relevel(dds$Treat, ref = "A")

dds <- DESeq(dds) 

resultsNames(dds)
# "Intercept" "Subject1"  "Subject2"  "Treat1"    "Treat2"    "Treat3"    "Treat4"    "Treat5" 

sessionInfo( )
# R version 3.6.0 (2019-04-26)
# Platform: x86_64-apple-darwin15.6.0 (64-bit)
# Running under: macOS  10.15.7

Thank you!

resultsNames DESeq2 coefficient • 162 views
ADD COMMENT
0
Entering edit mode
@mikelove
Last seen 11 hours ago
United States

Looks strange. And Subject is a simple factor also? What about class(samples)?

ADD COMMENT
0
Entering edit mode

Hi Michael,

Yes, Subject is also a simple factor. samples$Subject returns:

12 14 15 12 14 15 12 14 15 12 14 15 12 14 15 12 14 15 12 14 15 12 14 15 Levels: 12 14 15

class(samples) returns

"data.frame"

The strange thing is that I had it working properly last week, but it stopped giving specific coefficients when I ran it this week. I am not sure what I did differently - nothing I don't think.

ADD REPLY
0
Entering edit mode

Check if your installation is still ok: BiocManager::valid()

ADD REPLY
0
Entering edit mode

I get:

Bioconductor version '3.9'

  • 115 packages out-of-date
  • 0 packages too new

Versions of Bioconductor packages used for the analysis: tximport_1.12.3 tximportData_1.12.0 DESeq2_1.24.0 apeglm_1.6.0 IHW_1.12.0 rhdf5_2.28.1

I did the package update recommended by Bioconductor and re-ran DESeq2 but am still getting the weird coefficient names.

ADD REPLY
0
Entering edit mode

I'm not sure, could you try with the same data on a different machine?

Also what do you get with:

model.matrix(~Treat, data=samples)
ADD REPLY
0
Entering edit mode

I have 8 conditions with 3 biol reps each. I releveled the treatment factor so that every sample is compared to A (control).

As a (slight) aside, how would I specify the comparison B vs C if Treat2 is B vs A and Treat 2 is C vs A? Is this correct? contrast = c("Treat", "Treat3", "Treat2")

A_rep1 1 -1 -1 -1 -1 -1 -1 -1

A_rep2 1 -1 -1 -1 -1 -1 -1 -1

A_rep3 1 -1 -1 -1 -1 -1 -1 -1

B_rep1 1 1 -1 -1 -1 -1 -1 -1

B_rep2 1 1 -1 -1 -1 -1 -1 -1

B_rep3 1 1 -1 -1 -1 -1 -1 -1

C_rep1 1 0 2 -1 -1 -1 -1 -1

C_rep2 1 0 2 -1 -1 -1 -1 -1

C_rep3 1 0 2 -1 -1 -1 -1 -1

D_rep1 1 0 0 3 -1 -1 -1 -1

D_rep2 1 0 0 3 -1 -1 -1 -1

D_rep3 1 0 0 3 -1 -1 -1 -1

E_rep1 1 0 0 0 4 -1 -1 -1

E_rep2 1 0 0 0 4 -1 -1 -1

E_rep3 1 0 0 0 4 -1 -1 -1

F_rep1 1 0 0 0 0 5 -1 -1

F_rep2 1 0 0 0 0 5 -1 -1

F_rep3 1 0 0 0 0 5 -1 -1

G_rep1 1 0 0 0 0 0 6 -1

G_rep2 1 0 0 0 0 0 6 -1

G_rep3 1 0 0 0 0 0 6 -1

H_rep1 1 0 0 0 0 0 0 7

H_rep2 1 0 0 0 0 0 0 7

H_rep3 1 0 0 0 0 0 0 7

attr(,"assign")

[1] 0 1 1 1 1 1 1 1

attr(,"contrasts")

attr(,"contrasts")$Treat

[1] "contr.helmert"

ADD REPLY
0
Entering edit mode

I see -- the issue is that you've set a non-default setting for contrasts in your .Rprofile.

I'd recommend to take that out if you want default contrasts from DESeq2.

ADD REPLY
0
Entering edit mode

Ahh that's right, I used helmert contrasts for a different project and didn't realize it changed it in the .Rprofile. I set it back to default and now it works. Mystery solved! Thanks so much for your prompt help and for writing these tools!

ADD REPLY

Login before adding your answer.

Traffic: 403 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6