How to access lower triangle of NMM numpy array

Question

I have a numpy array of the shape arr.shape = N,M,M.

I want to access the lower triangles for each M,M array. I tried using

arr1 = arr[:,np.tril_indices(M,-1)]
arr1 = arr[:][np.tril_indices(M,-1)]

etc, with the kernel dying in the first case, while in the second case I get an error saying that:

   ---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-23-1b36c5b12706> in <module>
----> 1 arr1 = arr[:][np.tril_indices(M,-1)]

IndexError: index 6 is out of bounds for axis 0 with size 6

Where

N=6

To clarify I want to find all the elements in the lower triangle of each M,M array(N such instances) and save the result in a new array of the shape:

arr1.shape = (N,(M*(M-1))/2)

Edit:

While np.tril(arr) works, it results in an array

arr1 = np.tril(arr)
arr1.shape

#(N,M,M)

I want the resulting array to be of the specified shape, i.e. I dont want the upper parts of the arrays

Thank you

Please show what you tried and how it failed, i.e. be specific. Where does in not behave as expected? What does already work? As a new user, also take the [tour] and read [ask]. — Ulrich Eckhardt, Nov 10 '21 at 11:28
Alright, thank you for your feedback, I shall be more specific next time. — Aayush Desai, Nov 10 '21 at 11:54

Mortz · Answer 1 · 2021-11-10T12:36:50.063

import numpy as np 
a = np.random.rand(2, 5, 5) 
#array([[[0.28212197, 0.29827562, 0.05151153, 0.90448236, 0.07521404],
#        [0.38938978, 0.67007919, 0.83561652, 0.5950061 , 0.73563179],
#        [0.77515285, 0.31973392, 0.91861436, 0.87386527, 0.85917542],
#        [0.12588184, 0.09173029, 0.28577701, 0.4884228 , 0.07183555],
#        [0.68656271, 0.19941039, 0.07924489, 0.15046004, 0.91011737]],
#
#       [[0.18662788, 0.45745028, 0.14557573, 0.22425571, 0.14204739],
#        [0.44502694, 0.85773626, 0.78554919, 0.07306402, 0.14608384],
#        [0.70620254, 0.81497515, 0.09397011, 0.32053184, 0.255485  ],
#        [0.50139688, 0.51539848, 0.24719375, 0.80708819, 0.39685176],
#        [0.94052069, 0.53927081, 0.39567362, 0.06065674, 0.53479994]]])

np.tril(a) 
#array([[[0.28212197, 0.        , 0.        , 0.        , 0.        ],
#        [0.38938978, 0.67007919, 0.        , 0.        , 0.        ],
#        [0.77515285, 0.31973392, 0.91861436, 0.        , 0.        ],
#        [0.12588184, 0.09173029, 0.28577701, 0.4884228 , 0.        ],
#        [0.68656271, 0.19941039, 0.07924489, 0.15046004, 0.91011737]],
#
#       [[0.18662788, 0.        , 0.        , 0.        , 0.        ],
#        [0.44502694, 0.85773626, 0.        , 0.        , 0.        ],
#        [0.70620254, 0.81497515, 0.09397011, 0.        , 0.        ],
#        [0.50139688, 0.51539848, 0.24719375, 0.80708819, 0.        ],
#        [0.94052069, 0.53927081, 0.39567362, 0.06065674, 0.53479994]]])

If you want to remove the zeros and flatten it to a (2, 15) array (note that there are 10 zeros in each lower triangle array) -

a_no_zeros = np.array([el
for mat in a_lower
for row in mat
for el in row
if el > 0
]).reshape(2, 15)

