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

Submitter:  None
Submitted:  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 Fixed Release:  None
Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

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>
Group 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>
Group 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>
Group 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>
Group 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>
Group 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>
Group 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>
Group 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>
Group Member
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)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group 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 &lt; 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.13-f8d8.
    Corresponding source code