Tjyylli

Why did the Germans forbid the possession of pet pigeons in Rostov-on-Don in 1941?

Is ipsum/ipsa/ipse a third person pronoun, or can it serve other functions?

How is it possible for user's password to be changed after storage was encrypted? (on OS X, Android)

Where else does the Shulchan Aruch quote an authority by name?

What to wear for invited talk in Canada

Can the Produce Flame cantrip be used to grapple, or as an unarmed strike, in the right circumstances?

Extreme, but not acceptable situation and I can't start the work tomorrow morning

COUNT(*) or MAX(id) - which is faster?

LWC and complex parameters

Was there ever an axiom rendered a theorem?

Are objects structures and/or vice versa?

Copycat chess is back

Why was the "bread communication" in the arena of Catching Fire left out in the movie?

Is it wise to focus on putting odd beats on left when playing double bass drums?

Unbreakable Formation vs. Cry of the Carnarium

What does it exactly mean if a random variable follows a distribution

Is Social Media Science Fiction?

Does the average primeness of natural numbers tend to zero?

Need help identifying/translating a plaque in Tangier, Morocco

Domain expired, GoDaddy holds it and is asking more money

Calculate Levenshtein distance between two strings in Python

Why doesn't a const reference extend the life of a temporary object passed via a function?

"My colleague's body is amazing"

Why is the design of haulage companies so “special”?

How to efficiently unroll a matrix by value with numpy?

How to create mask with one at specified index?How to merge two dictionaries in a single expression?How do I check if a list is empty?How do I check whether a file exists without exceptions?How do I check if an array includes an object in JavaScript?How to append something to an array?How do I sort a dictionary by value?How do I determine whether an array contains a particular value in Java?How do I list all files of a directory?How do I remove a particular element from an array in JavaScript?How to use foreach with array in JavaScript?

.everyoneloves__top-leaderboard:empty,.everyoneloves__mid-leaderboard:empty,.everyoneloves__bot-mid-leaderboard:empty{ height:90px;width:728px;box-sizing:border-box;
}

7

I have a matrix M with values 0 through N within it. I'd like to unroll this matrix to create a new matrix A where each submatrix A[i, :, :] represents whether or not M == i.

The solution below uses a loop.

# Example Setup

import numpy as np



np.random.seed(0)

N = 5

M = np.random.randint(0, N, size=(5,5))



# Solution with Loop

A = np.zeros((N, M.shape[0], M.shape[1]))

for i in range(N):

    A[i, :, :] = M == i

This yields:

M

array([[4, 0, 3, 3, 3],

       [1, 3, 2, 4, 0],

       [0, 4, 2, 1, 0],

       [1, 1, 0, 1, 4],

       [3, 0, 3, 0, 2]])



M.shape

# (5, 5)





A 

