以经典例题理解经典算法：动态规划

一、从实例的角度理解动态规划

  什么是动态规划Dynamic Programming, DP？

  先来看看一种较为"正规"的解释：动态规划就是通过拆分问题，定义问题状态和状态之间的关系，使得问题能够以递推（或者说分治）的方式去解决。

  对上面的解释进行拆分，便不难得到动态规划的四个步骤：

• 拆分问题：能够通过动态规划解决的问题必然是可拆分的，这一步就是将一个复杂的问题拆分为多个可以逐步解决的子问题
• 定义状态：所谓状态其实就是每个子问题的解决方案，实战中我们需要找到合适的数据结构来存储状态
• 定义状态转移方程：拆分后的多个子问题其实是层层递进的，因此后一个子问题的解决一般都依赖于前一个子问题的解决方案，而所谓状态转移方程就是用来描述这个依赖关系的
• 进行状态推导：在定义好状态和状态转移方程之后，便可以通过进行状态的推导解决问题。这一步骤的关键是需要知道状态的初始值和状态推导的结束条件

  光看上面的文字描述其实还是很难弄清楚什么是动态规划，下面就结合一个动态规划的经典实例来尝试对其进行理解。

---

Example: （三角矩阵最长路径）假设有如下的三角矩阵：

data = [
	[6],
	[3, 8],
    [8, 1, 0],
	[2, 7, 3, 4]
]

现在要求找到从第一行到最后一行的最长路径长度，不要求找到具体的路径。

Solution: 认真观察不难发现三角矩阵从第一行到最后一行的路线其实是一种树状结构，如下图所示：

  怎么找到最长路径呢？粗暴一点的方法是遍历二叉树，找到最长的路径。但这种方法需要遍历所有的路径，不仅费时费力，而且还显得很low。

  怎么改善呢？我们先结合一点动态规划的思想。这个问题本质上就是求解一棵二叉树的最长路径，记为 $D(T)$。根据拆分问题的思想不难发现，其实二叉树的最长路径就等于其左子树最长路径和右子树最长路径中的最大值再加上这颗二叉树的根节点。用数学表示也就是：
$$D(T)=\max (D(T_{left}),D(T_{right}))+root$$

  仔细一看这不就是递归吗，于是乎，可以写一个solution如下：

def recursize_f(n_row, data, i, j):
    if i == n_row - 1:
        return data[i][j]
    else:
        return max(recursize_f(n_row, data, i + 1, j), recursize_f(n_row, data, i + 1, j + 1)) + data[i][j]

def solution1(n_row, data):
    return recursize_f(n_row, data, 0, 0)

计算结果显示最长路径长度为24。

  整个过程看起来似乎有动态规划内味了。问题已经拆分，所谓状态无非就是 $D(T)$ ，所谓状态转移方程也不就是：
D ( T ) = max ⁡ { D ( T l e f t ) ,   D ( T r i g h t ) } + r o o t D(T)=\max\{D(T_{left}),\ D(T_{right})\}+root D(T)=max{D(Tleft​), D(Tright​)}+root
那么下面就只剩下最后一个步骤：状态推导。

  首选确定状态的初始值，根据上面的递归过程不难发现初始状态就是二叉树的叶子结点。然后需要确定状态推导的结束条件，通过上面的问题拆分不难发现一共拆分出了三个级别的子问题。最低级别的子问题一共有2^(3-1)=4个，中间级别的子问题有2^(2-1)=2个，最高级别的子问题有2^(1-1)=1个，所以状态推导的结束条件就是遍历解决这4+2+1个子问题。

  确定了所有步骤，便可以实现该问题的动态规划solution，如下所示：

def solution2(n_row, data):
    assert (isinstance(data, list))
    dp = data[n_row - 1] # 初始状态
    # 遍历解决 4 + 2 + 1 个子问题
    for i in range(n_row - 2, -1, -1):
        dp_new = []
        for j in range(len(data[i])):
            assert (isinstance(dp, list))
            dp_new.append(max(dp[j], dp[j + 1]) + data[i][j])
        dp = dp_new # 更新状态
    return dp[0]

