Thank you very much for your answer, so if I understand well, even I don't
have interaction I need to use the model with the interaction term, isn't
it?

Aroa Suárez Vega
PhD student
Dpto. Producción Animal
Facultad de Veterinaria
Campus de Vegazana s/n
24071 Leon
Telef. 987291000 Ext. 5296


2014-04-04 17:59 GMT+02:00 Michael Love <michaelisaiahlove@gmail.com>:

>
> On Fri, Apr 4, 2014 at 11:30 AM, aroa [guest] <guest@bioconductor.org>wrote:
>
>>
>> We have and experiment to measure the differences in the milk
>> transcriptome for two breeds. We have RNA-seq samples for these two breeds
>> in 4 different time points, for each breed and time point we have three
>> replicates.
>>
>>         Day1    Day2    Day3    Day4
>> Breed1  3       3       3       3
>> Breed2  3       3       3       3
>>
>>
>> We want to use DESeq2 to extract the differential expressed (DE) genes
>> for each time point between the two breeds and the DE genes for each breed
>> comparing the different time points.
>> We have tested for the interaction (~breed+day+breed:day) and we cannot
>> find interaction between the breeds.
>> Now we are running the model like this:
>> design=data.frame(row.names = colnames(milk), breed = c("breed1",
>> "breed1",
>> "breed1","breed2","breed2","breed2","breed1","breed1","breed1","breed2","breed2","breed2","breed1","breed1","breed1","breed2","breed2","breed2","breed1","breed1","breed1","breed2","breed2","breed2"),
>> day =
>> c("D1","D1","D1","D1","D1","D1","D2","D2","D2","D2","D2","D2","D3","D3","D3","D3","D3","D3","D4","D4","D4","D4","D4","D4"))
>>
>
>
> side note: i would recommend putting the phenotypic data in a CSV file
> and reading it in with read.csv.
>
>
>
>
>> dds<- DESeqDataSetFromMatrix(countData= milk, colData= design, design= ~
>> breed + day)
>> dds$breed<- factor(dds$breed, levels=c("breed1","breed2"))
>> dds$day<- factor(dds$day, levels=c("D1","D2","D3","D4"))
>> dds<-DESeq(dds, betaPrior=FALSE)
>> resultsNames(dds)
>> [1] "Intercept"             "breed_breed2_vs_breed1" "day_D2_vs_D1"
>> [4] "day_D3_vs_D1"       "day_D4_vs_D1"
>>
>> We would like to know what is the meaning of the resultsNames?, we
>> understand them like this:
>> Intercept: breed1D1
>> breed_breed2_vs_breed1: breed2-breed1
>> day_D2_vs_D1: (D2-D1)breed1
>> day_D3_vs_D1: (D3-D1)breed1
>> day_D4_vs_D1: (D4-D1)breed1
>>
>>
> If you fit a model without interactions, then the last three terms here
> are, e.g. D2 - D1 across both breeds, not specifically for breed 1.
>
>
>
>> And how we can make the contrast to obtain the results that we want, that
>> it is comparing the different breeds in each time point and the different
>> time points in each breed?
>>
>
> To compare different breeds at each time point, you need to use the model
> with the interaction term. Interaction terms are how you encode the
> hypothesis, say, that the "breed 2 vs 1 on day 2" effect is not simply the
> sum (in log2FC space) of the main "breed 2 vs 1" effect and the main "day 2
> vs day 1" effect.
>
> In version >= 1.3, we have simplified the extraction of results for
> contrasts like this. This version will be released in 10 days, so I'd
> prefer to recommend you upgrade to this version. If you need to stick with
> v1.2 just email me and I can try to give details. For version 1.3, you
> would fit a model with an interaction term. Then the breed 2 vs 1 effect on
> day 2 can be done with a list() argument to contrasts (see ?results in
> DESeq2 version >= 1.3 for more details):
>
> results(dds, contrast=list(c("breed_breed2","breedbreed2.dayD2"),
>    c("breed_breed1","breedbreed1.dayD2"))
>
> Mike
>
>
>
> Thank you in advance.
>>
>>
>>  -- output of sessionInfo():
>>
>> R version 3.0.3 (2014-03-06)
>> Platform: i386-w64-mingw32/i386 (32-bit)
>>
>> locale:
>> [1] LC_COLLATE=French_France.1252  LC_CTYPE=French_France.1252
>> [3] LC_MONETARY=French_France.1252 LC_NUMERIC=C
>> [5] LC_TIME=French_France.1252
>>
>> attached base packages:
>> [1] parallel  stats     graphics  grDevices utils     datasets  methods
>> base
>>
>> other attached packages:
>>  [1] gplots_2.12.1             RColorBrewer_1.0-5        DESeq2_1.2.10
>>  [4] RcppArmadillo_0.4.100.2.1 Rcpp_0.11.1
>> GenomicRanges_1.14.4
>>  [7] XVector_0.2.0             IRanges_1.20.7
>>  BiocGenerics_0.8.0
>> [10] BiocInstaller_1.12.0
>>
>> loaded via a namespace (and not attached):
>>  [1] annotate_1.40.1      AnnotationDbi_1.24.0 Biobase_2.22.0
>> bitops_1.0-6
>>  [5] caTools_1.16         DBI_0.2-7            gdata_2.13.2
>> genefilter_1.44.0
>>  [9] grid_3.0.3           gtools_3.3.1         KernSmooth_2.23-10
>> lattice_0.20-27
>> [13] locfit_1.5-9.1       RSQLite_0.11.4       splines_3.0.3
>>  stats4_3.0.3
>> [17] survival_2.37-7      tools_3.0.3          XML_3.98-1.1
>> xtable_1.7-3
>>
>> --
>> Sent via the guest posting facility at bioconductor.org.
>>
>> _______________________________________________
>> Bioconductor mailing list
>> Bioconductor@r-project.org
>> https://stat.ethz.ch/mailman/listinfo/bioconductor
>> Search the archives:
>> http://news.gmane.org/gmane.science.biology.informatics.conductor
>>
>
>

	[[alternative HTML version deleted]]

