LeetCode 解题:统计全 1 子矩阵

上周末试着参加了一下 LeetCode 的周赛。

文章是 7/10 写的,所以这里指的是 7/5 的周赛

90 分钟,总共四题,分别一题简单,两题中等,一题困难。在我看来,第 1、2 题都算简单,第 3 题的中等我很多年前做过;只有第 4 题最优解要用到比较少见的算法,不过还是有网友暴力解通过。

就是这样的难度,第 2 题一开始没反应过来,先做后面,然后在第 3 题耗费了太多时间(约 1H),留给 2 、4 的时间都不太够,最后只 AC 1、3。相比之下,最快的参加者 17 分钟全部 AC。自尊心受到暴击,只能用刚刚恢复刷题,太生疏给自己台阶。做的时候感觉很顺畅,回头一看时间,都不知道怎么跑掉的。是中年人反应迟钝了吗?

不管怎样,3、4 题还是值得说一下的。今天先讲第 3 题,对应 1504 号。

英文:https://leetcode.com/problems/count-submatrices-with-all-ones/

中文:https://leetcode-cn.com/problems/count-submatrices-with-all-ones/

题意

一句话概括

给定由 0, 1 组成的矩阵,求 元素全是 1 的子矩阵 的数量。

例子

1
2
3
4
5
6
7
8
9
10
11
12
13
14
输入:
1,0,1
1,1,0
1,1,0
输出:13
解释:
有 6 个 1x1 的矩形。
有 2 个 1x2 的矩形。
有 3 个 2x1 的矩形。
有 1 个 2x2 的矩形。
有 1 个 3x1 的矩形。
矩形数目总共 = 6 + 2 + 3 + 1 + 1 = 13

更多例子看测试用例

约束

  • 行列最小 1,最大 150
  • 矩阵元素只有 0 和 1

为了后续方便讨论,约定行数为 R,列数为 C。

测试代码先行

如何快速生成单元测试

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
package main
import "testing"
func Test_numSubmat(t *testing.T) {
tests := []struct {
name string
mat [][]int
want int
}{
{
name: "min0",
mat: [][]int{
{0},
},
want: 0,
},
{
name: "min1",
mat: [][]int{
{1},
},
want: 1,
},
{
name: "1",
mat: [][]int{
{1, 0, 1},
{1, 1, 0},
{1, 1, 0},
},
want: 13,
},
{
name: "2",
mat: [][]int{
{0, 1, 1, 0},
{0, 1, 1, 1},
{1, 1, 1, 0},
},
want: 24,
},
{
name: "3",
mat: [][]int{
{1, 1, 1, 1, 1, 1},
},
want: 21,
},
{
name: "4",
mat: [][]int{
{1, 0, 1},
{0, 1, 0},
{1, 0, 1},
},
want: 5,
},
// 现在 LeetCode 周赛一次错误提交只罚时 5 分钟,代价不算大
// 简单添加一些能马上想到的 case 即可,折腾超过 5 分钟不值
// 能通过简单的 case 可以先提交一次,不行再针对强化
}
for _, tt := range tests {
t.Run(tt.name, func(t *testing.T) {
if got := numSubmat(tt.mat); got != tt.want {
t.Errorf("numSubmat() = %v, want %v", got, tt.want)
}
})
}
}

动态规划

暴力做法非常简单,一句话带过:首先 O(RC) 遍历左上角选点,然后 O(RC) 遍历右下角选点,再用 O(RC) 判断选定的矩阵是否全 1,时间复杂度为 O((RC)^3) ,空间复杂度 O(1)。(提醒:R 是行数,C 是列数)也许数据弱的话能 AC(我没试),比赛时间不够也许可以蒙一下,但没有讨论价值。

这道题明显有大量重复子问题,稍有经验都应该马上往动态规划上靠。关于动态规划的基础常识,可以先看之前的文章:《LeetCode 题解: Regular Expression Matching》,以后考虑围绕基础知识点写一些文章。

想到 DP 容易,难点是定义子问题,状态如何转移。多年前那次和这次,都第一时间知道是 DP 的题,仍然耗费不少时间想具体的状态转移方程。(虽然做过,相似的矩阵题太多,过了一阵才认出来,又花了一点时间回想有哪些方案。)

思考方向

最明显的重复运算,是对候选矩阵做全 1 判断。能否记录某些中间状态,根据少数点是否为 1, 结合这些中间状态得到下一个状态,避免前面判断过的 位置 再来一遍?

根据以往做 DP 的经验,一般是 dp[r][c] 代表到 r 行 c 列为止的某个结果。例如,我们定义两个矩阵 width[r][c]height[r][c] 代表以 (r, c) 为右下角的矩阵的最大宽度和最大高度,从左上到右下运算。有没有可能从前面的状态,推算出以(r, c) 为右下角的极大矩阵呢?

