Search
Question: coefficients meaning in time course analysis with DESeq2
0
gravatar for Assa Yeroslaviz
3.0 years ago by
Assa Yeroslaviz1.4k
Munich, Germany
Assa Yeroslaviz1.4k wrote:

Hi,

I decided to split DE genes over time with DESeq2 into several parts, as I figured out that it is too long. I hope it was a good idea.

I'm doing a time course analysis with deseq2, I have 8 different points and multiple replica.

when I am looking at resultsNames(dds) in my analysis I see

>resultsNames(dds)

# [1] "Intercept"      "replica_2_vs_1" "replica_3_vs_1" "time_16h_vs_0h" "time_1d_vs_0h"
# [6] "time_24h_vs_0h" "time_30h_vs_0h" "time_48h_vs_0h" "time_72h_vs_0h" "time_90h_vs_0h"

but when I do results(dds), I get only the last time point

> res
log2 fold change (MLE): time 90h vs 0h
LRT p-value: '~ replica + time' vs '~ replica'
DataFrame with 17558 rows and 6 columns
...

Does it means that with this coefficient I get all the genes which are changing over all TPs?

How can I interpret the other coefficients? Are the changes over time up to a specific TP (e.g. 16h_vs_0h) means all genes changing from 0h to 16h?

Is it possible to extract with this kind of analysis also all pair-wise comparisons (like 1d vs. 16h, or 24h vs. 48h), or do i need to run a "normal" wald test for the pair-wise comparisons?

 

the script I run is below:

#reading in the phenotypic data from a file
phenotype <- read.delim2("phenotypeData.txt", colClasses = c(rep("factor", 3)))
# create the DESeq object from a matrix.
dds<-DESeqDataSetFromMatrix(countData=countTable, colData=phenotype, design= ~ replica + time )
dds = DESeq(dds, test="LRT", reduced=~replica)

res<-results(dds)
ADD COMMENTlink modified 3.0 years ago by Michael Love19k • written 3.0 years ago by Assa Yeroslaviz1.4k
0
gravatar for Steve Lianoglou
3.0 years ago by
Denali
Steve Lianoglou12k wrote:

The results function takes a "contrast" or "name" argument so that you can specify the result you are interested in, cf the help page via ?results

ADD COMMENTlink written 3.0 years ago by Steve Lianoglou12k
0
gravatar for Assa Yeroslaviz
3.0 years ago by
Assa Yeroslaviz1.4k
Munich, Germany
Assa Yeroslaviz1.4k wrote:

@Steve - thanks for the response. I did read the help page and the user guide, but I'm still not sure if I understand it correctly.

From the list I get from resultsNames() I can use in the results() command with the name= parameter. Everything else I can constract with the contracts= parameter?

Does this means, that the contrast=c("time", "16", "0") and name="time_16h_vs_0h" will give me the same results?

Back from to my question from above. So if I understand it now correctly, I can just use contrast=c("time", "24", "48") to get all genes which are changed between the TP 48 hours to 24 hours, even though it is not in th resultNames() object.

Is it better to use here the Wald test, when comparing only two TP agianst each other?

Assa

ADD COMMENTlink modified 3.0 years ago • written 3.0 years ago by Assa Yeroslaviz1.4k
On Nov 5, 2015 3:56 AM, "Assa Yeroslaviz [bioc]" <noreply@bioconductor.org> wrote: > > Activity on a post you are following on support.bioconductor.org > > User Assa Yeroslaviz wrote Answer: coefficients meaning in time course analysis with DESeq2: > > @Steve - thanks for the response. I did read the help page and the user guide, but I'm still not sure if I understand it correctly. > > From the list I get from resultsNames() I can use in the results() command with the name= parameter. Everything else I can constract with the contracts= parameter? > > Does this means, that the contrast=c("time", "16", "0") and name="time_16h_vs_0h" will give me the same results? > Yes. > Assa > > ________________________________ > > Post tags: deseq2, timecourse, time course, design and contrast matrix, coefficient of variation > > You may reply via email or visit A: coefficients meaning in time course analysis with DESeq2
ADD REPLYlink written 3.0 years ago by Michael Love19k

Thanks Michael, I thought so.

Does it make than more sense to run the results() command for a pair-wise comparison using the Wald test?

I have run both and see, that the Wald statistics gives me in general lower adj. p-values. If I understand it correctly, the Wald statistics is doing an LFC shrinkage, while the LRT takes multiple parameters into account, when calculating the (adj.) p-values.

But which one of the two is "better" for a pair-wise comparison?

 

Thanks

Assa

 

ADD REPLYlink written 3.0 years ago by Assa Yeroslaviz1.4k
0
gravatar for Michael Love
3.0 years ago by
Michael Love19k
United States
Michael Love19k wrote:

Here's the relevant section from ?results

"For analyses using the likelihood ratio test (using nbinomLRT), the p-values are determined solely
by the difference in deviance between the full and reduced model formula. A log2 fold change is
included, which can be controlled using the name argument, or by default this will be the estimated
coefficient for the last element of resultsNames(object)."

There is also a section in the vignette which might help explain things:

"3.5 Likelihood ratio test

DESeq2 offers two kinds of hypothesis tests: the Wald test, where we use the estimated standard error of
a log2 fold change to test if it is equal to zero, and the likelihood ratio test (LRT)..."

For the likelihood ratio test results table, you can just disregard the LFC column and think of this as giving you a single p-value for a set of fold changes.

So for DESeq2:

Wald test for testing specific fold changes or linear combinations (condition specific changes at a particular time, using 'name' and telling results that you want test="Wald")

LRT for testing whether sets of fold changes are all equal to zero (any condition specific changes over time)

You can also take a look at the workflow which has a time series example:

www.bioconductor.org/help/workflows/rnaseqGene/#time

I might suggest you talk to a local statistician if you have further questions on the difference between a Wald test and a likelihood ratio test, or on how to interpret the various results tables. The coefficients given by DESeq2 here are standard for any kind of a linear model.

ADD COMMENTlink modified 3.0 years ago • written 3.0 years ago by Michael Love19k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.2.0
Traffic: 207 users visited in the last hour