计算结果也是24。

  我们输入一个更复杂的矩阵，并统计两种solution的计算次数：

count = 0

def recursize_f(n_row, data, i, j):
    global count
    count += 1
    if i == n_row - 1:
        return data[i][j]
    else:
        return max(recursize_f(n_row, data, i + 1, j), recursize_f(n_row, data, i + 1, j + 1)) + data[i][j]

def solution1(n_row, data):
    return recursize_f(n_row, data, 0, 0)

def solution2(n_row, data):
    global count
    assert (isinstance(data, list))
    dp = data[len(data) - 1]
    for i in range(n_row - 2, -1, -1):
        dp_new = []
        for j in range(len(data[i])):
            count += 1
            assert (isinstance(dp, list))
            dp_new.append(max(dp[j], dp[j + 1]) + data[i][j])
        dp = dp_new
    return dp[0]

if __name__ == "__main__":
    data = [
        [6],
        [3, 8],
        [8, 1, 0],
        [2, 7, 3, 4],
        [6, 7, 3, 2, 4],
        [9, 2, 1, 7, 8, 5]
    ]
    print(solution1(6, data), count)
    count = 0
    print(solution2(6, data), count)

结果如下：

可见，当问题复杂时，动态规划能够大大降低计算量，性能远远优于常规的解法。

二、动态规划经典例题解析(持续更新)

2.1 最长上升子序列

问题描述：从一个一维的数字序列中找到一个最长的上升子序列，这个序列不要求是连续的
输入：[10, 9, 2, 5, 3, 7, 101, 18]
输出：4
解释：最长的上升子序列是[2, 3, 7, 101]，它的长度是4

  同样结合动态规划的思想，我们将原问题转换如下：记长度为 $m$ 的序列为 $seq$ ，序列中索引为 $n$ 的元素为 $s_n$ ，序列以 $s_n$ 结束的最长上升子序列长度为 $S(seq,n)$ ，求 $S(seq,m)$ 。那么不难得到如下状态转移方程：
S ( s e q ,   n ) = max ⁡ { S ( s e q ,   i )  if  s i < s n ,   i < n } + 1 S(seq,\ n)=\max\{S(seq,\ i)\ \text{if}\ s_i<s_n,\ i<n\}+1 S(seq, n)=max{S(seq, i) if si​<sn​, i<n}+1
  根据上面的状态转移方程，可以实现如下递归算法：

def recursive_f(seq, n):
    lengths = [0]
    for i in range(0, n):
        if seq[i] < seq[n]:
            lengths.append(recursive_f(seq, i))
    return 1 + max(lengths)

def solution1(seq):
    return recursive_f(seq, len(seq) - 1)

计算结果为4。

  由于以 $s_n$ 结束的最长上升子串长度至少为1，所以初始状态为1。从 $S(seq,1)$ 开始，当逐步推导到 $S(seq,m)$ 时结束。由此便可以将上面的递归算法转换为动态规划：

def solution2(seq):
    dp = [1 for s in seq]
    for i in range(1, len(seq)):
        lens = [0]
        for j in range(0, i):
            if seq[j] < seq[i]:
                lens.append(dp[j])
        dp[i] += max(lens)
    return max(dp)

  对比二者的运算量：

if __name__ == "__main__":
    seq = [10, 9, 2, 5, 4, 7, 101, 18]
    print(solution1(seq), count)
    count = 0
    print(solution2(seq), count)

结果如下：

2.2 乘积最大子数组

问题描述：从一个一维输出中找到乘积最大的子数组，输出这个子输出的乘积
输入：[2, 3, -2, 4]
输出：6
解释：乘积最大子数组为[2, 3]，乘积为6

  这题最直观的解法就是遍历所有子数组，然后通过比较得到乘积最大子数组的乘积。solution如下：