绿色表示已经处理过的点,白色表示待处理的点。下同。

很难。我们找极大矩阵,是因为知道极大矩阵的大小,就很容易计算它有多少个子矩阵。但极大矩阵可以有很多个,光凭最大宽度和最大高度两个信息,无法计算出这么多个极大矩阵。那还要记录多少信息才够?

极大矩阵是指,这个矩阵在边界上碰到了 0,没办法变得更大了。但它不一定是最大的矩阵。

记录最大宽和最大高之所以很难用得上,是因为它们很可能不是在同一个矩阵上记录的,也很难用来计算那么多个极大矩阵。

换一个角度看为什么这么难。假定可以用 O(1) 完成状态转移,计算出 (r, c) 点为右下角的所有极大矩阵,由于右下角只有 O(RC) 种情况,换言之时间复杂度从 O((RC)^3) 降到了 O(RC),优化幅度非常大 。但多个极大矩阵真的能用 O(1) 时间就计算出来吗?要借用哪些中间结果?

事实上,这个方向是可行的。但我们先从简单一点的开始。

压缩 + DP

二维的情况复杂,是因为长宽有多种组合;如果长宽不能同时达到最大,就会产生大量极大矩阵。

那就 从一维开始找规律

一维问题:全 1 子段数

一维的情况,问题简化成『统计全 1 子段数』,其中任意长度不为 0 的连续区间 ,都是子段(interval)。假定数组长度为 N,Brute-Force 的做法是,O(N) 选左端点,O(N) 选右端点,O(N) 判断是否全 1,时间复杂度 O(N^3)。

这时用 DP,状态转移方程很好写。假定有一维数组 A,里面同样是 0 或 1。定义 L[i] 为 以第 i 个元素为右边界的子段的最大长度(换言之,连续多少个 1),很容易得到

