0

I have one problem while try to computing the 1-norm of a sparse matrix. I am using the function scipy.sparse.linalg.onenormest but it gives me an error because the operator can act only onto square matrix.

Here a code example:

from scipy import sparse

row = array([0,2,2,0,1,2])
col = array([0,0,1,2,2,2])
data = array([1,2,3,4,5,6])

A = sparse.csc_matrix( (data,(row,col)), shape=(5,3) )

onenormest(A)

this is the error:

Traceback (most recent call last):
  File "<ipython console>", line 1, in <module>
  File "C:\Python27\lib\site-packages\scipy\sparse\linalg\_onenormest.py", line 76, in onenormest
    raise ValueError('expected the operator to act like a square matrix')
ValueError: expected the operator to act like a square matrix

The operator onenormest works if I define A as a square matrix, but this is not what I want.

Anyone knows how to calculate the 1-norm of a sparse non-square matrix?

Alessandro Vianello
  • 447
  • 2
  • 6
  • 9

2 Answers2

1

I think that you want numpy.linalg.norm instead;

from numpy import linalg
from scipy import sparse


row = array([0,2,2,0,1,2])
col = array([0,0,1,2,2,2])
data = array([1,2,3,4,5,6])

A = sparse.csc_matrix( (data,(row,col)), shape=(5,3) )

print linalg.norm(A.todense(), ord=1) #15

It does not work to call A.data, since .data of a sparse matrix object is just the data - it appears as a vector instead.

If your sparse matrix is only small, then this is fine. If it is large, then obviously this is a problem. In which case, you can write your own routine.

If you are only interested in the L^1-norm, and casting to dense is not possible, then you could do it via something like this:

def sparseL1Norm = lambda A: max([numpy.abs(A).getcol(i).sum() for i in range(A.shape[1])])
will
  • 10,260
  • 6
  • 46
  • 69
1

This finds the L1-norm of each column:

from scipy import sparse
import numpy as np

row = np.array([0,2,2,0,1,2])
col = np.array([0,0,1,2,2,2])
data = np.array([1,2,3,-4,-5,-6]) # made negative to exercise abs
A = sparse.csc_matrix( (data,(row,col)), shape=(5,3) )
print(abs(A).sum(axis=0))

yields

[[ 3  3 15]]

You could then take the max to find the L1-norm of the matrix:

print(abs(A).sum(axis=0).max())
# 15

abs(A) is a sparse matrix:

In [29]: abs(A)
Out[29]: 
<5x3 sparse matrix of type '<type 'numpy.int64'>'
    with 6 stored elements in Compressed Sparse Column format>

and sum and max are methods of the sparse matrix, so abs(A).sum(axis=0).max() computes the L1-norm without densifying the matrix.

Note: Most NumPy functions (such a np.abs) are not designed to work with sparse matrices. Although np.abs(A) returns the correct result, it arrives there through an indirect route. The more direct route is to use abs(A) which calls A.__abs__(). Thanks to pv. for point this out.

Community
  • 1
  • 1
unutbu
  • 842,883
  • 184
  • 1,785
  • 1,677
  • Thanks! This seems to work. But what I would like to have is the same behavior of the Matlab norm(X,p) function, which gives as a result 15. 15 is the L1 norm of the third column of A, which is the higher L1 norm (first column norm = 4, second column norm = 3, third column norm = 15) – Alessandro Vianello Nov 06 '14 at 14:29
  • This does not work, since `A.data` for a sparse matrix is just `[1 2 3 4 5 6]` – will Nov 06 '14 at 14:36
  • Yes exactly. So the result is then 21. But if for example I define A = sparse.csc_matrix( (data,(row,col)), shape=(5,5) ) and then I calculate onenormest(A) the result is 15, so is working for square sparse matrix, but not with non-square matrix :( – Alessandro Vianello Nov 06 '14 at 14:45
  • `np.abs(A).sum(axis=0).max()` is much neater than my solution :/ – will Nov 06 '14 at 14:49
  • 1
    Do `abs(A).sum(axis=0).max()` instead. `np.abs` goes via some indirect route to call `A.__abs__()`, as you can verify with `pdb.runcall(np.abs, A)`, whereas `abs(A)` calls it directly. In general, Numpy routines and sparse matrices do not work together as of currently. – pv. Nov 06 '14 at 15:46
  • @pv.: Thanks very much for the correction and information. – unutbu Nov 06 '14 at 15:56
  • `A.__abs__().sum(axis=0).max()` makes it clear that you are using the `abs` defined for the sparse matrix. – hpaulj Nov 07 '14 at 21:40
  • @hpaulj: True, but the general rule is to ["never call a magic method directly (but always indirectly through a built-in) unless you know exactly why you need to do that"](http://stackoverflow.com/a/2481631/190597). – unutbu Nov 07 '14 at 21:58
  • In this case there seems to an ambiguity, or at least misunderstandings, as to whether `abs` or `np.abs` is the desired call. – hpaulj Nov 07 '14 at 23:02
  • @hpaulj: I know my opinion on matters of style have shifted and probably will continue to change as time goes on. But currently I believe ideally, code should be written in a way an expert would deem correct but with enough bread crumbs (e.g. in the form of comments) that an intelligent novice would understand. Therefore I'm reluctant to change `abs(A)` to `A.__abs__` since experts like Martelli and pv. probably wouldn't do that. – unutbu Nov 09 '14 at 14:45