Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How to calculate sum of abs of all off-diagonal elements of a numpy array?

I'm new to numpy and I want to calculate the sum of abs of all off-diagonal elements of a numpy array. off-diagonal elements of a matrix are all elements of a matrix except ones that are in the main diagonal of it.

I want to calculate the sum of abs of them so I can implement the Jacobi eigenvalue algorithm of this lecture

so, to compute it, I think this code will work:

import numpy as np
off_diagonal_sum = 0
for i in range(n): # n is the dimension of our square matrix
    # mat is our matrix
    off_diagonal_sum = off_diagonal_sum + np.sum(np.abs(mat[i, (i + 1):n]))
    off_diagonal_sum = off_diagonal_sum + np.sum(np.abs(mat[i, 0:(i - 1)]))

but as I'm new to numpy, I think there should be a simpler and shorter way to compute that. do you have any idea?

thanks in advance.

like image 747
Peyman Avatar asked Jan 28 '23 00:01

Peyman


2 Answers

There is a nice option in NumPy called diag_indices which returns you the indices of diagonal elements of a 2-d array. Using this, you can get the sum of diagonal elements and subtract it from the sum of complete array np.sum(arr) to get the sum of off diagonal elements without using any explicit for-loops. To get the absolute sum, you simply use np.abs to get absolute values of each element and then perform the task as follows.

Example

import numpy as np
arr = np.random.randint(1, 20, size=(3,3)) # Define a random 2d array
print (arr)
dia = np.diag_indices(3) # indices of diagonal elements
dia_sum = sum(arr[dia]) # sum of diagonal elements
off_dia_sum = np.sum(arr) - dia_sum # subtract the diagonal sum from total array sum
print (off_dia_sum)

[[12 19 10]
 [ 3 13 18]
 [16 16  6]]
82

Alternate 1

You can also use simply np.trace to get the sum of diagonal elements and then subtract it from the total array sum to get the sum of the off-diagonal elements.

off_dia_sum = np.sum(arr) - np.trace(arr)

Alternate 2

Using np.diagonal to get the diagonal elements and then take the sum and subtract from the total sum as

dia_sum = sum(np.diagonal(arr))
off_dia_sum = np.sum(arr) - dia_sum

Alternate 3

Using list comprehension you can do the following where you only store the elements in the list if it is off diagonal, which means if both indices i and j are not equal.

size = len(arr)
off_dia_sum = sum([arr[i][j] for i in range(size) for j in range(size) if i!=j])
like image 161
Sheldore Avatar answered Jan 29 '23 13:01

Sheldore


Subtracting the trace from the full sum is fast and most of the time perfectly fine.

If, however, off-diagonal elements are much (many orders of magnitude) smaller than diagonal elements, results can become imprecise.

Therefore here are two numpy solutions that are not quite as fast but safer.

# make an example with very small off-diag elements

from scipy.spatial.distance import squareform

N = 100
A = np.identity(N) + 1e-8*np.random.random((N, N))
A = np.abs(A+A.T) # A is symmetric, also we take the absolute value right here

# method 1
np.bincount(np.identity(N, int).ravel(), A.ravel())[0]
# 9.934386601640464e-05

# method 2 - this only works with symmetric A:
2*squareform(A, checks=False).sum()
# 9.934386601640431e-05

# and the "subtract trace" method
A.sum()-A.trace()
# 9.934386602594714e-05

# with double precision all looks ok in this example
# so for the sake of demonstration let's go to single precision
A = A.astype(np.single)

np.bincount(np.identity(N, int).ravel(), A.ravel())[0]
# 9.93438660777668e-05
2*squareform(A, checks=False).sum()
# 9.93438734440133e-05

# so far so good, but ...

A.sum()-A.trace()
# 7.6293945e-05

# Ouch!
like image 33
Paul Panzer Avatar answered Jan 29 '23 14:01

Paul Panzer