def solution1(seq):
    max_result = max(seq)
    for i in range(len(seq) - 1):
        result = seq[i]
        for j in range(i + 1, len(seq)):
            result *= seq[j]
            max_result = max(max_result, result)
    return max_result

  考虑动态规划的话，则需要先将原问题转换为如下形式：记长度为 $m$ 的序列为 $seq$ ，序列中索引为 $n$ 的元素为 $s_n$ ，序列以 $s_n$ 结束的乘积最大子数组乘积为 $S(seq,n)$ ，求 $S(seq,m)$ 。

  容易想到的是每次推导存储 $S(seq,n)$ 状态，那么状态转移方程如下：
S ( s e q ,   n ) = max ⁡ { s n S ( s e q ,   n − 1 ) ,   s n } S(seq,\ n)=\max\{s_nS(seq,\ n-1),\ s_n\} S(seq, n)=max{sn​S(seq, n−1), sn​}
初始状态就是序列 $seq$ 的第一个元素，推导结束条件则是将状态推导至 $S(seq,m)$ 。最终可得动态规划算法实现如下：

def solution2(seq):
    res = local_min = local_max = seq[0]
    for i in range(1, len(seq)):
        local_max = max(local_max * seq[i], seq[i])
        res = max(local_max, res)
    return res

计算结果为4，看似没有问题。

  我们修改一下输入，将序列修改为[2, 3, -2, 4, -1]，那么此时的乘积最大子序列乘积应当为2*3*-2*4*-1=48。然而上面动态规划算法的输出结果却仍然是6。稍微分析一下不难发现其原因，当状态推导至 $S(seq,3)$ 时，由于 $s_3=-2$ 为负数，所以 $S(seq,3)=S(seq,2)=6$ ，丢弃了-2，并导致后续的推导状态都等于 $S(seq,2)=6$ 。然而事实上，如果保留了-2，并在推导 $S(seq,5)$ 时乘上-1， $S(seq,5)$ 能够更大。总结一下其原因，就是上述算法忽略了最小值为负数时，乘上一个负数可能变成最大值的情况。

  为了解决这个问题，我们存储状态时不仅需要存储最大值 $S(seq,n{)}_{max}$ ，还需要存储最小值 $S(seq,n{)}_{min}$ 。那么新的状态转移方程如下：
{ S ( s e q ,   n ) m a x = max ⁡ { s n S ( s e q ,   n − 1 ) m a x ,   s n S ( s e q ,   n − 1 ) m i n ,   s n } S ( s e q ,   n ) m i n = min ⁡ { s n S ( s e q ,   n − 1 ) m a x ,   s n S ( s e q ,   n − 1 ) m i n ,   s n } \begin{cases} S(seq,\ n)_{max}=\max\{s_nS(seq,\ n-1)_{max},\ s_nS(seq,\ n-1)_{min},\ s_n\}\\ S(seq,\ n)_{min}=\min\{s_nS(seq,\ n-1)_{max},\ s_nS(seq,\ n-1)_{min},\ s_n\} \end{cases} {S(seq, n)max​=max{sn​S(seq, n−1)max​, sn​S(seq, n−1)min​, sn​}S(seq, n)min​=min{sn​S(seq, n−1)max​, sn​S(seq, n−1)min​, sn​}​
初始状态和推导结束条件不变，则可以重新实现动态规划算法如下：

def solution2(seq):
    res = local_min = local_max = seq[0]
    for i in range(1, len(seq)):
    	tmp = local_max
        local_max = max(tmp * seq[i], local_min * seq[i], seq[i])
        local_min = min(tmp * seq[i], local_min * seq[i], seq[i])
        res = max(local_max, res)
    return res

计算结果为48无误。

  对比solution1和动态规划算法的计算量：

if __name__ == "__main__":
    seq1 = [2, 3, -2, 4, -1]
    print("solution1", solution1(seq1), count)
    count = 0
    print("solution2", solution2(seq1), count)

结果如下：