LAPACK/BLAS error with edgeR GLM
5
0
Entering edit mode
@markrogers-8451
Last seen 6.8 years ago

I am trying to run edgeR with the GLM on some data that I've already used with the 'classic' edgeR analysis, as I would like to try it with RUVSeq.  However, I am getting the following errors when I attempt to use the edgeR estimateGLMCommonDisp function:

Error in mglmLevenberg(y, design = design, dispersion = dispersion, offset = offset,  :
BLAS/LAPACK routine 'DSYRK ' gave error code -10

If I use estimateCommonDisp instead, it works just fine, so I don't think it's a data issue.  I have reinstalled the LAPACK libraries several times and reinstalled edgeR after each time, but the error remains.  Error code -10 doesn't tell me much(!) so I'm about at my wit's end.

lapack/blas version: 3.5.0 (Nov. 2013)

edgeR version: 3.10.2

Linux Ubuntu 14.04

Any ideas at all??  Many thanks in advance!

glm edger lapack segfault ruvseq • 1.7k views
0
Entering edit mode

As Davide says, can you give us the design matrix that causes the problem, and a set of counts for one gene that can reproduce the error? DSYRK is called by DPOTRF for the Cholesky decomposition, so I guess that's where the problem is occurring within mglmLevenberg.

0
Entering edit mode

I've roughly narrowed it to the following code that I copied/pasted from the RUVSeq user's guide.  I should emphasize that this is not a problem with RUVSeq, as I have another example with a similar issue that doesn't use it at all.  However, since you're the one who responded...

library(RUVSeq)
library(zebrafishRNASeq)
library(RColorBrewer)

data(zfGenes)

filter   <- apply(zfGenes,1,function(x) length(x[x>5])>=2)
filtered <- zfGenes[filter,]
genes    <- rownames(filtered)[grep("^ENS", rownames(filtered))]
spikes   <- rownames(filtered)[grep("^ERCC", rownames(filtered))]
x        <- as.factor(rep(c("Ctl", "Trt"), each=3))
set      <- newSeqExpressionSet(as.matrix(filtered), phenoData=data.frame(x, row.names=colnames(filtered)))
uqset    <- betweenLaneNormalization(set, which="upper")

# Run RUV adjustment on upper-quantile
uqRUVSet <- RUVg(uqset, spikes, k=1)
pData(uqRUVSet)
colors   <- brewer.pal(3, "Set2")
plotRLE(uqRUVSet, outline=FALSE, ylim=c(-4, 4), col=colors[x], main='Old version')

# Run RUV adjustment on spike-ins
spikeRUVSet <- RUVg(set, spikes, k=1)
pData(spikeRUVSet)
plotRLE(spikeRUVSet, outline=FALSE, ylim=c(-4, 4), col=colors[x], main='New version')

# Now use the revised matrix for DE analysis (either one causes error below)
#design <- model.matrix(~x + W_1, data=pData(uqRUVSet))
design <- model.matrix(~x + W_1, data=pData(spikeRUVSet))
y      <- DGEList(counts=counts(set), group=x)
y      <- calcNormFactors(y, method="upperquartile")
# ERROR:
y      = estimateGLMCommonDisp(y, design)

0
Entering edit mode

I've copied/pasted the tail of my output below.  The W_1 values are slightly different, but I don't imagine it is the source of the problem.  I didn't do anything with LAPACK until I got this error the first time; then I thought it must be a library version issue, so I reinstalled the latest LAPACK before I reinstalled R and then the R packages.  However, I was unaware of these other LAPACK nuances (didn't even know about the --with-lapack flag) so I will read that documentation -- many thanks for the link.

