[BioC] design matrix paired data

[Bernd Klaus]

Hi Ninni,

yes I overlooked that line ;-)

Bernd

On Thu, 8 May 2014 13:38:18 +0200
Ninni Nahm <ninninahm at gmail.com> wrote:

> Hi Bernd!
> 
> Thanks again for your help!
> Maybe you did not see it, but I have this line here:
> 
> colnames(design) <- c("intercept", "patient1", "patient2", "treatment")
> 
> to rename the column names, so the
> 
> tt  <- topTable(fitNoAdj, coef= "treatment", adjust="fdr")
> 
> command should work :) (unless I misunderstood you)
> 
> Thank you :)
> 
> Ninni
> 
> 
> 
> On Thu, May 8, 2014 at 11:03 AM, Bernd Klaus <bernd.klaus at embl.de> wrote:
> 
> > Hi Ninni,
> >
> > looks good!
> >
> > However I would suggest removing the
> > "clutter" in front of the design matrix column
> > names.
> >
> > The command
> >
> > tt  <- topTable(fitNoAdj, coef= "treatment", adjust="fdr")
> >
> > will probably not work since the last column is called
> >
> > "factor(treatment, levels = c("control", "treatment"))treatment"
> >
> > and not
> >
> > "treatment". Another reason to rename the column names of
> > the design matrix ...
> >
> > Best wishes,
> >
> > Bernd
> >
> > On Thu,  8 May 2014 01:18:36 -0700 (PDT)
> > "Ninni Nahm \[guest\]" <guest at bioconductor.org> wrote:
> >
> > > Dear all,
> > >
> > > I made a design matrix according to the limma user guide. However, I
> > wanted to check with you, whether it is correct or not. I feel like I made
> > a mistake, but I do not know for sure.
> > >
> > >
> > > I have 3 patients and 2 conditions:
> > >
> > >
> > > df <- data.frame(patient = c("p1","p2","p3", "p1","p2","p3"),
> > >               treatment = c("control", "control", "control",
> > "treatment", "treatment", "treatment")
> > >               )
> > >
> > > df
> > >   patient treatment
> > > 1      p1   control
> > > 2      p2   control
> > > 3      p3   control
> > > 4      p1 treatment
> > > 5      p2 treatment
> > > 6      p3 treatment
> > >
> > >
> > > design <-  model.matrix(~factor(patient) + factor(treatment, levels =
> > c("control","treatment")), data= df)
> > >
> > >
> > >
> > >
> > > design
> > >   (Intercept) factor(patient)p2 factor(patient)p3
> > > 1           1                 0                 0
> > > 2           1                 1                 0
> > > 3           1                 0                 1
> > > 4           1                 0                 0
> > > 5           1                 1                 0
> > > 6           1                 0                 1
> > >   factor(treatment, levels = c("control", "treatment"))treatment
> > > 1                                                              0
> > > 2                                                              0
> > > 3                                                              0
> > > 4                                                              1
> > > 5                                                              1
> > > 6                                                              1
> > > attr(,"assign")
> > > [1] 0 1 1 2
> > > attr(,"contrasts")
> > > attr(,"contrasts")$`factor(patient)`
> > > [1] "contr.treatment"
> > >
> > > attr(,"contrasts")$`factor(treatment, levels = c("control",
> > "treatment"))`
> > > [1] "contr.treatment"
> > >
> > > colnames(design) <- c("intercept", "patient1", "patient2", "treatment")
> > >
> > >
> > > fit <- lmFit(eset, design)
> > > fit <- eBayes(fitNoAdj)
> > > tt  <- topTable(fitNoAdj, coef= "treatment", adjust="fdr")
> > >
> > >
> > > Is this correct?
> > >
> > > Thank you in advance!
> > >
> > > Ninni
> > >
> > >
> > >
> > >  -- output of sessionInfo():
> > >
> > > > sessionInfo()
> > > R version 3.1.0 (2014-04-10)
> > > Platform: x86_64-pc-linux-gnu (64-bit)
> > >
> > > locale:
> > >  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
> > >  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
> > >  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
> > >  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
> > >  [9] LC_ADDRESS=C               LC_TELEPHONE=C
> > > [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
> > >
> > > attached base packages:
> > > [1] parallel  stats     graphics  grDevices utils     datasets  methods
> > > [8] base
> > >
> > > other attached packages:
> > >  [1] annotate_1.42.0                       pd.hugene.2.0.st_3.8.1
> > >  [3] oligo_1.28.0                          Biostrings_2.32.0
> > >  [5] XVector_0.4.0                         IRanges_1.22.6
> > >  [7] oligoClasses_1.26.0                   RColorBrewer_1.0-5
> > >  [9] xtable_1.7-3
> >  hugene20sttranscriptcluster.db_2.14.0
> > > [11] org.Hs.eg.db_2.14.0                   RSQLite_0.11.4
> > > [13] DBI_0.2-7                             AnnotationDbi_1.26.0
> > > [15] GenomeInfoDb_1.0.2                    Biobase_2.24.0
> > > [17] BiocGenerics_0.10.0                   limma_3.20.1
> > >
> > > loaded via a namespace (and not attached):
> > >  [1] affxparser_1.36.0     affyio_1.32.0         BiocInstaller_1.14.2
> > >  [4] bit_1.1-12            codetools_0.2-8       ff_2.2-13
> > >  [7] foreach_1.4.2         GenomicRanges_1.16.2  iterators_1.0.7
> > > [10] preprocessCore_1.26.0 splines_3.1.0         stats4_3.1.0
> > > [13] XML_3.98-1.1          zlibbioc_1.10.0
> > >
> > > --
> > > Sent via the guest posting facility at bioconductor.org.
> > >
> > > _______________________________________________
> > > Bioconductor mailing list
> > > Bioconductor at r-project.org
> > > https://stat.ethz.ch/mailman/listinfo/bioconductor
> > > Search the archives:
> > http://news.gmane.org/gmane.science.biology.informatics.conductor
> >