Question

DESeq2: Extracting Averaged Counts Information from dds and rld

0

Entering edit mode

vanbelj ▴ 30

@vanbelj-21216

Last seen 10 months ago

United States

I'd like to extract count information from my DESeqDataSet that has been averaged by Treatment. I'd like to do the same for my rlog normalized counts. I know how to access the counts data in each case, but I'm unsure how to group samples by Treatment before generating averages. I'd like to be able to recycle this regardless of Treatment names.

Thanks for any help.

>samples
| Identifier   | Treatment |
|--------------|-----------|
| Control_Rep1 | Control   |
| Ethanol_Rep1 | Ethanol   |
| Control_Rep2 | Control   |
| Ethanol_Rep2 | Ethanol   |
| Control_Rep3 | Control   |
| Ethanol_Rep3 | Ethanol   |

control_name <- "Control" 
dds$Treatment <- relevel(dds$Treatment, ref = control_name)
dds <- DESeqDataSetFromTximport(txi,
                                colData = samples,
                                design = ~ Treatment)

>head(counts(dds_filtered), 3)
|        | Control_Rep1 | Ethanol_Rep1 | Control_Rep2 | Ethanol_Rep2 | Control_Rep3 | Ethanol_Rep3 |
|--------|--------------|--------------|--------------|--------------|--------------|--------------|
| Q0120  | 22           | 29           | 25           | 39           | 13           | 23           |
| R0010W | 3694         | 6205         | 3322         | 7110         | 4985         | 10513        |
| R0020C | 3024         | 3564         | 2799         | 4191         | 5030         | 6214         |

>Desired_Table
|        | Control    | Ethanol    |
|--------|------------|------------|
| Q0120  | 20         | 30.3333333 |
| R0010W | 4000.33333 | 7942.66667 |
| R0020C | 3617.66667 | 4656.33333 |

DESeq2 • 1.6k views

ADD COMMENT • link 2.6 years ago vanbelj ▴ 30

score 2 · Accepted Answer · 2022-05-11

2

Entering edit mode

ATpoint ★ 4.6k

@atpoint-13662

Last seen 7 hours ago

Germany

There are many ways, but the simplest is probably to use the avereps() function from limma, see ?limma::avereps.

ADD COMMENT • link 2.6 years ago ATpoint ★ 4.6k

0

Entering edit mode

Looks like that is for an S3 class, whereas the DESeq2 datasets are all S4

ADD REPLY • link 2.6 years ago vanbelj ▴ 30

0

Entering edit mode

It accepts a count matrix, so take the transformed (or whatever) counts and feed it into that function. As said, ?avereps:

(...)
x       a matrix-like object, usually a matrix, MAList or EList object.

Don't concern yourself with these class definitions, as an end user you almost never have to care about it.

ADD REPLY • link 2.6 years ago ATpoint ★ 4.6k

score 1 · Accepted Answer · 2022-05-17

Below is a solution to extract and average counts across replicates from a dds using Base R. If you're familiar with limma you may prefer ATpoint's solution above.

dds_df <- as.data.frame(counts(dds, normalized=FALSE))
dds_counts <- sapply(split.default(dds_df, samples$Treatment[match(names(dds_df), samples$Identifier)]),
                     rowMeans, na.rm = FALSE)

> head(dds_counts, 3)
            Control  Ethanol
YAL068C   0.0000000 3.666667
YAL067W-A 0.3333333 0.000000
YAL067C   5.6666667 8.333333