Modifying block matrices in Python

Question

I would like to take a matrix and modify blocks of it. For example, with a 4x4 matrix the {1,2},{1,2} block is to the top left quadrant ([0,1;4,5] below). The {4,1},{4,1} block is the top left quadrant if we rearrange the matrix so the 4th row/column is in position 1 and the 1st in position 2.

Let's made such a 4x4 matrix:

a = np.arange(16).reshape(4, 4)
print(a)

## [[ 0  1  2  3]
##  [ 4  5  6  7]
##  [ 8  9 10 11]
##  [12 13 14 15]]

Now one way of selecting the block, where I specify which rows/columns I want beforehand, is as follows:

C=[3,0]
a[[[C[0],C[0]],[C[1],C[1]]],[[C[0],C[1]],[C[0],C[1]]]]

## array([[15, 12],
##        [ 3,  0]])

Here's another way:

a[C,:][:,C]

## array([[15, 12],
##        [ 3,  0]])

Yet, if I have a 2x2 array, call it b, setting

a[C,:][:,C]=b

doesn't work but

a[[[C[0],C[0]],[C[1],C[1]]],[[C[0],C[1]],[C[0],C[1]]]]=b

does.

Why is this? And is this second way the most efficient possible? Thanks!

If C was a numpy array then a[C[:,None], C] would work as well. — Ashwini Chaudhary
– Ashwini Chaudhary, Commented Feb 18, 2015 at 0:55
Ashwini: when I try that I get 'TypeError: list indices must be integers, not tuple'? — mww
– mww, Commented Feb 18, 2015 at 1:11
@user65257 I did say: "If C was a numpy array", so it won't work with a list. — Ashwini Chaudhary
– Ashwini Chaudhary, Commented Feb 18, 2015 at 1:14

hpaulj · Accepted Answer · 2015-02-18 03:35:42Z

The relevant section from the numpy docs is http://docs.scipy.org/doc/numpy/reference/arrays.indexing.html#purely-integer-array-indexing Advanced array indexing.

Adapting that example to your case:

In [213]: rows=np.array([[C[0],C[0]],[C[1],C[1]]])
In [214]: cols=np.array([[C[0],C[1]],[C[0],C[1]]])

In [215]: rows
array([[3, 3],
       [0, 0]])

In [216]: cols
array([[3, 0],
       [3, 0]])

In [217]: a[rows,cols]
array([[15, 12],
       [ 3,  0]])

due to broadcasting, you don't need to repeat duplicate indices, thus:

a[[[3],[0]],[3,0]]

does just fine. np.ix_ is a convenience function to produce just such a pair:

np.ix_(C,C) 
(array([[3],
        [0]]), 
 array([[3, 0]]))

thus a short answer is:

a[np.ix_(C,C)]

A related function is meshgrid, which constructs full indexing arrays:

a[np.meshgrid(C,C,indexing='ij')]

np.meshgrid(C,C,indexing='ij') is the same as your [rows, cols]. See the functions doc for the significance of the 'ij' parameter.

np.meshgrid(C,C,indexing='ij',sparse=True) produces the same pair of arrays as np.ix_.

I don't think there's a serious difference in computational speed. Obviously some require less typing on your part.

a[:,C][C,:] works for viewing values, but not for modifying them. The details have to do with which actions make views and which make copies. The simple answer is, use only one layer of indexing if you want to modify values.

The indexing documentation:

Thus, x[ind1,...,ind2,:] acts like x[ind1][...,ind2,:] under basic slicing.

Thus a[1][3] += 7 works. But the doc also warns

Warning The above is not true for advanced indexing.

Collectives™ on Stack Overflow

Modifying block matrices in Python

1 Answer 1

Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related