I'm trying to find an efficient way of applying list convolve to a 3d matrix of 0 and 1, where about 5% of the matrix are 1, but the 1s may either be randomly placed, or grouped in one large contiguous "island" (note that this relates to a separate question I asked here about more efficient ways to generate this matrix).

Essentially I want to convolve a set of small (3x3x3) kernels with this large matrix of 1s and 0s. I have an implementation (see below), but I am wondering if there is a smarter way to take advantage of the sparse nature of my matrix - the fact that the 1s are randomly placed makes this difficult though.

Ultimately my goal in doing this is to count the number of "three-neighbor interactions" (i.e. the number of 1-1-1, 1-1-0, 1-0-0, and 0-0-0 interactions) in the matrix, where these interactions could have some arbitrary geometry specified by 3x3x3 kernel. I have tried thinking if there would be a more efficient way to implement this via encoding the matrix as a graph of some sort, but am not familiar enough with that infrastructure in Mathematica to make much headway. Any help or advice to optimizing this calculation would be greatly appreciated.

The working code: given a large but sparse 3 dimensional matrix, M (which here I create via GrowRandArray), and a set of small 3x3x3 kernels (which here I just randomly pick from all possible 3x3x3 permutations of binary arrays with 3 1s), is there a more efficient way to implement "GetCounts" below:

GrowRandArray[percent_, size_] := Module[{locs, phaseb, phase},
  locs = Flatten[
    Table[{i, j, k}, {i, 1, size}, {j, 1, size}, {k, 1, size}], 2];
  phaseb = RandomSample[locs, Ceiling[percent*size^3]];
  phase = 
    phaseb -> ConstantArray[1, Length[phaseb]], {size, size, size}]
  ] (*this is just to generate the large matrix for this example *)

GetCounts[Mat_, kernels_, num_] := 
  Count[Flatten[ListConvolve[#, Mat, 1] & /@ kernels], #] & /@ 
   Range[0, num];

M = GrowRandArray[0.1, 30];
kernelsFull = 
  ArrayReshape[#, {3, 3, 3}] & /@ 
   Permutations[{1, 1, 1}~Join~ConstantArray[0, 3^3 - 3]];
kernels = RandomSample[kernelsFull, 20];
GetCounts[M, kernels, 3]

Making ListConvolve output a SparseArray and using Counts or Tally on the "NonzeroValues" of sparse array:

ClearAll[nzvcounts, nzvtally]
nzvcounts[m_, k_] := Module[{
    counts = Counts[SparseArray[ListConvolve[#, m, 1] & /@ k]["NonzeroValues"]]}, 
  Lookup[counts, {0, 1, 2, 3}, -Total@counts + Length[k] Times @@ Dimensions[m]]]

nzvtally[m_, k_] := Module[{t = Length[k] Times @@ Dimensions[m],
    nzv = SparseArray[ListConvolve[#, m, 1] & /@ k ]["NonzeroValues"]},
  Flatten[{t - Total[#], #}] & @ SortBy[Tally[nzv], First][[All, 2]]]

Both are slightly faster than Carl's gc for M and kernels:

nzvcounts[M, kernels] // RepeatedTiming

{0.145, {393820, 130892, 14756, 532}}

nzvtally[M, kernels] // RepeatedTiming

{0.143, {393820, 130892, 14756, 532}}

Lookup[gc[M, kernels], {0, 1, 2, 3}] // RepeatedTiming

{0.18, {393820, 130892, 14756, 532}}

GetCounts[M, kernels, 3] // RepeatedTiming

{0.68, {393820, 130892, 14756, 532}}


At a minimum, it would help if you only ListConvolve a kernel once instead of once for each number between 0 and num:

gc[m_, k_] := Counts@Flatten[ListConvolve[#,m,1]&/@k]


r1 = GetCounts[M, kernels, 3]; //RepeatedTiming
r2 = gc[M, kernels]; //RepeatedTiming

r1 == Lookup[r2, {0, 1, 2, 3}]

{0.557, Null}

{0.16, Null}


  • $\begingroup$ Ah yes... overlooked the obvious on that one I suppose $\endgroup$ – KHAAAAAAAAN Jun 15 '19 at 7:35

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.