Function Repository Resource:

SubmatrixReplace

Source Notebook

Efficiently replace rectangular submatrices of an input data matrix

Contributed by: Bradley Klee

ResourceFunction["SubmatrixReplace"][data,rules]

applies rules to all relevant rectangular submatrices of the data matrix, by default taking the first sorted value when replacements collide on an element.

Details

Each rule in rules must map between submatrices of fixed rectangular dimensions m×n, but m and n may vary across rules.
The replacement procedure groups rules by submatrix dimension, and applies them to individual submatrices using ReplaceAll.
When application of rules over all submatrices leads to two or more plausible image values for one particular matrix element, a collision is said to occur.
Collisions occur in two different cases:
Applying only one replacement rule may cause a collision within a region where two submatrices overlap.
Applying many distinct replacement rules may cause a collision within a region where two submatrices overlap.
Collisions do not occur in only one submatrix, because the underlying mechanism involves ReplaceAll.
ResourceFunction["SubmatrixReplace"] takes the following options:
"Cyclic"{True,True}determines cyclic or terminal boundary conditions for row and columns
CollisionFunctionFirst[#2]&determines how to resolve collisions
Method"Structured"specifies how to handle submatrices
The CollisionFunction acts locally on elements by returning an output value from three sequential inputs:
#1a single initial value
#2a list containing the stack of candidate replacement values
#3a list containing row and column indices
Choosing Method "Structured" uses the resource function BlockSubmatrices to build an optimized data structure.
Choosing Method "Simple" uses raw output from Partition.

Examples

Basic Examples (3) 

Replace the central diagonal of an identity matrix:

In[1]:=
ResourceFunction["SubmatrixReplace"][IdentityMatrix[6],
  {IdentityMatrix[2] -> Reverse[IdentityMatrix[2]]}] // MatrixForm
Out[1]=

Use a functional notation to trim the output of DiamondMatrix:

In[2]:=
ArrayPlot[ResourceFunction["SubmatrixReplace"][
  DiamondMatrix[10], { {{x_, w_}, {y_, z_}
     } :> If[MatchQ[Count[{w, x, y, z}, 0], 2 | 3],
     {{0, 0}, {0, 0}} , {{x, w}, {y, z}}]}],
 ColorRules -> {0 -> White, 1 -> Black, 2 -> Red}]
Out[2]=

Notice collision stacking when highlighting pattern fragments within a HadamardMatrix output:

In[3]:=
ArrayPlot[
 ResourceFunction["SubmatrixReplace"][4 HadamardMatrix[16], {
   {{1, 1}, {1, 1}} -> {{2, 2}, {2, 2}},
   {{1, 1}, {1, -1}} -> {{3, 3}, {3, -1}},
   {{1}} -> 4}],
 ColorRules -> {1 -> Black, -1 -> White, 2 -> Red, 3 -> Blue, 4 -> Green}]
Out[3]=

Change resolution of collisions by setting the CollisionFunction option:

In[4]:=
ArrayPlot[
 ResourceFunction["SubmatrixReplace"][4 HadamardMatrix[16], {
   {{1, 1}, {1, 1}} -> {{2, 2}, {2, 2}},
   {{1, 1}, {1, -1}} -> {{3, 3}, {3, -1}}, {{1}} -> 0},
  "CollisionFunction" -> (Max[ #2] &) ], ColorRules -> {1 -> Black, -1 -> White, 2 -> Red, 3 -> Blue, 0 -> Green}]
Out[4]=

Scope (3) 

Test occurrence of certain blocks in a pseudorandom matrix pattern:

In[5]:=
With[{data = RandomInteger[{1, 4}, {101, 101}]},
 ArrayPlot[
  ResourceFunction["SubmatrixReplace"][data, Join[Map[# -> # + 4 &,
     Mod[{{2, 1}, {3, 4}} + #, 4, 1] & /@ Range[4]], {
     {{x_, x_, x_, x_}} :> ConstantArray[0, {1, 4}],
     {{x_}, {x_}, {x_}, {x_}} :> ConstantArray[0, {4, 1}]}],
   "Cyclic" -> {False, True},
   "CollisionFunction" -> (Max[ #2 ] &) ],
  ColorRules -> {0 -> Black, 1 -> LightMagenta,
    2 -> LightBlue, 3 -> LightGreen,
    4 -> LightOrange, 5 -> Magenta,
    6 -> Blue, 7 -> Green, 8 -> Orange
    }, PixelConstrained -> 4]]
Out[5]=

Restrict search range to a central patch:

In[6]:=
With[{data = RandomInteger[{1, 4}, {101, 101}]},
 ArrayPlot[ResourceFunction["SubmatrixReplace"][data,
   Join[Map[# -> # + 4 &,
     Mod[{{2, 1}, {3, 4}} + #, 4, 1] & /@ Range[4]], {
     {{x_, x_, x_, x_}} :> ConstantArray[0, {1, 4}],
     {{x_}, {x_}, {x_}, {x_}} :> ConstantArray[0, {4, 1}]}],
   "Cyclic" -> {False, True},
   "CollisionFunction" -> (Max[ If[
        (#3 - {50, 50}) . (#3 - {50, 50}) < 25^2,
        #2, {#1}] ] &) ],
  ColorRules -> {0 -> Black, 1 -> LightMagenta,
    2 -> LightBlue, 3 -> LightGreen,
    4 -> LightOrange, 5 -> Magenta,
    6 -> Blue, 7 -> Green, 8 -> Orange
    }, PixelConstrained -> 4]]
Out[6]=

Timing test Method Simple against the task of highlighting Rule 110's periodic background:

In[7]:=
With[{data = CellularAutomaton[110,
    SeedRandom["TimeTest110"]; RandomInteger[1, 101], 100]},
 Column@Reverse[RepeatedTiming[
    ArrayPlot[ResourceFunction["SubmatrixReplace"][data,
      Map[# -> (# /. {0 -> 2}) &,
       {{{1, 1, 1, 1}, {1, 0, 0, 0}, {1, 0, 0, 1}, {1, 0, 1, 1}},
        {{1, 1, 1}, {1, 0, 1}, {1, 1, 1}}}],
      "Method" -> "Simple",
      "Cyclic" -> {False, True},
      "CollisionFunction" -> (Min[#2] &) ],
     ColorRules -> {1 -> Black, 0 -> White, 2 -> Lighter[Blue, 0.6]},
     PixelConstrained -> 4]]]]
Out[7]=

Repeat timing test with default Method Structured:

In[8]:=
With[{data = CellularAutomaton[110,
    SeedRandom["TimeTest110"]; RandomInteger[1, 101], 100]},
 Column@Reverse[RepeatedTiming[
    ArrayPlot[ResourceFunction["SubmatrixReplace"][data,
      Map[# -> (# /. {0 -> 2}) &,
       {{{1, 1, 1, 1}, {1, 0, 0, 0}, {1, 0, 0, 1}, {1, 0, 1, 1}},
        {{1, 1, 1}, {1, 0, 1}, {1, 1, 1}}}],
      "Method" -> "Structured",
      "Cyclic" -> {False, True},
      "CollisionFunction" -> (Min[ #2] &) ],
     ColorRules -> {1 -> Black, 0 -> White, 2 -> Lighter[Blue, 0.6]},
     PixelConstrained -> 4]]]]
Out[8]=

Compare complex-valued collision alternatives by inspecting elements:

In[9]:=
With[{collideRules = MapThread[Rule[#1, #2*#1] &, {
     {2 FourierMatrix[4][[1 ;; 3, 2 ;; 4]],
      2 FourierMatrix[4][[2 ;; 4, 2 ;; 4]]},
     {I, -I}}]}, Grid[{
   {Text@"Input Data", Text@"Collision Alternatives", SpanFromLeft},
   MatrixForm /@ {2 FourierMatrix[4],
     ResourceFunction["SubmatrixReplace"][2 FourierMatrix[4], collideRules],
     ResourceFunction["SubmatrixReplace"][2 FourierMatrix[4], collideRules,
      "CollisionFunction" -> (Last[#2] &)]},
   Style[#, Gray] & /@ {"", First, Last},
   {Text@"Sorted Range:", Complement[{1, -1, I, -I}, {}],
    SpanFromLeft}}, Frame -> True, FrameStyle -> Lighter[Gray],
  Dividers -> {{}, {-2 -> Lighter[Gray]}},
  Spacings -> {2, 1}]]
Out[9]=

Even using only one rule, a collision can still occur:

In[10]:=
With[{collideRule = {IdentityMatrix[2] -> {{2, 0}, {0, 3}}}}, Grid[{
   {Text@"Input Data", Text@"Collision Alternatives", SpanFromLeft},
   MatrixForm /@ {IdentityMatrix[2],
     ResourceFunction["SubmatrixReplace"][IdentityMatrix[2], collideRule],
     ResourceFunction["SubmatrixReplace"][IdentityMatrix[2], collideRule,
      "CollisionFunction" -> (Last[#2] &)]},
   Style[#, Gray] & /@ {"", First, Last},
   {Text@"Sorted Range:", Complement[{2, 3}, {}],
    SpanFromLeft}}, Frame -> True, FrameStyle -> Lighter[Gray],
  Dividers -> {{}, {-2 -> Lighter[Gray]}},
  Spacings -> {2, 1}]]
Out[10]=

Neat Examples (1) 

Highlight persistent structures in a 2D Cellular Automaton output:

In[11]:=
With[{CADeloneData = MapThread[
    If[SameQ[##], #1, 1] &, CellularAutomaton[
      {129733, {2, {{2, 2, 2}, {2, 1, 2}, {2, 2, 2}}}, {1, 1}},
      {{{1, 1}, {1, 1}}, 0}, {{1, 64 + 32 + 16, 2}}][[-16 ;; -1]], 2]},
 ArrayPlot[ResourceFunction["SubmatrixReplace"][CADeloneData,
   MapIndexed[#1 -> (#1 /. {0 -> (#2[[1]] + 1)}) &,
    {{{1, 1, 1, 1}, {1, 0, 0, 1}, {1, 0, 0, 1}, {1, 1, 1, 1}},
     {{1, 0, 1}, {0, 0, 0}, {1, 0, 1}},
     {{1, 1, 1}, {1, 0, 1}, {1, 1, 1}}}],
   "CollisionFunction" -> (Min[#2] &) ,
   "Cyclic" -> {False, False}],
  ColorRules -> {0 -> Black, 1 -> Black,
    2 -> Lighter[Blue, .6],
    3 -> Lighter[Red, .6],
    4 -> Lighter[Green, .6]},
  PixelConstrained -> 2]]
Out[11]=

Publisher

Brad Klee

Version History

  • 1.0.0 – 18 March 2022

Related Resources

Author Notes

Unproven in theory whether method "Structured" should always outperform "Simple", but use case tests can even show factors of 2, so it looks like we're onto something here.

License Information