Composite Index Updates for Numpy Matrices

Question

I am trying to update a series of specific rows and columns of a numpy matrix. Here is an example:

import numpy as np A=np.zeros((8,8)) rows=[0, 1, 5] columns=[2, 3] #(What I am trying to achieve) The following does not update A A[rows][:,columns]+=1 #while this just does for i in rows: A[i][columns]+=1 

Expected Result:

 In [1]:print(A) Out[1]: array([[ 0., 0., 1., 1., 0., 0., 0., 0.], [ 0., 0., 1., 1., 0., 0., 0., 0.], [ 0., 0., 0., 0., 0., 0., 0., 0.], [ 0., 0., 0., 0., 0., 0., 0., 0.], [ 0., 0., 0., 0., 0., 0., 0., 0.], [ 0., 0., 1., 1., 0., 0., 0., 0.], [ 0., 0., 0., 0., 0., 0., 0., 0.], [ 0., 0., 0., 0., 0., 0., 0., 0.]]) 

Is there a way to perform multiple and quick updates at the same time without a loop?

Answer 1

Indexing by arrays or index lists belongs to the category of fancy indexing and always generates a copy of the corresponding part of the array, and how to give you an idea of ​​the source data.

Assigning or modifying a place as a result of one fancy indexing operation works fine (for example, A[rows] += 1 ). However, changing the result of chained indexing expressions, for example. A[rows][:, columns] will not work if the first expression uses fantastic indexing. The reason is because A[rows] creates a copy, so adding 1 to A[rows][:, columns] will not affect A

You can avoid this problem if you do all your indexing in one go. At some point, you probably tried something like this:

 In [38]: A[rows, columns] += 1 --------------------------------------------------------------------------- ValueError Traceback (most recent call last) <ipython-input-38-7d2300f59d51> in <module>() ----> 1 A[rows, columns] += 1 ValueError: shape mismatch: objects cannot be broadcast to a single shape 

Here numpy outputs your index lists to arrays, then tries to broadcast them into the same form, which fails because they are incompatible sizes. To do this successfully, you can make rows into an array (n, 1) and columns into an array (1, m) , where m and n are the number of separate row / column indices. Thus, they can be expanded in size 1 to size (n, m) and used for indexing in A :

 r = np.array(rows) c = np.array(columns) A[r[:, None], c[None, :]] += 1 print(A) # [[ 0. 0. 1. 1. 0. 0. 0. 0.] # [ 0. 0. 1. 1. 0. 0. 0. 0.] # [ 0. 0. 0. 0. 0. 0. 0. 0.] # [ 0. 0. 0. 0. 0. 0. 0. 0.] # [ 0. 0. 0. 0. 0. 0. 0. 0.] # [ 0. 0. 1. 1. 0. 0. 0. 0.] # [ 0. 0. 0. 0. 0. 0. 0. 0.] # [ 0. 0. 0. 0. 0. 0. 0. 0.]] 

Indexing with None is the same as using np.newaxis - it np.newaxis a new size of size 1.

Actually there is a convenient function np.ix_ for generating multidimensional indices from a sequence of 1D indices:

 A = np.zeros((8, 8)) A[np.ix_(rows, columns)] += 1 

Answer 2

rows should be a column vector e.g.

 rows=[[0],[1],[5]] cols=[2,3] A[rows,cols]+=1 

Two-stage indexing sometimes works, A[rows][:,cols] , but not always. In particular, this is not the case when rows not a slice. A[rows] now a copy, so changing it does not change A

There are various ways to index this block. plonser's use of the product works, although I rarely saw that it was used with numpy .

np.ix_ is a convenient tool for this:

 In [70]: np.ix_([0,1,5],[2,3]) Out[70]: (array([[0], [1], [5]]), array([[2, 3]])) 

np.newaxis also turns a row vector into a column:

 rows=np.array([0,1,5]) cols=np.array([2,3]) A[rows[:,None],cols] 

This pairing of column and row vectors works because numpy passes them to create an array of indices (3,2) .

itertools.product produces the same set of indices, but as a list (or generator) of tuples

 In [80]: list(itertools.product([0,1,5],[2,3])) Out[80]: [(0, 2), (0, 3), (1, 2), (1, 3), (5, 2), (5, 3)] In [84]: tuple(np.array(list(itertools.product([0,1,5],[2,3]))).T) Out[84]: (array([0, 0, 1, 1, 5, 5]), array([2, 3, 2, 3, 2, 3])) 

Answer 3

You can do this by first getting all the coordinates using itertool.product

 import itertools coord = np.array(list(itertools.product(rows,columns))) 

and then using advanced indexing

 A[tuple(coord.T)] += 1 

As mentioned in other answers, advanced indexing usually gives a copy. However, when handled correctly, it still works as an argument here .