estimateGLMCommonDisp
Error in mglmLevenberg(y, design = design, dispersion = dispersion, offset = offset,  :
BLAS/LAPACK routine 'DSYRK ' gave error code -10

> design
(Intercept) xTrt        W_1
Ctl1            1    0 -0.0390975
Ctl3            1    0  0.4088907
Ctl5            1    0  0.4414666
Trt9            1    1 -0.2138095
Trt11           1    1 -0.7527070
Trt13           1    1  0.1552567
attr(,"assign")
[1] 0 1 2
attr(,"contrasts")
attr(,"contrasts")$x [1] "contr.treatment" ADD REPLY 1 Entering edit mode Aaron Lun ★ 27k @alun Last seen 8 hours ago The city by the bay We have been unable to reproduce the problem with estimateGLMCommonDisp on our end. Just to be sure we're talking about the same thing, is this what the design matrix looks like on your machine? > design (Intercept) xTrt W_1 Ctl1 1 0 -0.04539413 Ctl3 1 0 0.50347642 Ctl5 1 0 0.40575319 Trt9 1 1 -0.30773479 Trt11 1 1 -0.68455406 Trt13 1 1 0.12845337 attr(,"assign") [1] 0 1 2 attr(,"contrasts") attr(,"contrasts")$x
[1] "contr.treatment"

On my system (CentOS 6.6, edgeR 3.10.2, LAPACK 3.5.0), the design matrix constructed with spikeRUVSet gives a common dispersion of 0.74, while that with uqRUVSet gives a common dispersion of 0.76.

Also, your original post suggests that you're installing LAPACK separately from R. Be careful with how you do this - perhaps a reading of the R documentation here may be helpful. In particular:

Please do bear in mind that using --with-lapack is ‘definitely not recommended’: it is provided only because it is necessary on some platforms and because some users want to experiment with claimed performance improvements. Reporting problems where it is used unnecessarily will simply irritate the R helpers.

0
Entering edit mode
Hi Aaron, As Gordon requested, I'm sending you the sessionInfo() output and the design matrix as attachments. I still haven't tried to address possible LAPACK install issues as I've got a meeting this morning; I hope to return to it this afternoon. Many thanks to you all! Mark -- Dr. Mark Rogers Research Associate Intelligent Systems Laboratory Department of Engineering Mathematics 0117 33 15328 On 25 July 2015 at 10:49, Aaron Lun [bioc] <noreply@bioconductor.org> wrote: > Activity on a post you are following on support.bioconductor.org > > User Aaron Lun <https: support.bioconductor.org="" u="" 6732=""/> wrote Comment: > LAPACK/BLAS error with edgeR GLM > <https: support.bioconductor.org="" p="" 70249="" #70385="">: > > We have been unable to reproduce the problem with on our end. Just to be > sure we're talking about the same thing, is this what the design matrix > looks like on your machine? > > > design > (Intercept) xTrt W_1 > Ctl1 1 0 -0.04539413 > Ctl3 1 0 0.50347642 > Ctl5 1 0 0.40575319 > Trt9 1 1 -0.30773479 > Trt11 1 1 -0.68455406 > Trt13 1 1 0.12845337 > attr(,"assign") > [1] 0 1 2 > attr(,"contrasts") > attr(,"contrasts")$x > [1] "contr.treatment" > > On my system, the design matrix constructed with spikeRUVSet gives a > common dispersion of 0.74, while that with uqRUVSet gives a common > dispersion of 0.76. > > Also, your original post suggests that you're installing LAPACK separately > from R. Be careful with how you do this - perhaps a reading of the R > documentation here > <https: cran.r-project.org="" doc="" manuals="" r-admin.html#lapack=""> may be > helpful. In particular: > > *Please do bear in mind that using --with-lapack is ‘definitely not > recommended’: it is provided only because it is necessary on some platforms > and because some users want to experiment with claimed performance > improvements. Reporting problems where it is used unnecessarily will simply > irritate the R helpers.* > > ------------------------------ > > Post tags: glm, edger, lapack, segfault, ruvseq > > You may reply via email or visit > C: LAPACK/BLAS error with edgeR GLM > ADD REPLY 1 Entering edit mode No problems on my end: > load("ydesign.RData") > require(edgeR) Loading required package: edgeR Loading required package: limma > y <- estimateGLMCommonDisp(y, design) > y$common
[1] 0.7447239

Here's the output of sessionInfo(), which is mostly the same as your setup:

R version 3.2.1 (2015-06-18)
Platform: x86_64-unknown-linux-gnu (64-bit)
Running under: CentOS release 6.4 (Final)

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
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

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

other attached packages:
[1] edgeR_3.10.2  limma_3.24.14

My best guess is that it's a problem with your edgeR or R installations. See if you get better results after re-installing either or both of them. For R itself, I usually compile it fresh from the CRAN source, but depending on how you've set it up, you might need to apt-get remove it or something equivalent before attempting re-installation.

0
Entering edit mode

Attachments don't seem to come through when you send it via the support site. It would probably be better to send it to me directly off-list, or to put up a link to DropBox (or whatever site you like to use for file sharing).

0
Entering edit mode
davide risso ▴ 890
@davide-risso-5075
Last seen 7 months ago

Hi Mark,

How does the design matrix look like?

My guess is that RUV is actually capturing the biology (i.e. Your genes are not truly negative controls) and this causes a (almost) colinearity in the design. Can you post your design matrix here?

Best,
Davide

0
Entering edit mode
@gordon-smyth
Last seen 1 hour ago
WEHI, Melbourne, Australia

Dear Mark,

As Aaron has told you, we are unable to reproduce the error you report, nor have we seen this error before. Can I suggest that you

1. Show us the entire output from sessionInfo() and

2. Save y and design to a binary file and email this file to Aaron. In other words, type

save(y,design,file="ydesign.RData")

immediately before running estimateGLMCommonDisp() and then email the file 'ydesign.RData' to Aaron.

0
Entering edit mode
@markrogers-8451
Last seen 6.8 years ago

Scanning the system, it's not hard to believe that LAPACK installations are causing a problem.  Although I never explicitly installed it myself (until this recent problem arose), I see liblapack.a and liblapack.so files in several areas (see below).  It appears the libraries are used by CGAL, scipy, numpy and MATLAB, possibly others.

/etc/alternatives
/usr/lib
/usr/lib/atlas-base/atlas
/usr/lib/lapack
/usr/local/lib
/var/lib/dpkg/alternatives

I will try to track down the source of the problem (prefer not to uninstall these other packages, especially MATLAB), but if anyone sees an obvious workaround, I'd love to hear about it.

0
Entering edit mode

The R installation should come with its own LAPACK and BLAS libraries, to be found at:

list.files(file.path(R.home(), "lib"))

I don't think it would use the external libraries unless explicitly requested to do so with --with-lapack during compilation.

0
Entering edit mode
@markrogers-8451
Last seen 6.8 years ago

Aaron, it looks as if your trick to download the R tarball and compile it myself worked.  It's unclear to me why the synaptic package manager yields a different executable; but for now I'm just happy to have it working.  Many thanks for your help!