Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How to fill upper triangle of numpy array with zeros in place?

Tags:

python

numpy

What is the best way to fill in the lower triangle of a numpy array with zeros in place so that I don't have to do the following:

a=np.random.random((5,5))
a = np.triu(a)

since np.triu returns a copy, not a view. Preferable this would require no list indexing as well since I am working with large arrays.

like image 477
methane Avatar asked May 23 '14 23:05

methane


2 Answers

Digging into the internals of triu you'll find that it just multiplies the input by the output of tri.

So you can just multiply the array in-place by the output of tri:

>>> a = np.random.random((5, 5))
>>> a *= np.tri(*a.shape)
>>> a
array([[ 0.46026582,  0.        ,  0.        ,  0.        ,  0.        ],
       [ 0.76234296,  0.5298908 ,  0.        ,  0.        ,  0.        ],
       [ 0.08797149,  0.14881991,  0.9302515 ,  0.        ,  0.        ],
       [ 0.54794779,  0.36896506,  0.92901552,  0.73747726,  0.        ],
       [ 0.62917827,  0.61674542,  0.44999905,  0.80970863,  0.41860336]])

Like triu, this still creates a second array (the output of tri), but at least it performs the operation itself in-place. The splat is a bit of a shortcut; consider basing your function on the full version of triu for something robust. But note that you can still specify a diagonal:

>>> a = np.random.random((5, 5))
>>> a *= np.tri(*a.shape, k=2)
>>> a
array([[ 0.25473126,  0.70156073,  0.0973933 ,  0.        ,  0.        ],
       [ 0.32859487,  0.58188318,  0.95288351,  0.85735005,  0.        ],
       [ 0.52591784,  0.75030515,  0.82458369,  0.55184033,  0.01341398],
       [ 0.90862183,  0.33983192,  0.46321589,  0.21080121,  0.31641934],
       [ 0.32322392,  0.25091433,  0.03980317,  0.29448128,  0.92288577]])

I now see that the question title and body describe opposite behaviors. Just in case, here's how you can fill the lower triangle with zeros. This requires you to specify the -1 diagonal:

>>> a = np.random.random((5, 5))
>>> a *= 1 - np.tri(*a.shape, k=-1)
>>> a
array([[0.6357091 , 0.33589809, 0.744803  , 0.55254798, 0.38021111],
       [0.        , 0.87316263, 0.98047459, 0.00881754, 0.44115527],
       [0.        , 0.        , 0.51317289, 0.16630385, 0.1470729 ],
       [0.        , 0.        , 0.        , 0.9239731 , 0.11928557],
       [0.        , 0.        , 0.        , 0.        , 0.1840326 ]])
like image 185
senderle Avatar answered Dec 22 '22 00:12

senderle


If speed and memory use are still a limitation and Cython is available, a short Cython function will do what you want. Here's a working version designed for a C-contiguous array with double precision values.

cimport cython
@cython.boundscheck(False)
@cython.wraparound(False)
cpdef make_lower_triangular(double[:,:] A, int k):
    """ Set all the entries of array A that lie above
    diagonal k to 0. """
    cdef int i, j
    for i in range(min(A.shape[0], A.shape[0] - k)):
        for j in range(max(0, i+k+1), A.shape[1]):
            A[i,j] = 0.

This should be significantly faster than any version that involves multiplying by a large temporary array.

like image 32
IanH Avatar answered Dec 21 '22 23:12

IanH