$$
\left{
\begin{array}{lr}
Li = L{i-1} + 1, & (A_i = 1)\
L_i = 0 & (A_i = 0)
\end{array}
\right.
$$

大白话就是,遇到 1,最大长度加 1;遇到 0,长度清零

再定义 C[i] 为到 第 i 个元素为止,累计的的字段数,状态转移方程为

$$
\left{
\begin{array}{lr}
Ci = C{i-1} + L_i, & (A_i = 1)\
Ci = C{i-1} & (A_i = 0)
\end{array}
\right.
$$

大白话就是,遇到 1,就增加『当前子段最大长度 L[i] 这么多个子段数』;遇到 0,统计数字不变。为什么是 L[i] ?因为 以 A[i] 为右端点的子段刚好是 L[i] 个,分别为长度为 1、长度为 2、长度为 … 长度为 L[i]

DP 跟 BF 相比,省略了重复运算,但思维不严谨的话,很容易出现重复或者遗漏。我们来确认一下,是否每个左端点和右端点的组合都包含了一次:对于右端点,每个 A[i] = 1 ,都访问且只访问了一次;在确定右端点之后,每个和右端点连续(中间没有 0)的左端点都刚好统计了一次(L[i] 个)。

对于每个元素 A[i],分别需要常数时间计算 L[i] 和 C[i] ,时间复杂度为 O(N)。由于 L 和 C 都只和当前元素和前一个状态相关,只需保留最后一个状态,没有必要保存数组,空间复杂度为 O(1)。

代码非常简单

1
2
3
4
5
6
7
8
9
10
11
12
13
func numInterval(A []int) int {
length := 0
count := 0
for _, elem := range A {
if elem == 1 {
length++
count += length
} else {
length = 0
}
}
return count
}

二维问题:压缩成一维

一维的情况简单,我还是尽量详细地讲,因为这是理解接下来内容的基础。

对于高度为 2 的矩阵,求全 1 子矩阵数。这些子矩阵,高度要么为 1,要么为 2。对于高度为 1 的子矩阵,情况退回到 一维,分别对每一行应用一维的算法即可 。而对于高度为 2 的子矩阵,可以『压缩』成高度为 1,再应用一维的算法。

以此类推,任意高度的子矩阵,也可以压缩成一行处理 。而且,R 行压缩成一行,并不需要每次从头开始压缩;只需要在前 R - 1 行的结果上,再合并一行

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
func numSubmat(mat [][]int) int {
count := 0
rowNum := len(mat)
colNum := len(mat[0])
// 压缩的结果只需要保留一行进行运算
// 用 bool 值保存节省空间
mymat := make([]bool, colNum)
// 遍历开始压缩的行数
for rowSt := 0; rowSt < rowNum; rowSt++ {
for col := 0; col < colNum; col++ {
mymat[col] = true
}
// 遍历结束的行数
for rowEd := rowSt; rowEd < rowNum; rowEd++ {
length := 0
for col := 0; col < colNum; col++ {
if mat[rowEd][col] == 0 {
// 同一列上任意一行为0,这个位置就断开了
mymat[col] = false
}
// 压缩后,应用一维的算法
if mymat[col] {
length++
count += length
} else {
length = 0
}
}
}
}
return count
}

O(R) 遍历开始行,O(R) 遍历结束行,O(C) 的时间压缩每一行并应用一维的算法,时间复杂度为 O(R^2 * C)。只用了一行的额外空间保存压缩行,空间复杂度为 O(C)。由于 R 比 C 遍历多了一层,如果行数比列数大得多,可以将行列交换运算。代码比较麻烦,不再展开。

回到正确性的问题,是否不重不漏?横向在一维时已经讨论过了。纵向,O(R^2) 把所有开始行和结束行的组合都遍历了。换言之,所有的左上角和右下角的组合,都刚好地访问了一遍。

单调栈 + DP

上述解法已经相当优秀,将六次方的 BF 解法,降了 3 个阶。(这里为了方便讨论,把 R、C 看做相近的数量级,都当成 N 处理,相当于从 O(N^6) 降到了 O(N^3) 时间复杂度)

当年第一次遇到这道题,我甚至都没能一下子想到这种办法,想通时已经不够时间实现和调试。还是事后同事 Lucy 跟我们讲了思路。之后很长一段时间,没有再细想,把这当成最优解。

直到这次,顺便逛了一下讨论区,居然发现还能 O(RC) !?

以下解法并非我的原创,而是看了 LeetCode 美国站讨论区后,根据网友提供的思路,重新实现和用自己的话解释一遍。主要为了分享和加深理解。

平时做的答案,其核心思想,其实也是学来的(主要是大学算法课)。毕竟常见,不需要特别声明。

这次的解法比较绕,不敢掠美,特此声明。

子问题:分割极大矩阵

我们回到一开始想到用 动态规划 时提出的问题:『能不能遍历为 1 的点作为右下角,结合少量中间状态,算出新的中间状态,并逐渐推导出答案?

这个问题难在子问题和中间状态的定义。怎么定义,中间状态定义多少个,才能应付最多时跟 R、C 同量级的极大矩阵?

我们把问题简化,先不想怎么判断和划分,假定已经找到了一个 r 行 c 列的极大全 1 矩阵,里面有多少个 包含右下角的点 的子矩阵?

答:r x c 个。

包含右下角(最后一行结束),子矩阵高度(开始行)有 r 种可能;包含右下角(最后一列结束),矩阵宽度(开始列)有 c 种可能。组成起来有 r x c 种 子矩阵。换言之,只要找到以当前点为右下角的极大子矩阵,就能直接算出以此为右下角的全部子矩阵数量

难点在于极大子矩阵根本不止一个,叠在一起可能是锯齿状的,也可能是驼峰状的,极难记录和统计。不过只要用递归分治的思想划分一下子问题,难点就变简单了。

我们定义 C[i][j] 为 『以点 (i, j) 为右下角的子矩阵数量』。mat[i][j] 为 0 时,这个点不能构成任何子矩阵,所以只讨论为 1 的情况。如图所示,以第一个白色点 (4, 5) 为右下角的极大子矩阵,一共有 3 种宽高组合(1x6, 2x4, 4x2),它们重叠在一起,构成了一个不规则的形状。将其纵向划分成 3 个互不重叠的 3 个矩阵 X, Y, Z。(不一定是 3 个,实际可能更多,最多时可以分成每列一个)

那么 C[4][5] 就等于 X (包含右下角)的子矩阵数(4x2),加上 C[4][3]。之所以这样算,是因为所有以 (4, 3) 为右下角的子矩阵,都可以列数加 2 ,变成以 (4, 5) 为右下角。通过这样划分子问题,无论有多少个极大子矩阵,我们只需要计算包含当前右下角最右(同时也是最高)的一个极大子矩阵,再引用之前的子问题即可

为什么最右的极大子矩阵一定是最高的呢?因为子矩阵的高度无法超过右下角所在列的高度。上图中 列4 的高度为 5 (蓝色的(0, 4)点),但是这个点无法跟右下角构成矩阵。

引入 高度记录 和 高度栈

不过问题还没解决:

  • 如何知道最右的极大子矩阵的宽高?(4 和 2 怎么得到)
  • 如何知道上一个子问题是哪个?(怎么得知要引用 (4, 3) 的子问题)

后面这个问题好回答,要引用的子问题,就是当前子问题减掉当前极大子矩阵宽度(2)的子问题。问题变成了:当前极大子矩阵的宽高如何获取

为此,引入两个记录:

  • 高度记录 h[i][j] ,记录每列到当前点的最大高度(连续 1 的个数)。
  • 每行的高度栈,记录单调递增的高度序列。

先说高度记录,计算方法跟上一个解法的最大长度一样,只是变成纵向。有了它,极大子矩阵的宽高都能获得:高,就是当前列的高度(例子中为 4);宽则稍复杂一点,需要找到左边第一个比当前列低的列(例子中为 列 3),两列之间宽度就是极大子矩阵的宽度(5 - 3 = 2)。

如果不考虑时间复杂度,高度记录就足以解决问题。但是每次都要搜索左边的高度记录,导致每行的复杂度为 O(C^2) ,乘以行数后整体为 O(R*C^2),则跟第一种解法相比并无改进。所以引入高度栈,维护一个严格单调递增的高度序列。每当拿到新一列的高度,先将栈里大于等于当前列的高度出栈,剩下的栈顶的列,就是左边第一个比当前列低的列。完成之后将当前列进栈,继续下一步。

实际实践中,栈里一般储存列号,方便计算宽度,高度根据列号到 高度数组查询。

由于高度一旦下降,之前更高的列就无法与后续的点构成矩阵,所以出栈的高度再也不会用到。每个列都只会进栈一次,也(最多)只会出栈一次,遍历一行的时间复杂度就降低到了 O(C)。

如果这部分感觉有点迷糊,跟其它难点混在一起难以理解,可以先尝试做 LeetCode 84 ,在一维的情况下,把高度栈搞懂。

上代码:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
func numSubmat(mat [][]int) int {
count := 0
rowNum := len(mat)
colNum := len(mat[0])
// 高度记录只用到最后一行,保留一行即可
heights := make([]int, colNum)
// DP 解,以当前点为右下角的子矩阵数,同样只用到最后一行,只保留一行
dpCount := make([]int, colNum)
// 高度栈,注意保存的是列号,列高还得到 heights 取
hStack := newStack()
// 行循环
for row := 0; row < rowNum; row++ {
// 只有一个高度栈循环利用,每行开始前摇清空
hStack.clear()
// 列循环
for col := 0; col < colNum; col++ {
// 计算当前列高度
if mat[row][col] == 1 {
heights[col]++
} else {
heights[col] = 0
}
// 将栈里高度 大于等于当前列 的列出栈
for hStack.size() > 0 && heights[hStack.top()] >= heights[col] {
hStack.pop()
}
if hStack.size() == 0 {
// 一直到最左,没有列比当前列 低,宽度直接为 列号 + 1 (列号从 0 开始)
dpCount[col] = heights[col] * (col + 1)
} else {
// 找到左边第一个比当前列低的列,宽度为列号的差,并加上该列的子矩阵数
preCol := hStack.top()
dpCount[col] = dpCount[preCol] + heights[col]*(col-preCol)
}
// 当前列进栈
hStack.push(col)
}
// 结束本行前,将本行结果累计,下一循环就会被覆盖
for col := 0; col < colNum; col++ {
count += dpCount[col]
}
}
return count
}
// 基于切片自行实现了一个简单的栈
type stack struct {
data []int
}
func newStack() *stack { return &stack{data: make([]int, 0)} }
func (s *stack) size() int { return len(s.data) }
func (s *stack) push(v int) { s.data = append(s.data, v) }
func (s *stack) top() int { return s.data[len(s.data)-1] }
func (s *stack) pop() { s.data = s.data[:len(s.data)-1] }
func (s *stack) clear() { s.data = s.data[:0] }

正确性

  • 不重:右下角不同的子矩阵肯定不是同一个子矩阵。
  • 不漏:每一个点都作为右下角遍历了一次,每次把以该点为右下角的子矩阵都计算在内。

复杂度

  • 时间复杂度:最外层每行遍历一次,内层同时做两件事:

    • 每个点作为右下角计算一次,计算过程为常数;
    • 每个列下标都进栈一次,最多出栈一次。(如果一直单调递增就不出栈);

      整体的复杂度为 O(RC)。

  • 空间复杂度:常数组辅助变量,每组变量大小为 C;单调栈也最多进栈 C 个列下标,复杂度为 O(C)。


知识共享 “署名-非商业性使用-相同方式共享” 4.0 (CC BY-NC-SA 4.0)”许可协议
本文为本人原创,采用知识共享 “署名-非商业性使用-相同方式共享” 4.0 (CC BY-NC-SA 4.0)”许可协议进行许可。
本作品可自由复制、传播及基于本作品进行演绎创作。如有以上需要,请留言告知,在文章开头明显位置加上署名(Jayce Chant)、原链接及许可协议信息,并明确指出修改(如有),不得用于商业用途。谢谢合作。
请点击查看协议的中文摘要。