# Longest path in 0/1 matrix

I have a matrix of 0's and 1's forming a number of disjoint paths:

I would like to find the lengths of the paths, and from that "spectrum," the longest length (in the above example: 27, starting at {1,14}). (The shortest length possible is 3, just from how I generate these paths. There are never trees or cycles—just paths.) I can do it by identifying start 1-cells as either on the boundary, or having three 0-neighbors and one 1-neighbor, and then tracing from the start cell along the path. This is quite clunky.

Can anyone see a slicker method? Efficiency is not my main concern at the moment. Thanks for your ideas!

Here's the matrix displayed above:

{{0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1,
0, 0, 0, 0, 0, 0, 0}, {1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0,
1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1}, {0, 1, 0, 0, 0, 0, 0,
1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0}, {0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0,
0, 1, 0, 0, 1, 0, 0, 1, 0}, {0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 0,
1, 1, 1, 1, 0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 1}, {0, 0, 0, 0, 1,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0,
0, 0}, {0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {1, 1, 0, 0, 1, 0, 1, 1, 0, 1, 1, 0,
0, 1, 0, 0, 1, 0, 1, 1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1}, {0, 1, 0, 0,
1, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1,
0, 0, 1, 0}, {0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0}, {1, 1, 1, 0, 1, 1, 0, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 0, 0, 1, 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, 1, 0, 0,
1, 0, 0, 0, 0}, {0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0,
0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0}, {1, 1, 0, 0, 1, 1, 0, 1, 0,
0, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 1}, {0,
0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1,
0, 0, 0, 0, 0, 1, 0}, {0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0,
0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0}, {0, 1, 1, 1, 1, 0, 0,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1,
0}, {0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
0, 1, 0, 0, 1, 0, 0, 0, 0}, {0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0,
1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0}, {0, 1, 1, 0, 1,
1, 1, 1, 0, 1, 1, 0, 1, 1, 0, 0, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1,
1, 0}, {0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0,
0, 0, 1, 0, 0, 1, 0, 0, 0, 0}, {0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0,
0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0}, {0, 1, 0, 1,
1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 0, 0, 1,
1, 1, 1, 1}, {0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0}, {0, 0, 0, 0, 1, 0, 0, 1, 0, 0,
1, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 1,
1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 0, 1, 1, 0,
1, 1, 1, 1, 1}, {0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1,
0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0}, {0, 1, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0}, {1,
1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1,
0, 1, 1, 0, 1, 1, 0}, {0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}}

• You could use the image processing functions to find the connected components and sort them by length. The pixel positions can easily be mapped back to the matrix index
– rm -rf
Jun 21, 2014 at 20:46
• Possibly food for thought: blog.wolfram.com/2010/12/21/… and linked posts. Jun 21, 2014 at 20:57
• How about providing the matrix?
– Öskå
Jun 21, 2014 at 21:00
• @Öskå: Thanks for the suggestion; included now. Jun 21, 2014 at 21:29
• The answers attest to an amazing inventiveness of this community! Jun 22, 2014 at 21:07

As I expressed in my comment above, it is possible (and easy) to use the image processing functions for this. Taking m to be the matrix above the following steps illustrate the idea:

img = Image@m;
ComponentMeasurements[img, "PerimeterCount"]
(* {1 -> 3, 2 -> 27, 3 -> 9, 4 -> 6, 5 -> 15, 6 -> 3, 7 -> 6, 8 -> 3, 9 -> 3, 10 -> 3,
11 -> 3, 12 -> 3, 13 -> 3, 14 -> 3, 15 -> 3, 16 -> 6, 17 -> 18, 18 -> 6, 19 -> 3,
20 -> 3, 21 -> 3, 22 -> 3, 23 -> 12, 24 -> 3, 25 -> 3, 26 -> 12, 27 -> 15, 28 -> 6,
29 -> 6, 30 -> 3, 31 -> 3, 32 -> 9, 33 -> 9, 34 -> 9, 35 -> 3, 36 -> 3, 37 -> 6,
38 -> 3, 39 -> 3, 40 -> 3, 41 -> 3, 42 -> 6, 43 -> 15, 44 -> 6, 45 -> 9, 46 -> 3,
47 -> 3, 48 -> 6, 49 -> 3} *)


From the above, the largest length is 27 (which you can also find programmatically) and the corresponding path is:

longestPath = SelectComponents[img, "PerimeterCount", # == 27 &]


To get the list of indices:

SparseArray[ImageData@longestPath]["NonzeroPositions"]
(* {{1, 14}, {2, 14}, {3, 14}, {4, 14}, {5, 14}, {5, 15}, {5, 16}, {5, 17}, {6, 17}, {7, 17},
{8, 7}, {8, 8}, {8, 17}, {9, 8}, {9, 17}, {10, 8}, {10, 17}, {11, 8}, {11, 9}, {11, 10},
{11, 11}, {11, 12}, {11, 13}, {11, 14}, {11, 15}, {11, 16}, {11, 17}} *)

• Wow! Beautiful! Jun 21, 2014 at 21:30
• To find the largest length you may use: Max[List @@@( ComponentMeasurements[img, "PerimeterCount"])[[;; , 2]]] Jun 22, 2014 at 8:08
• You can also use SelectComponents[img, "PerimeterCount", -1] to get the longest component, and SelectComponents[img, "PerimeterCount", -3] to get the longest 3.
– kglr
May 26, 2019 at 7:51

Here is a Graph-based solution inspired by this Q&A where mat is your given matrix.

binaryGraph[mat_] := Module[{pos, edge, dedge}, pos = Position[mat, 1];
edge = Select[Subsets[Range@Length@pos, {2}],
Last@# - First@# <= (Max@Dimensions@mat + 1) &];
dedge = DeleteDuplicates[UndirectedEdge @@@ (Extract[edge, #] & /@
With[{dist = N@(EuclideanDistance[pos[[#]], pos[[#2]]] & @@@ edge)},
Flatten[Position[dist, #] & /@ DeleteDuplicates@N@Select[dist, # == 1 &]]])];
Graph[dedge, VertexCoordinates -> Rule @@@ Thread[{Range@Length@pos, ({#2, -#1} & @@@ pos)}],
VertexLabels -> "Name", ImagePadding -> 20]]

g = binaryGraph@mat;
pos = Flatten@Position[#, Max@#] &@(VertexEccentricity[g, #] & /@ Range@Length@VertexList@g);
HighlightGraph[g, PathGraph[FindShortestPath[g, First@pos, Last@pos]],
ImageSize -> 350]
Length@FindShortestPath[g, First@pos, Last@pos]


If one wishes to go back on a matrix-like side, the exact path taken from the matrix mat can be given by posmax:

posmax = Position[mat, 1][[#]] & /@ FindShortestPath[g, First@pos, Last@pos];
newmat = mat;
(newmat[[Sequence @@ #]] = Red) & /@ posmax
ArrayPlot[newmat, Mesh -> All, MeshStyle -> Directive[{Thick, Gray}]]


• I love that display of the matrix & paths--Very clever! Jun 21, 2014 at 21:59
• @JosephO'Rourke Thank you :) I kinda like it too ;o)
– Öskå
Jun 21, 2014 at 22:00
– eldo
Jun 21, 2014 at 22:36
• @eldo thanks :) the elements were here.., I just needed to put them together :P
– Öskå
Jun 21, 2014 at 23:25

Assuming m is the test matrix; the following is a graph based solution:

nbours[x_, y_] := Module[{xr, yr, nts},
xr = Select[{x - 1, x, x + 1}, # >= 1 && # <= Length@m &];
yr = Select[{y - 1, y, y + 1}, # >= 1 && # <= Length@m &];
nts = Cases[
Table[{xrp, y}, {xrp, xr}]~Join~Table[{x, yrp}, {yrp, yr}],
Except[{x, y}]];
If[m[[y]][[x]] == 1, {x, y} \[DirectedEdge] # & /@
Select[nts, Part @@ {m}~Join~Reverse[#] == 1 &], {}]
];

map = Flatten[Table[nbours[x, y], {x, 1, Length@m}, {y, 1, Length@m}]];
g = Graph[map]; dist = GraphDistanceMatrix[g];
max = Max[Replace[dist, \[Infinity] -> -1, 2]];
VertexList[g][[#]] & /@ First@Position[dist, max]


Produces the pair of points that form the longest path of distance max.

• Nice use of various tools, e.g., GraphDistanceMatrix[]. Jun 22, 2014 at 1:46

Using a slightly modified version this answer in a related Q/A:

ClearAll[m, edges, v, vcs, g, cc, ap, ap2, pos];
m = Module[{i = 1}, mat /. 1 :> i++] (* where mat is the input matrix *);
v = ComponentMeasurements[m, "Label"][[All, 1]];
vcs = ComponentMeasurements[m, "Centroid"];
edges = UndirectedEdge @@@ DeleteDuplicates[
Sort /@ Flatten[
Thread /@ ComponentMeasurements[m, "Neighbors", CornerNeighbors -> False]]];

g = Graph[v, edges, VertexCoordinates -> vcs, ImageSize -> 400,
VertexSize -> .5, VertexStyle -> Blue, EdgeStyle -> Thickness[.01]];
cc = ConnectedComponents[g]


Note: ConnectedComponents are sorted by length in version 9.0.1.0. For earlier versions you need to use SortBy[cc, -Length[#]&] (thanks: Oska)

ap = ArrayPlot[m, Mesh -> All, ImageSize -> 400, MeshStyle -> Directive[{Gray, Thick}],
ColorRules -> Join[Thread[cc[[1]] -> Red], {0 -> White, _ -> Black}]];

ap2 = ArrayPlot[m, Mesh -> All, ImageSize-> 400, MeshStyle -> Directive[{Gray, Thick}],
ColorData[{"DeepSeaColors", "Reverse"}][Length[#]/Max[Length /@ cc]]] & /@ cc,
{0 -> White, _ -> Black}}, 2]] (* paths color-coded by length *);

Row[{HighlightGraph[g, cc[[1]], GraphHighlightStyle -> "Thick"], ap, ap2}]


To get the positions:

pos = First[Position[m, #]] & /@ # & /@ cc;
pos[[1]] (* positions of the longest-path elements *)
(* {{8, 7}, {8, 8}, {9, 8}, {10, 8}, {11, 8}, {11, 9}, {11, 10},
{11, 11}, {11, 12}, {11, 13}, {11, 14}, {11, 15}, {11, 16},
{11, 17}, {10, 17}, {9, 17}, {8, 17}, {7, 17}, {6, 17}, {5, 17},
{5, 16}, {5, 15}, {5, 14}, {4, 14}, {3, 14}, {2, 14}, {1, 14}} *)

• I was wondering how to tweak your code to make it work without the "diagonal neighbours", and it simply is CornerNeighbors -> False..! :D It's always obvious after seeing the answer :P
– Öskå
Jun 22, 2014 at 16:16
• @Öskå, it took me a while to remember CornerNeighbors:)
– kglr
Jun 22, 2014 at 16:21
• I love the 3rd display using DeepSeaColors! Jun 22, 2014 at 16:24
• @Öskå & Joseph, thank you both. Oska, I added a note on version differences.
– kglr
Jun 22, 2014 at 16:34

Not nearly as clean as @rm -rf's solution, but here is another image processing solution. Once again m is taken to be the matrix above.

mc=MorphologicalComponents@m;

Max@Tally[Join @@ mc][[2 ;;, 2]]
(*27*)


If one wishes to vizualize the exact path:

p1=ArrayPlot[mc /. Rule @@@ Rest[Tally[Join@@mc]], Mesh -> All];
(* paths color-coded by length *)

p2=ArrayPlot[mc /. {Last@Commonest[Join@@mc, 2] -> Red, a_ /; a != 0 -> 1}, Mesh -> All];
(* colors longest path red *)

GraphicsRow[{p1, p2}]


Here is a direct approach w/o the graphics/graph functions.

Algorithm:

• create an array with every '1' changed to a unique number (pad zeros around array to avoid dealing with edge issues )
• loop over nonzero entries. If any of the four neighbors (on further thought you may only need to look left & down) swap indices to match - this is a global replacement over the whole matrix.
• repeat the process until no unmatched neighbors are found (for this example it gets it on the first pass, but i think for some cases the iteration is needed )
• the last step reassigns the region numbers so they are sequential and sorted by region size (not really needed depending on what you want to do with the result )

>

  (m = ArrayPad[m0 Partition[ Range[ Times @@ Dimensions[m0]] ,
Last@Dimensions[m0]], 1];
nonzero = Position[m, x_ /; x > 0, {2}];
NestWhile[(Function[{p},
If[ Length@# > 0 ,
m = m /. m[[Sequence @@ p]] -> Min[#] ] &@
Select[Table[ m[[Sequence @@ (p + s) ]] ,
{s , {{0, 1}, {0, -1}, {-1, 0}, {1, 0}} }] ,
(# != 0 && # != m[[Sequence @@ p]]) &] ] /@ nonzero; # ) & ,
m , (#1 =!= #2 ) &, 2];
m = m[[2 ;; -2, 2 ;; -2]] /. MapIndexed[ #1[[1]] -> First@#2 &,
Reverse@
SortBy[Tally@Cases[ Flatten@m , Except[0]], #[[2]] & ] ]) //Timing // First


0.0468

(surprisingly faster than ComponentMeasurements )

The result is a matrix with identifier number for each group (sorted by group size)

 nonzero = # - {1, 1} & /@ nonzero;(*shift because we removed the padding*)
Graphics[{
ColorData[60, "ColorList"][[Mod[m[[Sequence @@ #]] - 1, 21] + 1]],
Disk[ {#[[2]], -#[[1]]}, .5]} & /@ nonzero]
Reverse@SortBy[Tally@Cases[ Flatten@m , x_ /; x > 0], #[[2]] & ]


 Reverse@SortBy[Tally@Cases[ Flatten@m , x_ /; x > 0], #[[2]] & ]


{{1, 27}, {2, 18}, {5, 15}, {4, 15}, {3, 15}, {7, 12}, ...

• I'm sure it has a name, but i've not formally studied the topic. At a glance It doesn't look exactly like the algorithms on that page, since I'm using mathematicas global pattern replace rather than iteratively flipping neighbors. Jun 23, 2014 at 20:03
• OK, good job. +1. Jun 23, 2014 at 20:08
• Impressive! I know this is supposed to be self-evident from the code, but I for one could use at least a high-level description of how this works... Jun 23, 2014 at 20:22
• @george2079 - Great piece of code, but I don't like the lollipop-colors (It's difficult for me to immediately see the longest path).
– eldo
Jun 23, 2014 at 22:24
• Try this Hue[2/3 (m[[Sequence @@ #]]-1)/(Last@Union@Flatten@m]-1) instead of ColorData[] Jun 23, 2014 at 22:40

There are already several really good answers. I simply want to point out that this can be cast as using Graph in a way that is explicitly quite efficient. The steps are as follows.

(1) Make positions corresponding to 1 into vertices. (2) Find immediate neighbors and make all pairs into edges. (3) Create a Graph object with the above data and VertexCoordinates corresponding to the vertex coordinates (this is not quite the tautology it sounds like). (4) Find the connected components.

If we start with the original binary matrix then code to do this is as below. The reverse-transpose in the first step is to get the correct orientation for coordinates to correspond to the original picture.

smat = SparseArray[Map[Reverse, Transpose[mat]]];
verts = Most[ArrayRules[smat]][[All, 1]];


To get edges efficiently one might use a Nearest formulation. Alternatively one can use a loop approach that is effectively linear time in the number of vertices. With the slightly less efficient Nearest approach we will afterward remove duplicates where first vertex is lexicographically greater than second (since we'll be retaining the oppositely-ordered equivalents).

nf = Nearest[verts -> Thread[v[verts]]];
nbors = Map[nf[#, {5, 1}] &, verts];
edges = Flatten[
v -> Identity];
edges2 = Select[
edges, #[[1, 1]] <= #[[2, 1]] && #[[1, 2]] <= #[[2, 2]] &];


Now we make this graph.

gr = Graph[edges];


Place the vertices where they belong. (Surely there must be a better way to do this in the graph creation step...)

PropertyValue[gr, VertexCoordinates] = VertexList[gr];


Extract the connected components.

paths = WeaklyConnectedComponents[gr];
pathlens = Map[Length, paths];
maxlen = Max[pathlens]

(* Out[279]= 27 *)


I will make no effort to compete with the far nicer pictures that have already been shown*.

HighlightGraph[gr,
Subgraph[gr, SelectFirst[paths, Length[#] == maxlen &]] ]


*Mostly because I'd lose. And badly.

• Thanks for the annotations of your code! Jun 24, 2014 at 23:26

On the Wikipedia page for connected-component labeling a couple of algorithms are explained that probably are similar to what the component image processing functions are doing. Although I'm not really satisfied with the code, I'll go ahead and post an implementation of the two-pass algorithm, since it seems to be pretty popular for this task.

listNeighbors[matrix_, pos_] := DeleteCases[
Extract[ArrayPad[matrix, 1], pos + {1, 1} + # & /@ {
{1, 1}, {1, -1}, {1, 0}, {-1, 1}, {-1, -1}, {-1, 0}, {0,
1}, {0, -1}
}], 0]

firstPass[matrix_] := Module[{linked, labels, nextLabel, neighbors, rows, cols},
labels = ConstantArray[0, {rows, cols} = Dimensions[matrix]];
nextLabel = 1;

Table[
If[
Extract[matrix, {i, j}] != 0,
neighbors = listNeighbors[labels, {i, j}];
If[
Length[neighbors] == 0,
labels[[i, j]] = nextLabel;
nextLabel++;,
labels[[i, j]] = Min[neighbors];
Do[
{label, neighbors}
]
]
],
{i, rows}, {j, cols}
];

secondPass[firstPass[matrix]] // Colorize

The longest path can be retrieved by Tally as in some of the other answers.