Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with HTTPS or Subversion.

Download ZIP
Browse files

ENH: add 'return_rank' keyword to pseudo-inverse

  • Loading branch information...
commit 6cc12fe08595435af1da6ef7173ac1d32eb00a51 1 parent 1b51c8d
@jakevdp jakevdp authored pv committed
Showing with 37 additions and 7 deletions.
  1. +37 −7 scipy/linalg/basic.py
View
44 scipy/linalg/basic.py
@@ -2,6 +2,7 @@
# Author: Pearu Peterson, March 2002
#
# w/ additions by Travis Oliphant, March 2002
+# and Jake Vanderplas, August 2012
__all__ = ['solve', 'solve_triangular', 'solveh_banded', 'solve_banded',
'inv', 'det', 'lstsq', 'pinv', 'pinv2', 'pinvh']
@@ -502,7 +503,7 @@ def lstsq(a, b, cond=None, overwrite_a=False, overwrite_b=False):
return x, resids, rank, s
-def pinv(a, cond=None, rcond=None):
+def pinv(a, cond=None, rcond=None, return_rank=False):
"""
Compute the (Moore-Penrose) pseudo-inverse of a matrix.
@@ -517,11 +518,15 @@ def pinv(a, cond=None, rcond=None):
Cutoff for 'small' singular values in the least-squares solver.
Singular values smaller than ``rcond * largest_singular_value``
are considered zero.
+ return_rank : bool, optional
+ if True, return the effective rank of the matrix
Returns
-------
B : array, shape (N, M)
The pseudo-inverse of matrix `a`.
+ rank : int
+ The effective rank of the matrix. Returned if return_rank == True
Raises
------
@@ -542,10 +547,16 @@ def pinv(a, cond=None, rcond=None):
b = np.identity(a.shape[0], dtype=a.dtype)
if rcond is not None:
cond = rcond
- return lstsq(a, b, cond=cond)[0]
+
+ x, resids, rank, s = lstsq(a, b, cond=cond)
+
+ if return_rank:
+ return x, rank
+ else:
+ return x
-def pinv2(a, cond=None, rcond=None):
+def pinv2(a, cond=None, rcond=None, return_rank=False):
"""
Compute the (Moore-Penrose) pseudo-inverse of a matrix.
@@ -562,11 +573,15 @@ def pinv2(a, cond=None, rcond=None):
Singular values smaller than ``rcond*largest_singular_value``
are considered zero.
If None or -1, suitable machine precision is used.
+ return_rank : bool, optional
+ if True, return the effective rank of the matrix
Returns
-------
B : array, shape (N, M)
The pseudo-inverse of matrix `a`.
+ rank : int
+ The effective rank of the matrix. Returned if return_rank == True
Raises
------
@@ -597,10 +612,15 @@ def pinv2(a, cond=None, rcond=None):
psigma_diag = np.zeros_like(s)
psigma_diag[above_cutoff] = 1.0 / s[above_cutoff]
- return np.transpose(np.conjugate(np.dot(u * psigma_diag, vh)))
+ B = np.transpose(np.conjugate(np.dot(u * psigma_diag, vh)))
+
+ if return_rank:
+ return B, np.sum(above_cutoff)
+ else:
+ return B
-def pinvh(a, cond=None, rcond=None, lower=True):
+def pinvh(a, cond=None, rcond=None, lower=True, return_rank=False):
"""Compute the (Moore-Penrose) pseudo-inverse of a hermetian matrix.
Calculate a generalized inverse of a symmetric matrix using its
@@ -619,10 +639,15 @@ def pinvh(a, cond=None, rcond=None, lower=True):
lower : boolean
Whether the pertinent array data is taken from the lower or upper
triangle of a. (Default: lower)
+ return_rank : bool, optional
+ if True, return the effective rank of the matrix
Returns
-------
B : array, shape (N, N)
+ The pseudo-inverse of matrix `a`.
+ rank : int
+ The effective rank of the matrix. Returned if return_rank == True
Raises
------
@@ -651,9 +676,14 @@ def pinvh(a, cond=None, rcond=None, lower=True):
factor = {'f': 1E3, 'd': 1E6}
cond = factor[t] * np.finfo(t).eps
- # unlike svd case, eigh can lead to negative eigenvalues
+ # For hermitian matrices, singular values equal abs(eigenvalues)
above_cutoff = (abs(s) > cond * np.max(abs(s)))
psigma_diag = np.zeros_like(s)
psigma_diag[above_cutoff] = 1.0 / s[above_cutoff]
- return np.dot(u * psigma_diag, np.conjugate(u).T)
+ B = np.dot(u * psigma_diag, np.conjugate(u).T)
+
+ if return_rank:
+ return B, np.sum(above_cutoff)
+ else:
+ return B
Please sign in to comment.
Something went wrong with that request. Please try again.