Finding maximum size sub-matrix of all 1's in a matrix having 1's and 0's

Suppose you are given an mXn bitmap, represented by an array M[1..m,1.. n] whose entries are all 0 or 1. A all-one block is a subarray of the form M[i .. i0, j .. j0] in which every bit is equal to 1. Describe and analyze an efficient algorithm to find an all-one b开发者_如何转开发lock in M with maximum area

I am trying to make a dynamic programming solution. But my recursive algorithm runs in O(n^n) time, and even after memoization I cannot think of bringing it down below O(n^4). Can someone help me find a more efficient solution?

An O(N) (number of elements) solution:

A
1 1 0 0 1 0
0 1 1 1 1 1
1 1 1 1 1 0
0 0 1 1 0 0 

Generate an array C where each element represents the number of 1s above and including it, up until the first 0.

C
1 1 0 0 1 0
0 2 1 1 2 1
1 3 2 2 3 0
0 0 3 3 0 0 

We want to find the row R, and left, right indices l , r that maximizes (r-l+1)*min(C[R][l..r]). Here is an algorithm to inspect each row in O(cols) time:

Maintain a stack of pairs (h, i), where C[R][i-1] < h ≤ C[R][i]. At any position cur, we should have h=min(C[R][i..cur]) for all pairs (h, i) on the stack.

For each element:

• If h_cur>h_top
• • Push (h, i).
• Else:
• • While h_cur<h_top:
• • • Pop the top of the stack.
• • • Check whether it would make a new best, i.e. (i_cur-i_pop)*h_pop > best.
• • If h_cur>h_top
• • • Push (h, i_lastpopped).

An example of this in execution for the third row in our example:

  i =0      1      2      3      4      5
C[i]=1      3      2      2      3      0
                                 (3, 4)
 S=         (3, 1) (2, 1) (2, 1) (2, 1)
     (1, 0) (1, 0) (1, 0) (1, 0) (1, 0)    
     (0,-1) (0,-1) (0,-1) (0,-1) (0,-1) (0,-1) 

i=0, C[i]=1) Push (1, 0).
i=1, C[i]=3) Push (3, 1).
i=2, C[i]=2) Pop (3, 1). Check whether (2-1)*3=3 is a new best.
        The last i popped was 1, so push (2, 1).
i=3, C[i]=2) h_cur=h_top so do nothing.
i=4, C[i]=3) Push (3, 4).
i=5, C[i]=0) Pop (3, 4). Check whether (5-4)*3=3 is a new best.
        Pop (2, 1). Check whether (5-1)*2=8 is a new best.
        Pop (1, 0). Check whether (5-0)*1=5 is a new best.
        End. (Okay, we should probably add an extra term C[cols]=0 on the end for good measure).

Here's an O(numCols*numLines^2) algorithm. Let S[i][j] = sum of the first i elements of column j.

I will work the algorithm on this example:

M
1 1 0 0 1 0
0 1 1 1 0 1
1 1 1 1 0 0
0 0 1 1 0 0 

We have:

S
1 1 0 0 1 0
1 2 1 1 1 1
2 3 2 2 1 1 
2 3 3 3 1 1

Now consider the problem of finding the maximum subarray of all ones in a one-dimensional array. This can be solved using this simple algorithm:

append 0 to the end of your array
max = 0, temp = 0
for i = 1 to array.size do
  if array[i] = 1 then
    ++temp
  else
    if temp > max then
      max = temp
    temp = 0

For example, if you have this 1d array:

1 2 3 4 5 6
1 1 0 1 1 1

you'd do this:

First append a 0:

1 2 3 4 5 6 7
1 1 0 1 1 1 0

Now, notice that whenever you hit a 0, you know where a sequence of contiguous ones ends. Therefore, if you keep a running total (temp variable) of the current number of ones, you can compare that total with the maximum so far (max variable) when you hit a zero, and then reset the running total. This will give you the maximum length of a contiguous sequence of ones in the variable max.

Now you can use this subalgorithm to find the solution for your problem. First of all append a 0 column to your matrix. Then compute S.

Then:

max = 0
for i = 1 to M.numLines do
  for j = i to M.numLines do
    temp = 0
    for k = 1 to M.numCols do
      if S[j][k] - S[i-1][k] = j - i + 1 then
        temp += j - i + 1
      else
        if temp > max then
          max = temp
        temp = 0

