From e680e5f9c570211a3e631569205abd0f0914ec1b Mon Sep 17 00:00:00 2001 From: Lars Buitinck Date: Tue, 26 Feb 2013 18:31:23 +0100 Subject: [PATCH 1/2] ENH: sparse: min and max methods on matrices with .data, except dia axis keyword not supported (yet). --- scipy/sparse/bsr.py | 4 +-- scipy/sparse/compressed.py | 4 +-- scipy/sparse/coo.py | 4 +-- scipy/sparse/data.py | 43 +++++++++++++++++++++++++++++++++ scipy/sparse/tests/test_base.py | 42 ++++++++++++++++++++++++++++++++ 5 files changed, 91 insertions(+), 6 deletions(-) diff --git a/scipy/sparse/bsr.py b/scipy/sparse/bsr.py index 6529b40efb14..d06103d7b6e7 100644 --- a/scipy/sparse/bsr.py +++ b/scipy/sparse/bsr.py @@ -10,7 +10,7 @@ import numpy as np -from .data import _data_matrix +from .data import _data_matrix, _minmax_mixin from .compressed import _cs_matrix from .base import isspmatrix, _formats from .sputils import isshape, getdtype, to_native, upcast @@ -18,7 +18,7 @@ from .sparsetools import bsr_matvec, bsr_matvecs, csr_matmat_pass1, \ bsr_matmat_pass2, bsr_transpose, bsr_sort_indices -class bsr_matrix(_cs_matrix): +class bsr_matrix(_cs_matrix, _minmax_mixin): """Block Sparse Row matrix This can be instantiated in several ways: diff --git a/scipy/sparse/compressed.py b/scipy/sparse/compressed.py index d750d79cf6d5..2a6f52aa4f61 100644 --- a/scipy/sparse/compressed.py +++ b/scipy/sparse/compressed.py @@ -10,13 +10,13 @@ from scipy.lib.six.moves import xrange from .base import spmatrix, isspmatrix, SparseEfficiencyWarning -from .data import _data_matrix +from .data import _data_matrix, _minmax_mixin from . import sparsetools from .sputils import upcast, upcast_char, to_native, isdense, isshape, \ getdtype, isscalarlike, isintlike -class _cs_matrix(_data_matrix): +class _cs_matrix(_data_matrix, _minmax_mixin): """base matrix class for compressed row and column oriented matrices""" def __init__(self, arg1, shape=None, dtype=None, copy=False): diff --git a/scipy/sparse/coo.py b/scipy/sparse/coo.py index e1afdeee400e..c2f8ae7eb089 100644 --- a/scipy/sparse/coo.py +++ b/scipy/sparse/coo.py @@ -13,10 +13,10 @@ from .sparsetools import coo_tocsr, coo_todense, coo_matvec from .base import isspmatrix -from .data import _data_matrix +from .data import _data_matrix, _minmax_mixin from .sputils import upcast, upcast_char, to_native, isshape, getdtype, isintlike -class coo_matrix(_data_matrix): +class coo_matrix(_data_matrix, _minmax_mixin): """ A sparse matrix in COOrdinate format. diff --git a/scipy/sparse/data.py b/scipy/sparse/data.py index 7f72f208e14b..76162e647548 100644 --- a/scipy/sparse/data.py +++ b/scipy/sparse/data.py @@ -92,3 +92,46 @@ def method(self): return method setattr(_data_matrix, name, _create_method(npfunc)) + + +class _minmax_mixin(object): + """Mixin for min and max methods. + + These are not implemented for dia_matrix, hence the separate class. + """ + + def max(self): + """Maximum of the elements of this matrix. + + This takes all elements into account, not just the non-zero ones. + + Returns + ------- + amax : self.dtype + Maximum element. + """ + zero = self.dtype.type(0) + if self.nnz == 0: + return zero + mx = np.max(self.data) + if self.nnz != np.product(self.shape): + mx = max(zero, mx) + return mx + + def min(self): + """Minimum of the elements of this matrix. + + This takes all elements into account, not just the non-zero ones. + + Returns + ------- + amin : self.dtype + Minimum element. + """ + zero = self.dtype.type(0) + if self.nnz == 0: + return zero + mn = np.min(self.data) + if self.nnz != np.product(self.shape): + mn = min(zero, mn) + return mn diff --git a/scipy/sparse/tests/test_base.py b/scipy/sparse/tests/test_base.py index efe8d8f9956a..c1acb75abebb 100644 --- a/scipy/sparse/tests/test_base.py +++ b/scipy/sparse/tests/test_base.py @@ -1575,6 +1575,26 @@ def test_eliminate_zeros(self): assert_array_equal(asp.data,[1, 2, 3]) assert_array_equal(asp.todense(),bsp.todense()) + def test_minmax(self): + for dtype in [np.float32, np.float64, np.int32, np.int64]: + D = np.arange(20, dtype=dtype).reshape(5,4) + + X = csr_matrix(D) + assert_equal(X.min(), 0) + assert_equal(X.max(), 19) + assert_equal(X.min().dtype, dtype) + assert_equal(X.max().dtype, dtype) + + D *= -1 + X = csr_matrix(D) + assert_equal(X.min(), -19) + assert_equal(X.max(), 0) + + D += 5 + X = csr_matrix(D) + assert_equal(X.min(), -14) + assert_equal(X.max(), 5) + def test_ufuncs(self): X = csr_matrix(np.arange(20).reshape(4, 5) / 20.) for f in ["sin", "tan", "arcsin", "arctan", "sinh", "tanh", @@ -1685,6 +1705,21 @@ def test_sort_indices(self): assert_array_equal(asp.indices,[1, 2, 7, 4, 5]) assert_array_equal(asp.todense(),bsp.todense()) + def test_minmax(self): + # try a fully dense matrix + X = csc_matrix(np.arange(1, 10).reshape(3, 3)) + assert_equal(X.min(), 1) + assert_equal(X.min().dtype, X.dtype) + + X = -X + assert_equal(X.max(), -1) + + # and a fully sparse matrix + Z = csc_matrix(np.zeros(1)) + assert_equal(Z.min(), 0) + assert_equal(Z.max(), 0) + assert_equal(Z.max().dtype, Z.dtype) + def test_ufuncs(self): X = csc_matrix(np.arange(21).reshape(7, 3) / 21.) for f in ["sin", "tan", "arcsin", "arctan", "sinh", "tanh", @@ -2110,6 +2145,13 @@ def test_eliminate_zeros(self): assert_array_equal(asp.nnz, 3*4) assert_array_equal(asp.todense(),bsp.todense()) + def test_bsr(self): + D = np.arange(20, dtype=float).reshape(5,4) + D[0:2, :] = 0 + X = bsr_matrix(D) + assert_equal(X.min(), 0) + assert_equal(X.max(), 19) + def test_bsr_matvec(self): A = bsr_matrix( arange(2*3*4*5).reshape(2*4,3*5), blocksize=(4,5) ) x = arange(A.shape[1]).reshape(-1,1) From 3807b6dc879fd16527386008cb47c3f860e52f16 Mon Sep 17 00:00:00 2001 From: Pauli Virtanen Date: Mon, 11 Mar 2013 22:42:53 +0200 Subject: [PATCH 2/2] TST: sparse: refactor minmax test --- scipy/sparse/tests/test_base.py | 97 +++++++++++++++++---------------- 1 file changed, 51 insertions(+), 46 deletions(-) diff --git a/scipy/sparse/tests/test_base.py b/scipy/sparse/tests/test_base.py index c1acb75abebb..54e8a662cbe0 100644 --- a/scipy/sparse/tests/test_base.py +++ b/scipy/sparse/tests/test_base.py @@ -1410,6 +1410,49 @@ def test_mu(self): assert_equal(S1.dtype,D1.dtype) +class _TestMinMax(object): + def test_minmax(self): + for dtype in [np.float32, np.float64, np.int32, np.int64]: + D = np.arange(20, dtype=dtype).reshape(5,4) + + X = self.spmatrix(D) + assert_equal(X.min(), 0) + assert_equal(X.max(), 19) + assert_equal(X.min().dtype, dtype) + assert_equal(X.max().dtype, dtype) + + D *= -1 + X = self.spmatrix(D) + assert_equal(X.min(), -19) + assert_equal(X.max(), 0) + + D += 5 + X = self.spmatrix(D) + assert_equal(X.min(), -14) + assert_equal(X.max(), 5) + + # try a fully dense matrix + X = self.spmatrix(np.arange(1, 10).reshape(3, 3)) + assert_equal(X.min(), 1) + assert_equal(X.min().dtype, X.dtype) + + X = -X + assert_equal(X.max(), -1) + + # and a fully sparse matrix + Z = self.spmatrix(np.zeros(1)) + assert_equal(Z.min(), 0) + assert_equal(Z.max(), 0) + assert_equal(Z.max().dtype, Z.dtype) + + # another test + D = np.arange(20, dtype=float).reshape(5,4) + D[0:2, :] = 0 + X = self.spmatrix(D) + assert_equal(X.min(), 0) + assert_equal(X.max(), 19) + + #------------------------------------------------------------------------------ # Tailored base class for generic tests #------------------------------------------------------------------------------ @@ -1443,7 +1486,8 @@ def wrapper(*a, **kw): def sparse_test_class(getset=True, slicing=True, slicing_assign=True, fancy_indexing=True, fancy_assign=True, - fancy_multidim_indexing=True, fancy_multidim_assign=True): + fancy_multidim_indexing=True, fancy_multidim_assign=True, + minmax=True): """ Construct a base class, optionally converting some of the tests in the suite to check that the feature is not implemented. @@ -1462,6 +1506,7 @@ def sparse_test_class(getset=True, slicing=True, slicing_assign=True, fancy_indexing and fancy_multidim_indexing), _possibly_unimplemented(_TestFancyMultidimAssign, fancy_multidim_assign and fancy_assign), + _possibly_unimplemented(_TestMinMax, minmax), TestCase) # check that test names do not clash @@ -1575,26 +1620,6 @@ def test_eliminate_zeros(self): assert_array_equal(asp.data,[1, 2, 3]) assert_array_equal(asp.todense(),bsp.todense()) - def test_minmax(self): - for dtype in [np.float32, np.float64, np.int32, np.int64]: - D = np.arange(20, dtype=dtype).reshape(5,4) - - X = csr_matrix(D) - assert_equal(X.min(), 0) - assert_equal(X.max(), 19) - assert_equal(X.min().dtype, dtype) - assert_equal(X.max().dtype, dtype) - - D *= -1 - X = csr_matrix(D) - assert_equal(X.min(), -19) - assert_equal(X.max(), 0) - - D += 5 - X = csr_matrix(D) - assert_equal(X.min(), -14) - assert_equal(X.max(), 5) - def test_ufuncs(self): X = csr_matrix(np.arange(20).reshape(4, 5) / 20.) for f in ["sin", "tan", "arcsin", "arctan", "sinh", "tanh", @@ -1705,21 +1730,6 @@ def test_sort_indices(self): assert_array_equal(asp.indices,[1, 2, 7, 4, 5]) assert_array_equal(asp.todense(),bsp.todense()) - def test_minmax(self): - # try a fully dense matrix - X = csc_matrix(np.arange(1, 10).reshape(3, 3)) - assert_equal(X.min(), 1) - assert_equal(X.min().dtype, X.dtype) - - X = -X - assert_equal(X.max(), -1) - - # and a fully sparse matrix - Z = csc_matrix(np.zeros(1)) - assert_equal(Z.min(), 0) - assert_equal(Z.max(), 0) - assert_equal(Z.max().dtype, Z.dtype) - def test_ufuncs(self): X = csc_matrix(np.arange(21).reshape(7, 3) / 21.) for f in ["sin", "tan", "arcsin", "arctan", "sinh", "tanh", @@ -1773,7 +1783,8 @@ def test_slicing_3(self): class TestDOK(sparse_test_class(slicing=False, slicing_assign=False, fancy_indexing=False, - fancy_assign=False)): + fancy_assign=False, + minmax=False)): spmatrix = dok_matrix def test_mult(self): @@ -1908,7 +1919,7 @@ def test_fancy_assign_slice(self): def test_fancy_indexing_multidim_set(self): pass -class TestLIL(sparse_test_class()): +class TestLIL(sparse_test_class(minmax=False)): spmatrix = lil_matrix def test_dot(self): @@ -2069,7 +2080,8 @@ def test_constructor4(self): class TestDIA(sparse_test_class(getset=False, slicing=False, slicing_assign=False, - fancy_indexing=False, fancy_assign=False)): + fancy_indexing=False, fancy_assign=False, + minmax=False)): spmatrix = dia_matrix def test_constructor1(self): @@ -2145,13 +2157,6 @@ def test_eliminate_zeros(self): assert_array_equal(asp.nnz, 3*4) assert_array_equal(asp.todense(),bsp.todense()) - def test_bsr(self): - D = np.arange(20, dtype=float).reshape(5,4) - D[0:2, :] = 0 - X = bsr_matrix(D) - assert_equal(X.min(), 0) - assert_equal(X.max(), 19) - def test_bsr_matvec(self): A = bsr_matrix( arange(2*3*4*5).reshape(2*4,3*5), blocksize=(4,5) ) x = arange(A.shape[1]).reshape(-1,1)