array([[[0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [1, 0, 0, 0, 1],

        [0, 0, 1, 0, 0],

        [0, 1, 0, 1, 0]],

       ...

       [[1, 0, 0, 0, 0],

        [0, 0, 0, 1, 0],

        [0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [0, 0, 0, 0, 0]]])



A.shape

# (5, 5, 5)

Is there a faster way, or a way to do it in a single numpy operation?

python arrays numpy

share|improve this question

edited Apr 5 at 22:24

coldspeed

141k25156241

asked Apr 5 at 22:06

seveibarseveibar

1,29911225

• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  It would be better if you explain it in detail.
  
  – Marios Nikolaou
  Apr 5 at 22:16
  

  
  
  
  


• 
  
  
  

  
  2
  

  
  
  
  
  
  

  
  
  @MariosNikolaou just copy/paste his code and print(M);print(A)...I edited it for you though
  
  – Reedinationer
  Apr 5 at 22:18
  
  
  

  
  
  
  


• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  @Reedinationer i did it.
  
  – Marios Nikolaou
  Apr 5 at 22:21
  

  
  
  
  


• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  I would not recommend pasting output for this code as the input is randomised without a seed.
  
  – coldspeed
  Apr 5 at 22:21
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  i == M compare int with array 5x5 ? and then save it in A?
  
  – Marios Nikolaou
  Apr 5 at 22:26
  

  
  
  
  

 | 
show 2 more comments

7

I have a matrix M with values 0 through N within it. I'd like to unroll this matrix to create a new matrix A where each submatrix A[i, :, :] represents whether or not M == i.

The solution below uses a loop.

# Example Setup

import numpy as np



np.random.seed(0)

N = 5

M = np.random.randint(0, N, size=(5,5))



# Solution with Loop

A = np.zeros((N, M.shape[0], M.shape[1]))

for i in range(N):

    A[i, :, :] = M == i

This yields:

M

array([[4, 0, 3, 3, 3],

       [1, 3, 2, 4, 0],

       [0, 4, 2, 1, 0],

       [1, 1, 0, 1, 4],

       [3, 0, 3, 0, 2]])



M.shape

# (5, 5)





A 

array([[[0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [1, 0, 0, 0, 1],

        [0, 0, 1, 0, 0],

        [0, 1, 0, 1, 0]],

       ...

       [[1, 0, 0, 0, 0],

        [0, 0, 0, 1, 0],

        [0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [0, 0, 0, 0, 0]]])



A.shape

# (5, 5, 5)

Is there a faster way, or a way to do it in a single numpy operation?

python arrays numpy

share|improve this question

edited Apr 5 at 22:24

coldspeed

141k25156241

asked Apr 5 at 22:06

seveibarseveibar

1,29911225

• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  It would be better if you explain it in detail.
  
  – Marios Nikolaou
  Apr 5 at 22:16
  

  
  
  
  


• 
  
  
  

  
  2
  

  
  
  
  
  
  

  
  
  @MariosNikolaou just copy/paste his code and print(M);print(A)...I edited it for you though
  
  – Reedinationer
  Apr 5 at 22:18
  
  
  

  
  
  
  


• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  @Reedinationer i did it.
  
  – Marios Nikolaou
  Apr 5 at 22:21
  

  
  
  
  


• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  I would not recommend pasting output for this code as the input is randomised without a seed.
  
  – coldspeed
  Apr 5 at 22:21
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  i == M compare int with array 5x5 ? and then save it in A?
  
  – Marios Nikolaou
  Apr 5 at 22:26
  

  
  
  
  

 | 
show 2 more comments

I have a matrix M with values 0 through N within it. I'd like to unroll this matrix to create a new matrix A where each submatrix A[i, :, :] represents whether or not M == i.

The solution below uses a loop.

# Example Setup

import numpy as np



np.random.seed(0)

N = 5

M = np.random.randint(0, N, size=(5,5))



# Solution with Loop

A = np.zeros((N, M.shape[0], M.shape[1]))

for i in range(N):

    A[i, :, :] = M == i

This yields:

M

array([[4, 0, 3, 3, 3],

       [1, 3, 2, 4, 0],

       [0, 4, 2, 1, 0],

       [1, 1, 0, 1, 4],

       [3, 0, 3, 0, 2]])



M.shape

# (5, 5)





A 

array([[[0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [1, 0, 0, 0, 1],

        [0, 0, 1, 0, 0],

        [0, 1, 0, 1, 0]],

       ...

       [[1, 0, 0, 0, 0],

        [0, 0, 0, 1, 0],

        [0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [0, 0, 0, 0, 0]]])



A.shape

# (5, 5, 5)

Is there a faster way, or a way to do it in a single numpy operation?

python arrays numpy

share|improve this question

edited Apr 5 at 22:24

coldspeed

141k25156241

asked Apr 5 at 22:06

seveibarseveibar

1,29911225

# Example Setup

import numpy as np



np.random.seed(0)

N = 5

M = np.random.randint(0, N, size=(5,5))



# Solution with Loop

A = np.zeros((N, M.shape[0], M.shape[1]))

for i in range(N):

    A[i, :, :] = M == i

M

array([[4, 0, 3, 3, 3],

       [1, 3, 2, 4, 0],

       [0, 4, 2, 1, 0],

       [1, 1, 0, 1, 4],

       [3, 0, 3, 0, 2]])



M.shape

# (5, 5)





A 

array([[[0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [1, 0, 0, 0, 1],

        [0, 0, 1, 0, 0],

        [0, 1, 0, 1, 0]],

       ...

       [[1, 0, 0, 0, 0],

        [0, 0, 0, 1, 0],

        [0, 1, 0, 0, 0],

        [0, 0, 0, 0, 1],

        [0, 0, 0, 0, 0]]])



A.shape

# (5, 5, 5)

share|improve this question

edited Apr 5 at 22:24

coldspeed

141k25156241

asked Apr 5 at 22:06

seveibarseveibar

1,29911225

share|improve this question

edited Apr 5 at 22:24

coldspeed

141k25156241

asked Apr 5 at 22:06

seveibarseveibar

1,29911225

edited Apr 5 at 22:24

coldspeed

141k25156241

edited Apr 5 at 22:24

coldspeed

141k25156241

asked Apr 5 at 22:06

seveibarseveibar

1,29911225

asked Apr 5 at 22:06

seveibarseveibar

1,29911225

3 Answers
3

active

oldest

votes

6

Broadcasted comparison is your friend:

B = (M[None, :] == np.arange(N)[:, None, None]).view(np.int8)



 np.array_equal(A, B)

# True

The idea is to expand the dimensions in such a way that the comparison can be broadcasted in the manner desired.

---

As pointed out by @Alex Riley in the comments, you can use np.equal.outer to avoid having to do the indexing stuff yourself,

B = np.equal.outer(np.arange(N), M).view(np.int8)



np.array_equal(A, B)

# True

share|improve this answer

edited Apr 5 at 22:45

answered Apr 5 at 22:19

coldspeedcoldspeed

141k25156241

• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  Good answer - just to point out there's a superfluous newaxis in your indexing for M (resulting in a 4D array). You could use M[None, :] instead to get the 3D array. An alternative to avoid fiddly indexing is to use np.equal.outer(np.arange(N), M).view(np.int8).
  
  – Alex Riley
  Apr 5 at 22:34
  
  
  

  
  
  
  


• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  @AlexRiley Thanks for that! And the outer solution is quite neat.
  
  – coldspeed
  Apr 5 at 22:45
  

  
  
  
  

add a comment | 

6

You can make use of some broadcasting here:

P = np.arange(N)

Y = np.broadcast_to(P[:, None], M.shape)

T = np.equal(M, Y[:, None]).astype(int)

---

Alternative using indices:

X, Y = np.indices(M.shape)

Z = np.equal(M, X[:, None]).astype(int)

share|improve this answer

edited Apr 5 at 22:39

answered Apr 5 at 22:19

user3483203user3483203

31.9k82857

• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  This answer was really helpful towards my understanding of broadcasting, thank you!
  
  – seveibar
  Apr 6 at 0:58
  

  
  
  
  

add a comment | 

3

You can index into the identity matrix like so

 A = np.identity(N, int)[:, M]

or so

 A = np.identity(N, int)[M.T].T

Or use the new (v1.15.0) put_along_axis

A = np.zeros((N,5,5), int)

np.put_along_axis(A, M[None], 1, 0)

Note if N is much larger than 5 then creating an NxN identity matrix may be considered wasteful. We can mitigate this using stride tricks:

def read_only_identity(N, dtype=float):

    z = np.zeros(2*N-1, dtype)

    s, = z.strides

    z[N-1] = 1

    return np.lib.stride_tricks.as_strided(z[N-1:], (N, N), (-s, s))

share|improve this answer

edited Apr 6 at 1:08

answered Apr 5 at 23:45

Paul PanzerPaul Panzer

31.6k21845

• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  This is great, really interesting answer.
  
  – user3483203
  Apr 5 at 23:47
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  Hi Paul, this answer is really elegant, but the identity answers seem specific to the case where the N=M.shape[0]=M.shape[1]. Is the solution similarly elegant for N=/=M.shape[0]=/=M.shape[1]? Thanks for the answer, learning a lot!
  
  – seveibar
  Apr 6 at 0:55
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  @seveibar I've updated the answer. It is really just a matter of replacing the correct 5s with Ns.
  
  – Paul Panzer
  Apr 6 at 1:13
  

  
  
  
  

add a comment | 

Your Answer

StackExchange.ifUsing("editor", function () {
StackExchange.using("externalEditor", function () {
StackExchange.using("snippets", function () {
StackExchange.snippets.init();
});
});
}, "code-snippets");

StackExchange.ready(function() {
var channelOptions = {
tags: "".split(" "),
id: "1"
};
initTagRenderer("".split(" "), "".split(" "), channelOptions);

StackExchange.using("externalEditor", function() {
// Have to fire editor after snippets, if snippets enabled
if (StackExchange.settings.snippets.snippetsEnabled) {
StackExchange.using("snippets", function() {
createEditor();
});
}
else {
createEditor();
}
});

function createEditor() {
StackExchange.prepareEditor({
heartbeatType: 'answer',
autoActivateHeartbeat: false,
convertImagesToLinks: true,
noModals: true,
showLowRepImageUploadWarning: true,
reputationToPostImages: 10,
bindNavPrevention: true,
postfix: "",
imageUploader: {
brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
allowUrls: true
},
onDemand: true,
discardSelector: ".discard-answer"
,immediatelyShowMarkdownHelp:true
});


}
});

draft saved

draft discarded

Sign up or log in

StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});

Sign up using Google

Sign up using Facebook

Sign up using Email and Password

Post as a guest

Name

Email

Required, but never shown

StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fstackoverflow.com%2fquestions%2f55543949%2fhow-to-efficiently-unroll-a-matrix-by-value-with-numpy%23new-answer', 'question_page');
}
);

Post as a guest

Name

Email

Required, but never shown

6

Broadcasted comparison is your friend:

B = (M[None, :] == np.arange(N)[:, None, None]).view(np.int8)



 np.array_equal(A, B)

# True

The idea is to expand the dimensions in such a way that the comparison can be broadcasted in the manner desired.

---

As pointed out by @Alex Riley in the comments, you can use np.equal.outer to avoid having to do the indexing stuff yourself,

B = np.equal.outer(np.arange(N), M).view(np.int8)



np.array_equal(A, B)

# True

share|improve this answer

edited Apr 5 at 22:45

answered Apr 5 at 22:19

coldspeedcoldspeed

141k25156241

• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  Good answer - just to point out there's a superfluous newaxis in your indexing for M (resulting in a 4D array). You could use M[None, :] instead to get the 3D array. An alternative to avoid fiddly indexing is to use np.equal.outer(np.arange(N), M).view(np.int8).
  
  – Alex Riley
  Apr 5 at 22:34
  
  
  

  
  
  
  


• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  @AlexRiley Thanks for that! And the outer solution is quite neat.
  
  – coldspeed
  Apr 5 at 22:45
  

  
  
  
  

add a comment | 

6

Broadcasted comparison is your friend:

B = (M[None, :] == np.arange(N)[:, None, None]).view(np.int8)



 np.array_equal(A, B)

# True

The idea is to expand the dimensions in such a way that the comparison can be broadcasted in the manner desired.

---

As pointed out by @Alex Riley in the comments, you can use np.equal.outer to avoid having to do the indexing stuff yourself,

B = np.equal.outer(np.arange(N), M).view(np.int8)



np.array_equal(A, B)

# True

share|improve this answer

edited Apr 5 at 22:45

answered Apr 5 at 22:19

coldspeedcoldspeed

141k25156241

• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  Good answer - just to point out there's a superfluous newaxis in your indexing for M (resulting in a 4D array). You could use M[None, :] instead to get the 3D array. An alternative to avoid fiddly indexing is to use np.equal.outer(np.arange(N), M).view(np.int8).
  
  – Alex Riley
  Apr 5 at 22:34
  
  
  

  
  
  
  


• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  @AlexRiley Thanks for that! And the outer solution is quite neat.
  
  – coldspeed
  Apr 5 at 22:45
  

  
  
  
  

add a comment | 

Broadcasted comparison is your friend:

B = (M[None, :] == np.arange(N)[:, None, None]).view(np.int8)



 np.array_equal(A, B)

# True

The idea is to expand the dimensions in such a way that the comparison can be broadcasted in the manner desired.

---

As pointed out by @Alex Riley in the comments, you can use np.equal.outer to avoid having to do the indexing stuff yourself,

B = np.equal.outer(np.arange(N), M).view(np.int8)



np.array_equal(A, B)

# True

share|improve this answer

edited Apr 5 at 22:45

answered Apr 5 at 22:19

coldspeedcoldspeed

141k25156241

B = (M[None, :] == np.arange(N)[:, None, None]).view(np.int8)



 np.array_equal(A, B)

# True

B = np.equal.outer(np.arange(N), M).view(np.int8)



np.array_equal(A, B)

# True

share|improve this answer

edited Apr 5 at 22:45

answered Apr 5 at 22:19

coldspeedcoldspeed

141k25156241

answered Apr 5 at 22:19

coldspeedcoldspeed

141k25156241

answered Apr 5 at 22:19

coldspeedcoldspeed

141k25156241

6

You can make use of some broadcasting here:

P = np.arange(N)

Y = np.broadcast_to(P[:, None], M.shape)

T = np.equal(M, Y[:, None]).astype(int)

---

Alternative using indices:

X, Y = np.indices(M.shape)

Z = np.equal(M, X[:, None]).astype(int)

share|improve this answer

edited Apr 5 at 22:39

answered Apr 5 at 22:19

user3483203user3483203

31.9k82857

• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  This answer was really helpful towards my understanding of broadcasting, thank you!
  
  – seveibar
  Apr 6 at 0:58
  

  
  
  
  

add a comment | 

6

You can make use of some broadcasting here:

P = np.arange(N)

Y = np.broadcast_to(P[:, None], M.shape)

T = np.equal(M, Y[:, None]).astype(int)

---

Alternative using indices:

X, Y = np.indices(M.shape)

Z = np.equal(M, X[:, None]).astype(int)

share|improve this answer

edited Apr 5 at 22:39

answered Apr 5 at 22:19

user3483203user3483203

31.9k82857

• 
  
  
  

  
  

  
  
  
  
  
  

  
  
  This answer was really helpful towards my understanding of broadcasting, thank you!
  
  – seveibar
  Apr 6 at 0:58
  

  
  
  
  

add a comment | 

You can make use of some broadcasting here:

P = np.arange(N)

Y = np.broadcast_to(P[:, None], M.shape)

T = np.equal(M, Y[:, None]).astype(int)

---

Alternative using indices:

X, Y = np.indices(M.shape)

Z = np.equal(M, X[:, None]).astype(int)

share|improve this answer

edited Apr 5 at 22:39

answered Apr 5 at 22:19

user3483203user3483203

31.9k82857

P = np.arange(N)

Y = np.broadcast_to(P[:, None], M.shape)

T = np.equal(M, Y[:, None]).astype(int)

X, Y = np.indices(M.shape)

Z = np.equal(M, X[:, None]).astype(int)

share|improve this answer

edited Apr 5 at 22:39

answered Apr 5 at 22:19

user3483203user3483203

31.9k82857

answered Apr 5 at 22:19

user3483203user3483203

31.9k82857

answered Apr 5 at 22:19

user3483203user3483203

31.9k82857

3

You can index into the identity matrix like so

 A = np.identity(N, int)[:, M]

or so

 A = np.identity(N, int)[M.T].T

Or use the new (v1.15.0) put_along_axis

A = np.zeros((N,5,5), int)

np.put_along_axis(A, M[None], 1, 0)

Note if N is much larger than 5 then creating an NxN identity matrix may be considered wasteful. We can mitigate this using stride tricks:

def read_only_identity(N, dtype=float):

    z = np.zeros(2*N-1, dtype)

    s, = z.strides

    z[N-1] = 1

    return np.lib.stride_tricks.as_strided(z[N-1:], (N, N), (-s, s))

share|improve this answer

edited Apr 6 at 1:08

answered Apr 5 at 23:45

Paul PanzerPaul Panzer

31.6k21845

• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  This is great, really interesting answer.
  
  – user3483203
  Apr 5 at 23:47
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  Hi Paul, this answer is really elegant, but the identity answers seem specific to the case where the N=M.shape[0]=M.shape[1]. Is the solution similarly elegant for N=/=M.shape[0]=/=M.shape[1]? Thanks for the answer, learning a lot!
  
  – seveibar
  Apr 6 at 0:55
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  @seveibar I've updated the answer. It is really just a matter of replacing the correct 5s with Ns.
  
  – Paul Panzer
  Apr 6 at 1:13
  

  
  
  
  

add a comment | 

3

You can index into the identity matrix like so

 A = np.identity(N, int)[:, M]

or so

 A = np.identity(N, int)[M.T].T

Or use the new (v1.15.0) put_along_axis

A = np.zeros((N,5,5), int)

np.put_along_axis(A, M[None], 1, 0)

Note if N is much larger than 5 then creating an NxN identity matrix may be considered wasteful. We can mitigate this using stride tricks:

def read_only_identity(N, dtype=float):

    z = np.zeros(2*N-1, dtype)

    s, = z.strides

    z[N-1] = 1

    return np.lib.stride_tricks.as_strided(z[N-1:], (N, N), (-s, s))

share|improve this answer

edited Apr 6 at 1:08

answered Apr 5 at 23:45

Paul PanzerPaul Panzer

31.6k21845

• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  This is great, really interesting answer.
  
  – user3483203
  Apr 5 at 23:47
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  Hi Paul, this answer is really elegant, but the identity answers seem specific to the case where the N=M.shape[0]=M.shape[1]. Is the solution similarly elegant for N=/=M.shape[0]=/=M.shape[1]? Thanks for the answer, learning a lot!
  
  – seveibar
  Apr 6 at 0:55
  

  
  
  
  


• 
  
  
  

  
  1
  

  
  
  
  
  
  

  
  
  @seveibar I've updated the answer. It is really just a matter of replacing the correct 5s with Ns.
  
  – Paul Panzer
  Apr 6 at 1:13
  

  
  
  
  

add a comment | 

You can index into the identity matrix like so

 A = np.identity(N, int)[:, M]

or so

 A = np.identity(N, int)[M.T].T

Or use the new (v1.15.0) put_along_axis

A = np.zeros((N,5,5), int)

np.put_along_axis(A, M[None], 1, 0)

Note if N is much larger than 5 then creating an NxN identity matrix may be considered wasteful. We can mitigate this using stride tricks:

def read_only_identity(N, dtype=float):

    z = np.zeros(2*N-1, dtype)

    s, = z.strides

    z[N-1] = 1

    return np.lib.stride_tricks.as_strided(z[N-1:], (N, N), (-s, s))

share|improve this answer

edited Apr 6 at 1:08

answered Apr 5 at 23:45

Paul PanzerPaul Panzer

31.6k21845

 A = np.identity(N, int)[:, M]

 A = np.identity(N, int)[M.T].T

A = np.zeros((N,5,5), int)

np.put_along_axis(A, M[None], 1, 0)

def read_only_identity(N, dtype=float):

    z = np.zeros(2*N-1, dtype)

    s, = z.strides

    z[N-1] = 1

    return np.lib.stride_tricks.as_strided(z[N-1:], (N, N), (-s, s))

share|improve this answer

edited Apr 6 at 1:08

answered Apr 5 at 23:45

Paul PanzerPaul Panzer

31.6k21845

answered Apr 5 at 23:45

Paul PanzerPaul Panzer

31.6k21845

answered Apr 5 at 23:45

Paul PanzerPaul Panzer

31.6k21845