bugGNU Octave - Bugs: bug #44061, U and V result matrices from svds...

 
 

bug #44061: U and V result matrices from svds are not othogonal when rows < 9

Submitted by:  None
Submitted on:  Wed 21 Jan 2015 04:23:36 PM UTC  
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Confirmed Assigned to:  None
Originator Name:  Svend H. Sørensen Originator Email:  -email is unavailable-
Open/Closed:  Open Release:  dev
Operating System:  GNU/Linux

Add a New Comment (Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

( Jump to the original submission)

Sun 18 Dec 2016 08:27:56 PM UTC, comment #10: 

This behavior seems to be still present in Octave 4.2.0.

The script from comment #0 still gives the same results.

Hartmut <hardy>
Fri 06 Feb 2015 04:33:17 PM UTC, comment #9: 

@marco: I put in another cset to make it explicit that tolower is to be used only on char inputs (http://hg.savannah.gnu.org/hgweb/octave/rev/6e2729924601).  I had done my testing on the development branch where tolower has been changed to return numeric input and only operate on char inputs.

Rik <rik5>
Project Administrator
Thu 05 Feb 2015 09:01:09 AM UTC, comment #8: 

Rik,

I think your patch should be modified, since when tmp is a number, tolower (tmp) is not defined (at least, if I try tolower(4) I see a strange symbol).

Coming back to the original bug, I have an explanation (but not yet a fix). If you call A the original matrix and try

[V,s] = eigs([sparse(4,4),A;A',sparse(4,4)],8)

you see four eigenvectors relative to 0, two of them with the last four entries 0. When you ask for 4 singular values of A, close to 4, you get the eigenvectors relative to 2 and two of the eigenvectors relative to 0, precisely the two with the last four entries 0. Try

[V,d]=eigs([sparse(4,4),A;A',sparse(4,4)],4,4)

to see it. This is precisely the command used in svds.
This is the problem, since the last entries are used to build the orthogonal matrix V. Here we get these two eigenvectors because if rows(A)<9, eigs calls eig (see call_eig in eigs.m). If eigs would call ARPACK, since it is an iterative method based on a random initial vector, you would (almost) never see such two vectors, but rather two "random" vectors spanned by the four eingenvectors relative to 0. This is a case in which a random initialization helps you. A second problem is that these two vectors, although orthogonal, are not orthognal when split into two parts. Then, you need to orthogonalize the two parts. I attached a script to demonstrate this second point.

(file #32999)

Marco Caliari <caliari>
Project Member
Thu 05 Feb 2015 02:12:50 AM UTC, comment #7: 

See comment #2, I've already committed a change for the 'LA' error.  You'll need to build from a recent Mercurial tip in order to see the fix.

Rik <rik5>
Project Administrator
Wed 04 Feb 2015 07:54:50 AM UTC, comment #6: 

Andrew,

I recently upgraded to 3.8.2 and I see

V =

   1
   0

s =  1

which is now correct. With 3.8.1 I'm pretty sure that V was wrong, try [V,s]=eigs(rand(2),1). Anyway, the original submission and your test svds(rand(4)) fail as described.

Marco Caliari <caliari>
Project Member
Wed 04 Feb 2015 07:33:43 AM UTC, comment #5: 

Rik, yes Matlab does output orthonormal vectors:

U =

    0.5323   -0.4654   -0.1215    0.6966
    0.4654    0.5323   -0.6966   -0.1215
   -0.4654   -0.5323   -0.6966   -0.1215
   -0.5323    0.4654   -0.1215    0.6966

V =

   -0.9978   -0.0669    0.0000    0.0000
   -0.0669    0.9978    0.0000   -0.0000
   -0.0000    0.0000    0.0002    1.0000
   -0.0000    0.0000    1.0000   -0.0002

Marco, what is wrong with the value returned by eigs?  I get

V =

   1   0

s =  1

Any (non-zero) vector is an eigenvector of the identity matrix, isn't it?  That means that the returned  V  is (the transpose of) an eigenvector corresponding to the eigenvalue returned.

However, for  n<4,  I get a failure of

 [U, S, V] = svds(rand(4))
error: unrecognized value for sigma: LA
error: called from:
error:   /usr/share/octave/3.8.1/m/sparse/eigs.m at line 331, column 9
error:   /usr/share/octave/3.8.1/m/sparse/eigs.m at line 254, column 15
error:   /usr/share/octave/3.8.1/m/sparse/svds.m at line 182, column 17

in version 3.8.1

Lachlan Andrew <lachlan>
Thu 29 Jan 2015 12:53:57 PM UTC, comment #4: 

1)

[U,S,V] = svds(rand(n))

gives orthogonal matrices, both Matlab and Octave (U'*U=1, not U*U').

2) But for this particular matrix, it does not work. It seems a problem of eigs when computing only some of the eigenvalues

[V,s] = eigs(eye(2),1)

V is not the matrix of eigenvectors as described in the help and as it is in Matlab. The dimensions are even wrong (and it is not just the transposed).

Marco Caliari <caliari>
Project Member
Tue 27 Jan 2015 09:34:59 PM UTC, comment #3: 

This may be related to interpretation, and maybe something Octave should be doing better.  If you read the documentation for svds in Octave it does not say that U or V will be orthonormal matrices such that

U*U' = 1
V*V' = 1.

The only thing it states is that

A_approx = U*S*V'

And this does in fact work with the example:

full (A)
ans =

  -1  -1   0   0
  -1   1   0   0
   1  -1   0   0
   1   1   0   0

[U, S, V] = svds (A)

U*S*V' - A
ans =

   0   0   0   0
   0   0   0   0
   0   0   0   0
   0   0   0   0

Can someone with access to Matlab run the example matrix?  According to their documentation they are supposed to output orthonormal U and V matrices.

This doesn't appear to be a problem with this particular matrix.  I tried svds on several matrices created with rand() and U and V are never orthonormal.

Rik <rik5>
Project Administrator
Tue 27 Jan 2015 07:58:47 PM UTC, comment #2: 

@Jordi: See this cset from 3 days ago on stable which fixes the uppercase argument thing (http://hg.savannah.gnu.org/hgweb/octave/rev/08dbf53f34c1).  I haven't pulled it through gui-release and on to dev yet.  I figured there would be something else that would flush that cset down.

Rik <rik5>
Project Administrator
Mon 26 Jan 2015 04:43:42 PM UTC, comment #1: 

Something is definitely fishy here... if nothing else, I see a failure with svds (A), because eigs () does not accept uppercase arguments.

Jordi Gutiérrez Hermoso <jordigh>
Project Administrator
Wed 21 Jan 2015 04:23:36 PM UTC, original submission:  

octave:1> load "test_files/in/two_vars_neg_in.sparse"
octave:2> A = spconvert(two_vars_neg_in)
A =

Compressed Column Sparse (rows = 4, cols = 4, nnz = 8 [50%])

  (1, 1) -> -1
  (2, 1) -> -1
  (3, 1) ->  1
  (4, 1) ->  1
  (1, 2) -> -1
  (2, 2) ->  1
  (3, 2) -> -1
  (4, 2) ->  1

octave:3>  [U,S,V]= svds(A,rows(A),rows(A))
U =

  -0.50000   0.50000   1.00000   0.00000
   0.50000   0.50000   0.00000   1.00000
  -0.50000  -0.50000   0.00000   1.00000
   0.50000  -0.50000   1.00000   0.00000

S =

Diagonal Matrix

   2   0   0   0
   0   2   0   0
   0   0   0   0
   0   0   0   0

V =

  -0.00000  -1.00000   0.00000   0.00000
   1.00000   0.00000   0.00000   0.00000
   0.00000   0.00000   0.00000   0.00000
   0.00000   0.00000   0.00000   0.00000

octave:4> #[U,S,V] = svd (full(A) )
octave:4> ss = diag(S)
ss =

   2
   2
   0
   0

octave:5> nnull=sum(ss(:)==0)
nnull =  2
octave:6> V'*V  ## this should result en the identity matrix but does not !!
ans =

   1.00000   0.00000   0.00000   0.00000
   0.00000   1.00000   0.00000   0.00000
   0.00000   0.00000   0.00000   0.00000
   0.00000   0.00000   0.00000   0.00000



octave:7> U'*U  ## here again we expect to see the identity matrix
ans =

   1.00000   0.00000   0.00000   0.00000
   0.00000   1.00000   0.00000   0.00000
   0.00000   0.00000   2.00000   0.00000
   0.00000   0.00000   0.00000   2.00000

Anonymous

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

Attached Files
file #32999:  script.m added by caliari (370B - text/x-objcsrc)
file #32889:  two_vars_neg_in.sparse added by None (121B - application/octet-stream)
file #32890:  svd-2-var.oct added by None (178B - application/octet-stream)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by hardy (Posted a comment)
  • -email is unavailable- added by lachlan (Posted a comment)
  • -email is unavailable- added by caliari (Posted a comment)
  • -email is unavailable- added by rik5 (Updated the item)
  • -email is unavailable- added by jordigh (Posted a comment)
  • -email is unavailable- added by None (Submitted the item)
  •  

    Do you think this task is very important?
    If so, you can add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

     

     

    Follow 8 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2016-12-30 rik5 SummaryU and V result matrices from svds are not othogonal when M or N < 9 => U and V result matrices from svds are not othogonal when rows < 9
    2016-12-30 rik5 Release3.8.2 => dev
        SummaryU and V result matrices from svds are not othogonal => U and V result matrices from svds are not othogonal when M or N < 9
    2015-02-05 caliari Attached File- => Added script.m, #32999
    2015-01-27 rik5 SummaryU and V result matrises are not othogonal => U and V result matrices from svds are not othogonal
    2015-01-26 jordigh StatusNone => Confirmed
    2015-01-21 None Attached File- => Added two_vars_neg_in.sparse, #32889
        Attached File- => Added svd-2-var.oct, #32890

    Back to the top


    Powered by Savane 3.5