mariottidae
mariottidae

Reputation: 61

All combinations of a numpy 2d array filled with 0s and 1s

Given K, I need to have all the possibile combinations of K x 2 numpy matrices so that in each matrix there are all 0s except for two 1s in different rows and columns. Something like this for K = 5:

  1. [[1,0],[0,1],[0,0],[0,0][0,0]]
  2. [[1,0],[0,0],[0,1],[0,0][0,0]]
  3. [[1,0],[0,0],[0,0],[0,1][0,0]]
  4. [[1,0],[0,0],[0,0],[0,0][0,1]]
  5. [[0,0],[1,0],[0,1],[0,0][0,0]]
  6. [[0,0],[1,0],[0,0],[0,1][0,0]]
  7. ... and so on

So the resulting array should be a K x 2 x (K*(K-1)/2). I want to avoid loops since it's not an efficient way when K is big enough (in my specific case K = 300)

Upvotes: 3

Views: 1138

Answers (3)

bb1
bb1

Reputation: 7863

Another possibility (with rearranged axes for clearer output):

from itertools import combinations
import numpy as np

k = 4
x = list(combinations(range(k), 2))
out = np.zeros((n := len(x), k, 2), dtype=int)
out[np.c_[:n], x, [0, 1]] = 1
print(out)

It gives:

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

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

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

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

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

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

Upvotes: 1

CB Madsen
CB Madsen

Reputation: 632

This is an oddly specific question, but an interesting problem, I'd love to know what the context is?

You are looking for all permutations of a multiset enter image description here, python's itertools doesn't currently support this. So simplest solution is to use the multiset tools of the sympy library.

The following code took about ~2.5 minutes to run on my machine, so is fairly fast for a single thread. You're looking at 179700 unique permutations for K=300.

(I took inspiration from https://stackoverflow.com/a/40289807/10739860)

from collections import Counter
from math import factorial, prod

import numpy as np
from sympy.utilities.iterables import multiset_permutations
from tqdm import tqdm


def No_multiset_permutations(multiset: list) -> int:
    """Calculates the No. possible permutations given a multiset.
    See: https://en.wikipedia.org/wiki/Permutation#Permutations_of_multisets

    :param multiset: List representing a multiset.
    """
    value_counts = Counter(multiset).values()
    denominator = prod([factorial(val) for val in value_counts])
    return int(factorial(len(multiset)) / denominator)


def multiset_Kx2_permutations(K: int) -> np.ndarray:
    """This will generate all possible unique Kx2 permutations of an array
    withsize K where two values are 1 and the rest are 0.

    :param K: The size of the array.
    """
    # Construct number multiset, e.g. K=5 gives [1, 1, 0, 0, 0, 0, 0, 0, 0, 0]
    numbers = [1, 1] + [0] * (K - 1) * 2

    # Use sympy's multiset_permutations to get a multiset permutation generator
    generator = multiset_permutations(numbers)

    # Calculate the No. possible permutations
    number_of_perms = No_multiset_permutations(numbers)

    # Get all permutations, bonus progress bar is included :)
    unique_perms = [next(generator) for _ in tqdm(range(number_of_perms))]

    # Reshape each permutation to Kx2
    unique_perms = np.array(unique_perms, dtype=np.int8)
    return unique_perms.reshape(-1, K, 2)


if __name__ == "__main__":
    solution = multiset_Kx2_permutations(300)

Upvotes: 1

I can't think of an elegant solution but here's a not-so-elegant pure numpy one:

import numpy as np

def combination_matrices(K):
    # get combination indices
    i, j = np.indices((K, K))
    comb_indices = np.transpose((i < j).nonzero())  # (num_combs, 2) array where ones are
    num_combs = comb_indices.shape[0]  # K*(K-1)/2

    # create a matrix of the desired shape, first axis enumerates combinations
    matrices = np.zeros((num_combs, K, 2), dtype=int)
    # broadcasting assignment of ones
    comb_range, col_index = np.ogrid[:num_combs, :2]
    matrices[comb_range, comb_indices, col_index] = 1
    return matrices

This first uses the indices of a (K, K)-shaped array to find the index pairs for every combination (these are indices that encode the upper triangle of the array, excluding the diagonal). Then we use a bit tricky broadcasting assignment (heavy fancy indexing) to set each corresponding element of the pre-allocated output array to 1.

Note that I put the K*(K-1)/2-sized axis first, because this makes the most sense in numpy with C-contiguous memory layout. This way when you take the matrix for combination index 3, arr[3, ...] will be a contiguous chunk of memory of shape (K, 2) that's fast to work with in vectorised operations.

The output for K = 4:

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

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

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

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

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

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

Upvotes: 2

Related Questions