Model Design

"Sander [guest]" <[email protected]>
Newsgroups gmane.science.biology.informatics.conductor
Message-ID <[email protected]>
Dear edgeR maintainers,

I have some troubles setting up the correct design for my DE experiment. I now how my experiment design should look like and that is like this: 
> design
      subject1 subject2 subject3 subject4
S7309        1        0        0        0
S7310        1        0        0        0
S7311        0        1        0        0
S7312        0        1        0        0
S7313        0        0        1        0
S7314        0        0        1        0
S7315        0        0        0        1
S7316        0        0        0        1
attr(,"assign")
[1] 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$subject
[1] "contr.treatment"

After that i make my contrasts:
> makeContrasts(subject1,subject2,subject3,subject4, levels = design)
          Contrasts
Levels     subject1 subject2 subject3 subject4
  subject1        1        0        0        0
  subject2        0        1        0        0
  subject3        0        0        1        0
  subject4        0        0        0        1

This all looks right since i need to compare the subjects with each other. But when i run my analysis function which looks like this:

library( edgeR );
library( ggplot2 );
library( reshape );
library( FactoMineR );

analyse <- function( counts, design, contrast, name, style ) {
  counts <- counts[ rowSums( counts, na.rm = TRUE ) > 0, ];
  y <- DGEList( counts = counts, genes = rownames( counts ) );
  y <- calcNormFactors( y );
  y <- estimateGLMCommonDisp( y, design );
  y <- estimateGLMTrendedDisp( y, design, df = 5 );
  y <- estimateGLMTagwiseDisp( y, design );
  
  fit <- glmFit( y, design );
  lrt <- glmLRT( fit, contrast = contrast );
  de <- decideTestsDGE( lrt, p = 0.05, adjust = "BH" );
  cpmY <- cpm( y );
  
  daf <- designAsFactor( design );
  orderedDesign <- design[ order( daf, names( daf ) ), ];
  tab <- data.frame(
    row.names = rownames( cpmY ),
    genes = rownames( cpmY ),
    de = de,
    cpmY[ ,order( daf, names( daf ) ) ]
  );
  
  aRepTab <- topTags( lrt, n = nrow( counts ) )$table;
  aRepTab$rank <- 1:nrow( counts );
#  repTab <- tab[ match( aRepTab$genes, rownames( tab ) ), ];
  
  repTab <- merge( aRepTab, tab, by = "genes", sort = FALSE );
  repTab <- repTab[ order( repTab$rank ), ];
#  data.frame(
#    row.names = rownames( aRepTab ),
#    aRepTab,
#    tab[ match( aRepTab$genes, tab$genes ), ]
#  )
  
  list(
    name = name,
    y = y,
    fit = fit,
    lrt = lrt,
    de = de,
    tab = tab,
    style = style,
    repTab = repTab,
    orderedDesign = orderedDesign
  );
}

I got the following error: 

Error in mglmLevenberg(y, design = design, dispersion = dispersion, offset = offset,  : 
  BLAS/LAPACK routine 'DGEMM ' gave error code -13 
5 mglmLevenberg(y, design = design, dispersion = dispersion, offset = offset, 
    weights = weights, coef.start = start, maxit = 250) 
4 glmFit.default(glmfit$counts, design = design0, offset = glmfit$offset, 
    weights = glmfit$weights, dispersion = glmfit$dispersion, 
    prior.count = 0) 
3 glmFit(glmfit$counts, design = design0, offset = glmfit$offset, 
    weights = glmfit$weights, dispersion = glmfit$dispersion, 
    prior.count = 0) 
2 glmLRT(fit, contrast = contrast) at diffexpr.R#15
1 analyse(counts, design, contrast, countsId, style) 

I tried to use different models but i cannot succeed to avoid the error for this comparison. (other comparisons do succeed) Any hints will be very appreciated.

Thanks in advance!

Best regards,

Sander

 -- 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               LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8     LC_MONETARY=en_US.UTF-8   
 [6] LC_MESSAGES=en_US.UTF-8    LC_PAPER=en_US.UTF-8       LC_NAME=C                  LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] splines   stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] FactoMineR_1.26 reshape_0.8.5   ggplot2_1.0.0   reshape2_1.4    edgeR_3.6.7     limma_3.20.8   

loaded via a namespace (and not attached):
 [1] car_2.0-20           cluster_1.15.2       colorspace_1.2-4     digest_0.6.4         grid_3.1.0           gtable_0.1.2        
 [7] htmltools_0.2.4      lattice_0.20-29      leaps_2.9            MASS_7.3-33          munsell_0.4.2        nnet_7.3-8          
[13] plyr_1.8.1           proto_0.3-10         Rcpp_0.11.2          rmarkdown_0.2.49     scales_0.2.4         scatterplot3d_0.3-35
[19] stringr_0.6.2        tools_3.1.0          yaml_2.1.13 

--
Sent via the guest posting facility at bioconductor.org.

_______________________________________________
Bioconductor mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/bioconductor
Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor
lmpx.com only provides a reader for public news (NNTP) servers. It is not affiliated with the servers or forums shown here and is not responsible for the content of articles, which is written by their respective authors.