Basically, for each possible height of a subarray (there are O(numLines^2) possible heights), you find the one with maximum area having that height by applying the algorithm for the one-dimensional array (in O(numCols)).

Consider the following "picture":

   M
   1 1 0 0 1 0 0
i  0 1 1 1 0 1 0
j  1 1 1 1 0 0 0
   0 0 1 1 0 0 0

This means that we have the height j - i + 1 fixed. Now, take all the elements of the matrix that are between i and j inclusively:

0 1 1 1 0 1 0
1 1 1 1 0 0 0

Notice that this resembles the one-dimensional problem. Let's sum the columns and see what we get:

1 2 2 2 0 1 0

Now, the problem is reduced to the one-dimensional case, with the exception that we must find a subsequence of contiguous j - i + 1 (which is 2 in this case) values. This means that each column in our j - i + 1 "window" must be full of ones. We can check for this efficiently by using the S matrix.

To understand how S works, consider a one-dimensional case again: let s[i] = sum of the first i elements of the vector a. Then what is the sum of the subsequence a[i..j]? It's the sum of all the elements up to and including a[j], minus the sum of all those up to and including a[i-1], meaning s[j] - s[i-1]. The 2d case works the same, except we have an s for each column.

I hope this is clear, if you have any more questions please ask.

I don't know if this fits your needs, but I think there's also an O(numLines*numCols) algorithm, based on dynamic programming. I can't figure it out yet, except for the case where the subarray you're after is square. Someone might have better insight however, so wait a bit more.

Define a new matrix A wich will store in A[i,j] two values: the width and the height of the largest submatrix with the left upper corner at i,j, fill this matrix starting from the bottom right corner, by rows bottom to top. You'll find four cases: Perform these cases when given matrix at [i,j]=1

case 1: none of the right or bottom neighbour elements in the original matrix are equal to the current one, i.e: M[i,j] != M[i+1,j] and M[i,j] != M[i,j+1] being M the original matrix, in this case, the value of A[i,j] is 1x1

case 2: the neighbour element to the right is equal to the current one but the bottom one is different, the value of A[i,j].width is A[i+1,j].width+1 and A[i,j].height=1

case 3: the neighbour element to the bottom is equal but the right one is different, A[i,j].width=1, A[i,j].height=A[i,j+1].height+1

case 4: both neighbours are equal: Three rectangles are considered:

1. A[i,j].width=A[i,j+1].width+1; A[i,j].height=1;

2. A[i,j].height=A[i+1,j].height+1; a[i,j].width=1;

3. A[i,j].width = min(A[i+1,j].width+1,A[i,j+1].width) and A[i,j].height = min(A[i,j+1]+1,A[i+1,j])

The one with the max area in the above three cases will be considered to represent the rectangle at this position.

The size of the largest matrix that has the upper left corner at i,j is A[i,j].width*A[i,j].height so you can update the max value found while calculating the A[i,j]

the bottom row and the rightmost column elements are treated as if their neighbours to the bottom and to the right respectively are different.

Here is a O(N) implementation in C#.

The idea is to use a dynamic programming to build an accumulated Matrix that has the size of the biggest submatrix including the current cell itself.

    public static int LargestSquareMatrixOfOne(int[,] original_mat)
    {
        int[,] AccumulatedMatrix = new int[original_mat.GetLength(0), original_mat.GetLength(1)];
        AccumulatedMatrix[0, 0] = original_mat[0, 0];
        int biggestSize = 1;
        for (int i = 0; i < original_mat.GetLength(0); i++)
        {
            for (int j = 0; j < original_mat.GetLength(1); j++)
            {
                if (i > 0 && j > 0)
                {
                    if (original_mat[i, j] == 1)
                    {
                        AccumulatedMatrix[i, j] = Math.Min(AccumulatedMatrix[i - 1, j - 1], (Math.Min(AccumulatedMatrix[i - 1, j], AccumulatedMatrix[i, j - 1]))) + 1;
                        if (AccumulatedMatrix[i, j] > biggestSize)
                        {
                            biggestSize = AccumulatedMatrix[i, j];
                        }
                    }
                    else
                    {
                        AccumulatedMatrix[i, j] = 0;
                    }
                }
                else if ( (i > 0 && j == 0) || (j > 0 && i == 0))
                {
                    if (original_mat[i, j] == 1) { AccumulatedMatrix[i, j] = 1; }
                    else { AccumulatedMatrix[i, j] = 0; }
                }
            }
        }
        return biggestSize;
    }