Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1"""Sparse matrix norms. 

2 

3""" 

4import numpy as np 

5from scipy.sparse import issparse 

6 

7from numpy import Inf, sqrt, abs 

8 

9__all__ = ['norm'] 

10 

11 

12def _sparse_frobenius_norm(x): 

13 if np.issubdtype(x.dtype, np.complexfloating): 

14 sqnorm = abs(x).power(2).sum() 

15 else: 

16 sqnorm = x.power(2).sum() 

17 return sqrt(sqnorm) 

18 

19 

20def norm(x, ord=None, axis=None): 

21 """ 

22 Norm of a sparse matrix 

23 

24 This function is able to return one of seven different matrix norms, 

25 depending on the value of the ``ord`` parameter. 

26 

27 Parameters 

28 ---------- 

29 x : a sparse matrix 

30 Input sparse matrix. 

31 ord : {non-zero int, inf, -inf, 'fro'}, optional 

32 Order of the norm (see table under ``Notes``). inf means numpy's 

33 `inf` object. 

34 axis : {int, 2-tuple of ints, None}, optional 

35 If `axis` is an integer, it specifies the axis of `x` along which to 

36 compute the vector norms. If `axis` is a 2-tuple, it specifies the 

37 axes that hold 2-D matrices, and the matrix norms of these matrices 

38 are computed. If `axis` is None then either a vector norm (when `x` 

39 is 1-D) or a matrix norm (when `x` is 2-D) is returned. 

40 

41 Returns 

42 ------- 

43 n : float or ndarray 

44 

45 Notes 

46 ----- 

47 Some of the ord are not implemented because some associated functions like, 

48 _multi_svd_norm, are not yet available for sparse matrix. 

49 

50 This docstring is modified based on numpy.linalg.norm. 

51 https://github.com/numpy/numpy/blob/master/numpy/linalg/linalg.py 

52 

53 The following norms can be calculated: 

54 

55 ===== ============================ 

56 ord norm for sparse matrices 

57 ===== ============================ 

58 None Frobenius norm 

59 'fro' Frobenius norm 

60 inf max(sum(abs(x), axis=1)) 

61 -inf min(sum(abs(x), axis=1)) 

62 0 abs(x).sum(axis=axis) 

63 1 max(sum(abs(x), axis=0)) 

64 -1 min(sum(abs(x), axis=0)) 

65 2 Not implemented 

66 -2 Not implemented 

67 other Not implemented 

68 ===== ============================ 

69 

70 The Frobenius norm is given by [1]_: 

71 

72 :math:`||A||_F = [\\sum_{i,j} abs(a_{i,j})^2]^{1/2}` 

73 

74 References 

75 ---------- 

76 .. [1] G. H. Golub and C. F. Van Loan, *Matrix Computations*, 

77 Baltimore, MD, Johns Hopkins University Press, 1985, pg. 15 

78 

79 Examples 

80 -------- 

81 >>> from scipy.sparse import * 

82 >>> import numpy as np 

83 >>> from scipy.sparse.linalg import norm 

84 >>> a = np.arange(9) - 4 

85 >>> a 

86 array([-4, -3, -2, -1, 0, 1, 2, 3, 4]) 

87 >>> b = a.reshape((3, 3)) 

88 >>> b 

89 array([[-4, -3, -2], 

90 [-1, 0, 1], 

91 [ 2, 3, 4]]) 

92 

93 >>> b = csr_matrix(b) 

94 >>> norm(b) 

95 7.745966692414834 

96 >>> norm(b, 'fro') 

97 7.745966692414834 

98 >>> norm(b, np.inf) 

99 9 

100 >>> norm(b, -np.inf) 

101 2 

102 >>> norm(b, 1) 

103 7 

104 >>> norm(b, -1) 

105 6 

106 

107 """ 

108 if not issparse(x): 

109 raise TypeError("input is not sparse. use numpy.linalg.norm") 

110 

111 # Check the default case first and handle it immediately. 

112 if axis is None and ord in (None, 'fro', 'f'): 

113 return _sparse_frobenius_norm(x) 

114 

115 # Some norms require functions that are not implemented for all types. 

116 x = x.tocsr() 

117 

118 if axis is None: 

119 axis = (0, 1) 

120 elif not isinstance(axis, tuple): 

121 msg = "'axis' must be None, an integer or a tuple of integers" 

122 try: 

123 int_axis = int(axis) 

124 except TypeError: 

125 raise TypeError(msg) 

126 if axis != int_axis: 

127 raise TypeError(msg) 

128 axis = (int_axis,) 

129 

130 nd = 2 

131 if len(axis) == 2: 

132 row_axis, col_axis = axis 

133 if not (-nd <= row_axis < nd and -nd <= col_axis < nd): 

134 raise ValueError('Invalid axis %r for an array with shape %r' % 

135 (axis, x.shape)) 

136 if row_axis % nd == col_axis % nd: 

137 raise ValueError('Duplicate axes given.') 

138 if ord == 2: 

139 raise NotImplementedError 

140 #return _multi_svd_norm(x, row_axis, col_axis, amax) 

141 elif ord == -2: 

142 raise NotImplementedError 

143 #return _multi_svd_norm(x, row_axis, col_axis, amin) 

144 elif ord == 1: 

145 return abs(x).sum(axis=row_axis).max(axis=col_axis)[0,0] 

146 elif ord == Inf: 

147 return abs(x).sum(axis=col_axis).max(axis=row_axis)[0,0] 

148 elif ord == -1: 

149 return abs(x).sum(axis=row_axis).min(axis=col_axis)[0,0] 

150 elif ord == -Inf: 

151 return abs(x).sum(axis=col_axis).min(axis=row_axis)[0,0] 

152 elif ord in (None, 'f', 'fro'): 

153 # The axis order does not matter for this norm. 

154 return _sparse_frobenius_norm(x) 

155 else: 

156 raise ValueError("Invalid norm order for matrices.") 

157 elif len(axis) == 1: 

158 a, = axis 

159 if not (-nd <= a < nd): 

160 raise ValueError('Invalid axis %r for an array with shape %r' % 

161 (axis, x.shape)) 

162 if ord == Inf: 

163 M = abs(x).max(axis=a) 

164 elif ord == -Inf: 

165 M = abs(x).min(axis=a) 

166 elif ord == 0: 

167 # Zero norm 

168 M = (x != 0).sum(axis=a) 

169 elif ord == 1: 

170 # special case for speedup 

171 M = abs(x).sum(axis=a) 

172 elif ord in (2, None): 

173 M = sqrt(abs(x).power(2).sum(axis=a)) 

174 else: 

175 try: 

176 ord + 1 

177 except TypeError: 

178 raise ValueError('Invalid norm order for vectors.') 

179 M = np.power(abs(x).power(ord).sum(axis=a), 1 / ord) 

180 return M.A.ravel() 

181 else: 

182 raise ValueError("Improper number of dimensions to norm.")