I think @hpaulj comment was for my original answer which had `np.tril(a[:,])` instead of `np.tril(a)` — Mortz, Nov 10 '21 at 11:49
I want the resulting array to be of the shape (2,10) i.e. I dont want the upper half of the array. Is there a way to do that? — Aayush Desai, Nov 10 '21 at 11:59
Ah okay, I was wondering if there was a method without for loop, with a numpy function or something cause I am working on an array of the dimensions (6,9404,9404) — Aayush Desai, Nov 10 '21 at 13:36
Just to point out, a_no_zeros = a[a!=0].reshape(2,10) works as well. — Aayush Desai, Nov 10 '21 at 13:49

hpaulj · Accepted Answer · 2021-11-10T16:16:51.713

When working with the tri... set of functions it can be useful to examine the source code. They are all python, and based on np.tri.

Make a small sample array - to illustrate and verify the answer:

In [205]: arr = np.arange(18).reshape(2,3,3)  # arange(1,19) might be better
In [206]: arr
Out[206]: 
array([[[ 0,  1,  2],
        [ 3,  4,  5],
        [ 6,  7,  8]],

       [[ 9, 10, 11],
        [12, 13, 14],
        [15, 16, 17]]])

tril sets the upper triangle values to 0. It works in this case, but application to 3d arrays is not documented.

In [207]: np.tril(arr) 
Out[207]: 
array([[[ 0,  0,  0],
        [ 3,  4,  0],
        [ 6,  7,  8]],

       [[ 9,  0,  0],
        [12, 13,  0],
        [15, 16, 17]]])

But in the code if first constructs a boolean mask from the last 2 dimensions:

In [208]: mask = np.tri(*arr.shape[-2:], dtype=bool)
In [209]: mask
Out[209]: 
array([[ True, False, False],
       [ True,  True, False],
       [ True,  True,  True]])

and uses np.where to set some values to 0. This works in the 3d case by broadcasting. mask and arr match on the last 2 dimensions, so mask can broadcast to match:

In [210]: np.where(mask, arr, 0)
Out[210]: 
array([[[ 0,  0,  0],
        [ 3,  4,  0],
        [ 6,  7,  8]],

       [[ 9,  0,  0],
        [12, 13,  0],
        [15, 16, 17]]])

Your tril_indices is just the indices of this mask:

In [217]: np.nonzero(mask)    # aka np.where
Out[217]: (array([0, 1, 1, 2, 2, 2]), array([0, 0, 1, 0, 1, 2]))
In [218]: np.tril_indices(3)
Out[218]: (array([0, 1, 1, 2, 2, 2]), array([0, 0, 1, 0, 1, 2]))

They can't be used directly to index arr:

In [220]: arr[np.tril_indices(3)].shape
Traceback (most recent call last):
  File "<ipython-input-220-e26dc1f514cc>", line 1, in <module>
    arr[np.tril_indices(3)].shape
IndexError: index 2 is out of bounds for axis 0 with size 2

In [221]: arr[:,np.tril_indices(3)].shape
Out[221]: (2, 2, 6, 3)

But unpacking the two indexing arrays:

In [222]: I,J = np.tril_indices(3)
In [223]: I,J
Out[223]: (array([0, 1, 1, 2, 2, 2]), array([0, 0, 1, 0, 1, 2]))
In [224]: arr[:,I,J]
Out[224]: 
array([[ 0,  3,  4,  6,  7,  8],
       [ 9, 12, 13, 15, 16, 17]])

The boolean mask can also be used directly:

In [226]: arr[:,mask]
Out[226]: 
array([[ 0,  3,  4,  6,  7,  8],
       [ 9, 12, 13, 15, 16, 17]])

The base np.tri works by simply doing an outer >= on indices

In [231]: m = np.greater_equal.outer(np.arange(3),np.arange(3))
In [232]: m
Out[232]: 
array([[ True, False, False],
       [ True,  True, False],
       [ True,  True,  True]])
In [234]: np.arange(3)[:,None]>=np.arange(3)
Out[234]: 
array([[ True, False, False],
       [ True,  True, False],
       [ True,  True,  True]])

How to access lower triangle of N*M*M numpy array

2 Answers2

How to access lower triangle of NMM numpy array