[BioC] design matrix paired data

Bernd Klaus bernd.klaus at embl.de
Thu May 8 11:03:35 CEST 2014


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



More information about the Bioconductor mailing list