[R] Observations on SVD linpack errors, and a workaround

[Simon Wood]

I get the same problem with this matrix with R 2.6.0 on a 64 bit (Intel Xeon) 
machine running Suse 9, but the svd is fine on a 32 bit windows machine ( R 
2.6.0 again). 

`badx' has another weird property (64 bit linux version)... 
> x <- matrix(runif(prod(dim(badx))),nrow(badx),ncol(badx))
> system.time(qr(x)) ## x same dimensions as badx
   user  system elapsed
  1.209   0.031   1.241
> system.time(qr(badx))
   user  system elapsed
845.896   0.164 846.182
Something similar happens on 32 bit windows...  odd for a non-iterative 
algorithm (qr.R(qr(badx)) has no NA entries, in case you are wondering, and 
has the same singular values as badx, as it should).

best,
Simon

On Wednesday 17 October 2007 07:06, Art Owen wrote:
> Lately I'm getting this error quite a bit:
> Error in La.svd(x, nu, nv) : error code 1 from Lapack routine 'dgesdd'
>
> I'm running R 2.5.0 on a 64 bit Intel machine running Fedora (8 I think).
> Maybe the 64 bit platform is more fragile about declaring convergence.
> I'm seeing way more of these errors than I ever have before.
>
>  From R-Help I see that this issue comes up from time to time.
>
> I'm posting an observation that might help diagnose
> the problem, and a workaround that improves the odds of success.
>
> I have found that sometimes  svd(t(x)) will work when
>
> svd(x) fails.  For example:
>  > load("badx")
>  > svd(badx)$d
>
> Error in La.svd(x, nu, nv) : error code 1 from Lapack routine 'dgesdd'
>
>  > svd(t(badx))$d
>
>   [1] 1.572739e+02 9.614579e+01 7.719867e+01 7.127926e+01 6.490623e+01
>  .... stuff deleted ....
> [126] 8.889272e+00 8.738343e+00 8.447202e+00 8.290393e+00 1.338621e-11
> [131] 1.590829e-12 6.154970e-13
>
> badx was a residual matrix, hence the 3 small singular values.
> I put the output of save(badx,file="badx") on the web if anybody wants to
> play with it.  That matrix is 132 x 10270 entries and the file is
> over 10Mb.  As I write this, it seems to be giving firefox a very bad
> time loading it.  So proceed with caution (if at all) to the file badx
> in the web page stat.stanford.edu/~owen/
> There is also a smaller one, called badx2 which illustrates the much
> rarer case where a skinny matrix makes svd choke, while its wide
> transpose causes no trouble.  Also badx2 did not make firefox hang
> so it might be a safer one to look at.
>
> For now my workaround is to write a wrapper that first tries svd(x).
> If that fails it then tries svd(t(x)).  In about 800 svds the first case
> failed about 100 times.  But the combination never failed.
>
> A simplistic wrapper is listed below.    If SVD failures get very
> common for lots of people then a better solution would be to have
> the svd function itself try both ways.  Another option is to have the
> svd code try the Golub and Reinsch algorithm (or some other SVD)
> on those cases where the Lapack one fails.
>
>
> -Art Owen,  Dept Statistics, Stanford University
>
>
> ##
> ##   Wrapper function for SVD.  If svd(x) fails, try svd( t(x) ).
> ##   If both fail you're out of luck in the SVD department.
> ##   You might succeed by writing a third option based on
> ##   eigen().  That is numerically inferior to svd when the latter
> ##   works.
> ##   -Art Owen, October 2007
> ##
> svdwrapper = function( x, nu, nv, verbose=T ){
> #   Caution: I have not tested this much.
> #   It's here as an example for an R-Help discussion.
>   gotit = F
>   try( {svdx = svd(x,nu,nv); gotit=T}, silent = !verbose )
>   if( gotit )return(svdx)
>   try( {svdtx = svd(t(x),nv,nu); gotit=T}, silent = !verbose )
>   if( !gotit )stop("svd(x) and svd(t(x)) both failed.")
>   if( verbose )print("svd(x) failed but svd(t(x)) worked.")
>   temp    = svdtx$u
>   svdtx$u = svdtx$v
>   svdtx$v = temp
>   svdtx
> }
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html and provide commented, minimal,
> self-contained, reproducible code.

-- 
> Simon Wood, Mathematical Sciences, University of Bath, Bath, BA2 7AY UK
> +44 1225 386603  www.maths.bath.ac.uk/~sw283