pcarbo / varbvs

Large-scale Bayesian variable selection for R and MATLAB.
http://pcarbo.github.io/varbvs
GNU General Public License v3.0
42 stars 14 forks source link

An edge case with `varbvsmix` #15

Closed gaow closed 6 years ago

gaow commented 6 years ago

Click here for the data.

To reproduce:

> attach(readRDS('test20170816.rds'))    
> varbvs::varbvsmix(X, Z, y, sa = c(0,mixsd^2))
Error in solve.default(crossprod(Z), c(y %*% Z)) : 
  system is computationally singular: reciprocal condition number = 3.53577e-17

I'm sure it is problem of input data, but would be nice if varbvs can handle it more robustly. Thanks a lot!

pcarbo commented 6 years ago

@gaow This isn't so much a bug as a problem with your covariate matrix "Z", but I've handled it a bit more robustly in the new version of varbvs; please follow the instructions in the README to install varbvs 2.3-4 from the master branch.

Although varbvsmix will now work, crossprod(Z) is close to singular because several of your covariates are highly correlated; if it is realistic to do so, you should consider ways to adjust Z so that you don't have eigenvalues near zero. Perhaps you could apply the QR decomposition to find an independent basis of some of the columns.

library(varbvs)
attach(readRDS("test20170816.rds"))
R <- cor(Z)
v <- eigen(R)$values
image(Matrix(abs(R)))
plot(1:20,log10(v),pch = 20)

I'll consider this "solved", but let me know if you need some help with coming up with a better matrix Z.

pcarbo commented 6 years ago

Also, I should, in principle, run better checks of input matrix Z; e.g., compute eigenvalues and returning a warning when crossprod(Z) is close to singular. But for now I will consider this issue "closed".