【动态规划】MATLAB和Python实现-Part02_matlab动态规划

往期系列：
【动态规划】MATLAB和Python实现-Part01

---

一、从斐波那契数列到动态规划

1.1 回顾与引入

斐波那契数列的例子严格来说不是动态规划，因为动态规划一般用来求解优化问题，但在这个例子中我们实际上能看到动态规划的影子 ~

---

在 《算法导论》 中，这样描述动态规划：
动态规划通过组合子问题的解来求解原问题，一般来说，动态规划应用于 重叠子问题 的情况，即不同的子问题具有公共的子问题。
动态规划算法对每个子问题只求解一次，将其解保存在一个表格中，从而无需每次求解一个子问题时都重新计算，避免了这种不必要的计算。

---

实际上，动态规划有两种等价的实现方法，也即我们上一篇提到的两种方法：

1. 带有备忘录的递归算法；
2. 自底向上法。

实际中我们常用第二种方法。

1.2 动态规划中常见的概念

为了便于理解概念，我们仍然以斐波那契数列的例子进行举例说明，即使它不算严格的动态规划问题。

原问题和子问题
原问题就是你要求解的这个问题本身，子问题是和原问题相似但规模较小的问题。
（原问题本身就是子问题的最复杂的情形，即子问题的特例）
例如： 要求 $F(n)$，那么 $F(n)$ 就是原问题，$F(k),k\le n$ 就是子问题。

状态
状态就是子问题中会变化的某个量，可以把状态看成要求解问题的自变量。
例如： 我们要求解 $F(n)$，其中 $n$ 就是一个状态。

状态转移方程
能够表示状态之间转移关系的方程，一般利用关于状态的某个函数建立起来。
例如： $F(n)=F(n-1)+F(n-2),n\ge 2$为一个状态转移方程。

边界条件
一般为简单的初始条件，表现为程序的出口。
例如： $F(0)=0,F(1)=1$ 就称为边界条件。

DP数组（DP=动态规划）
DP数组也可以叫做“子问题数组”，因为DP数组中的每一个元素都对应一个子问题的结果，DP数组的下标一般就是该子问题对应的状态。
例如： 当用带有备忘录的递归算法时，我们定义向量（或列表）memo就可以看作为一个DP数组，数组下标从1到n+1（或列表索引从0到n），对应元素从 $F(0)$ 到 $F(n)$。

二、动态规划初体验

2.1 题目介绍——打家劫舍

题目描述：
你是一个小偷，现在有一排相邻的房子等着你去偷窃。这些房子装有相互连通的防盗系统，如果两间相邻的房子在同一晚上被小偷闯入，系统会自动报警。给定一个代表每个房子存放金额的正整数数组 ，计算你不出动警报装置的情况下，一夜之内能够偷窃到的最高金额。（不用考虑偷窃时间）

示例1：
输入：[1, 2, 3, 1, 3]
输出：7

解释：偷窃1，3，5号房子，可以获得最高金额为1+3+3=7

示例2：
输入：[2, 7, 2, 3, 8]
输出：15

解释：偷窃2，5号房子，可以获得最高金额为7+8=15

2.2 问题的分析

假设一共有 $n$ 间房子，第 $i$ 间房子中有金额为 $M_i$ 的钱财。
步骤一：定义原问题和子问题
原问题： 从全部 $n$ 间房子中能够偷到的最大金额；
子问题： 从前 $k(k\le n)$个房子中能够偷到的最大金额。

步骤二：定义状态
我们记 $f(k)$ 为小偷能从前 $n$ 间房子中偷到的最大金额，这里的 $k$ 就是子问题的一个状态，当 $k=n$ 时，$f(n)$ 就是原问题的解。

前提条件：不能偷两个相邻的房子。

步骤三：寻找状态转移方程
所谓的寻找状态转移方程，本质上就是寻找子问题之间的一个递推关系，具体分析如下：

• 当 $k=1$ 时：
  只有 1 间房子可偷，那么就有 $f(1)=M_1$

• 当 $k=2$ 时：
  有 2 间房子可偷，而且不能偷相邻的房子，那就选择金额较大的房子偷，即： f ( 2 ) = m a x { M 1 , M 2 } f(2)=max\{M_{1},M_{2}\} f(2)=max{M1​,M2​}

• 当 $k\ge 3$ 时：
  此时我们根据第 $k-1$ 个房子是否被偷的状态，分为两种状态讨论。
  当第 $k-1$ 间房子 已经被偷，那就不能再偷第 $k$ 间房子，即能获得的最高金额为 $f(k-1)$.
  当第 $k-1$ 间房子 还没有被偷，那就相当于只偷了前 $k-2$ 间房子，那就一定会偷第 $k$ 间房子，即能获得的最高金额为 $f(k-2)+M_k$.
  那么，采取哪种情况取决于两种情况下获取的最高金额，即 f ( k ) = m a x { f ( k − 1 ) , f ( k − 2 ) + M k } f(k)=max\{f(k-1),f(k-2)+M_{k}\} f(k)=max{f(k−1),f(k−2)+Mk​}

综上，最终确定状态转移方程为：
f ( k ) = { M 1 k = 1 m a x { M 1 , M 2 } k = 2 m a x { f ( k − 1 ) , f ( k − 2 ) + M k } k ≥ 3 f(k)=

$$\begin{cases} M_{1} & k=1 \\ max\{M_{1},M_{2}\} & k=2 \\ max\{f(k-1),f(k-2)+M_{k}\} & k\geq 3 \end{cases}$$ f(k)=⎩⎪⎨⎪⎧​M1​max{M1​,M2​}max{f(k−1),f(k−2)+Mk​}​k=1k=2k≥3​

2.3 问题的求解

假设房子中钱财金额依次为：
7，5，8，10，6，9，11，6，8，13，15，12，4，3，9

递归求解
MATLAB求解：

%% 函数
function f = djjs_dg(M)
    k = length(M);
    if k == 1
        f = M(1);
    elseif k == 2
        f = max(M(1),M(2));
    else
        M_1 = M(1:k-1);
        M_2 = M(1:k-2);
        f = max(djjs_dg(M_1), djjs_dg(M_2)+M(k));
end
1
2
3
4
5
6
7
8
9
10
11
12

%% 调用
clear;clc;
M = [7, 5, 8, 10, 6, 9, 11, 6, 8, 13, 15, 12, 4, 3, 9];
tic;
djjs_dg(M)
toc;
1
2
3
4
5
6

结果：

---

Python求解：

# 打家劫舍——递归算法
import time
def djjs_dg(M):
    k = len(M)
    if k==1:
        f = M[0]
    elif k==2:
        f = max(M[0], M[1])
    else:
        M_1 = M[0:-1]
        M_2 = M[0:-2]
        f = max(djjs_dg(M_1), (djjs_dg(M_2)+M[-1]))
    return f

M = [7, 5, 8, 10, 6, 9, 11, 6, 8, 13, 15, 12, 4, 3, 9]
start = time.time()
print('打家劫舍递归算法得到的最大金额为：')
print(djjs_dg(M))
end = time.time()
print('打家劫舍递归算法用时(s)：')
print('%.8f' % (end-start))
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21

结果：

动态规划求解
MATLAB求解：

%% 函数
function f = djjs_dp(M)
    n = length(M);
    if n == 1
        f = M(1);
    elseif n == 2
        f = max(M(1),M(2));
    else
        result = -1*ones(1, n);
        result(1) = M(1);
        result(2) = max(M(1),M(2));
        for i = 3:n
            result(i) = max(result(i-1), result(i-2)+M(i));
        end
        f = result(i);
    end
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

%% 调用
clear;clc;
M = [7, 5, 8, 10, 6, 9, 11, 6, 8, 13, 15, 12, 4, 3, 9];
tic;
djjs_dp(M)
toc;
1
2
3
4
5
6

结果：

---

Python求解：

# 打家劫舍——动态规划算法
import time
def djjs_dp(M):
    n = len(M)
    if n==1:
        f = M[0]
    elif n==2:
        f = max(M[0],M[1])
    else:
        result = [-1 for i in range(n)]
        result[0] = M[0]
        result[1] = max(M[0], M[1])
        for i in range(2, n):
            result[i] = max(result[i-1], (result[i-2]+M[i]))
        f = result[-1]
    return f

M = [7, 5, 8, 10, 6, 9, 11, 6, 8, 13, 15, 12, 4, 3, 9]
start = time.time()
print('打家劫舍动态规划算法得到的最大金额为：')
print(djjs_dp(M))
end = time.time()
print('打家劫舍动态规划算法用时(s)：')
print('%.8f' % (end-start))
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24

结果：

---

2.4 动态规划算法优化——状态压缩

状态压缩 就是用来减少空间复杂度的一种方法，其基本思想是很多时候我们并不需要始终持有全部的DP数组，而是仅需要其中的一部分。

比如：对于 打家劫舍 问题，我们在计算 $f(n)$ 时，实际上只用了 $f(n-1)$ 和 $f(n-2)$ 的结果，再往前的子问题的结果实际上已经用不到了。
因此我们可以利用两个变量保存前两个子问题的结果，这样既可以依次计算出所有子问题结果，又能节省计算机的空间。

根据这个思路，我们重新设计MATLAB 和 Python 的函数：

MATLAB：

function f = djjs_dp(M)
    n = length(M);
    if n == 1
        f = M(1);
    elseif n == 2
        f = max(M(1),M(2));
    else
        pre2 = M(1);
        pre1 = max(M(1),M(2));
        for i = 3:n
            cur = max(pre1, pre2+M(i));
            pre2 = pre1;
            pre1 = cur;
        end
        f = cur;
    end
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

Python：

def djjs_dp(M):
    n = len(M)
    if n==1:
        f = M[0]
    elif n==2:
        f = max(M[0],M[1])
    else:
        pre2 = M[0]
        pre1 = max(M[0], M[1])
        for i in range(2, n):
            cur = max(pre1, pre2+M[i])
            pre2 = pre1
            pre1 = cur
        f = cur
    return f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

---

2.5 输出偷窃房子的编号

首先，我们要注意：相同的金额可能对应着不同的偷窃方案。
比如：当金额为[1, 3, 2]时，有 1+2=3，即选第一个和第三个 等同于 只选第二个。
我们的程序只需要输出其中一种结果就可以了。

算法原理
给定了每间房子的金额向量（或列表）为：M；
通过之前的算法，我们可以得到 DP数组 为：result，其中第 k 个元素（起始坐标为1）的含义是偷窃前 k 间房子可以得到的最大金额。

我们可以容易得出以下 2个结论：

1. 对于任意的 i，都有 result(i)>=M(i)，当且仅当只偷窃一间房子时，等号成立；
2. 随着 i 增大，result(i) 单调递增。

接下来，我们需要通过 result来推出对应的编号数组 IND（IND 初始为空）：

• 首先，找到 result 中最大值第一次出现的位置，这个位置就是偷窃的最后一个房子。记这个位置下标为 ind，并将其添加入 IND 中。

反证：如果不是最后一个偷窃的房子，那么必然不是最大值。

• 之后，判断 result(ind) 和 M(ind) 的大小：因为根据上述结论，必定有 result(i)>=M(i)，因此我们可以分两种情况分别讨论：
• • 若 result(ind)==M(ind)，则说明只偷窃了这一间房子，那么可以直接返回 IND；
• • 若 result(ind)>M(ind)，则说明之前还偷窃了其他房子，那么我们就计算两者的差 diff=result(ind)-M(6)，表示前 ind-1 间房子中偷窃的金额；之后在 result 中找到第一次出现 diff 的位置，有 result(new_ind)=diff，因此将 ind=new_ind 添加进 IND 中；

---

为了更好地理解上述过程，我们举一个例子：
给定：M=[1, 4, 1, 2, 5, 6, 3]；
得到：result=[1, 4, 4, 6, 9, 12, 12]
步骤：

• 步骤1，找到 result 中的最大值 12，其第一次出现的下标为 ind=6，将其添加进 IND，得到 IND=[6]
• 步骤2，比较 result(6)=12 和 M(6)=6 的大小，发现 result(6)>M(6)，说明除了第 6 间房子，还偷了其他房子，计算 diff=result(6)-M(6)，得到 diff=6，表示在前 6-1=5 间房子中偷窃到的金额；之后在 result 中找到第一次出现 diff=6 的位置，有 result(4)=6，因此将 ind=4 添加进 IND 中；
• 此后，不断重复步骤2，直到 result(ind)==M(ind) 为止。

---

接下来，我们用 MATLAB 和 Python 分别实现：
MATLAB：
函数实现：

function [f, IND] = djjs_with_ind(M)
    n = length(M);
    if n == 1
        f = M(1);
    elseif n == 2
        [f, IND] = max(M(1), M(2));
    else
        result = zeros(1, n);
        result(1) = M(1);
        result(2) = max(M(1), M(2));
        for i = 3:n
            result(i) = max(result(i-1), result(i-2)+M(i));
        end
        f = result(i);
        IND = [];
        ind = find(result==result(end), 1);
        IND = [IND, ind];
        while result(ind) > M(ind)
            ind = find(result==(result(ind)-M(ind)), 1);
            IND = [IND, ind];
        end
        IND = IND(end:-1:1);
    end
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24

函数调用：

clear;clc;
M = [7, 5, 8, 10, 6, 9, 11, 6, 8, 13, 15, 12, 4, 3, 9];
tic;
[f, IND] = djjs_with_ind(M)
toc;
1
2
3
4
5

结果：

---

Python：
函数实现和调用：

# 打家劫舍——递归算法，并输出偷窃的房子编号
import time
def djjs_with_ind(M):
    n = len(M)
    if n==1:
        f = M[0]
        IND = 0
    elif n==2:
        f = max(M[0],M[1])
        IND = M.index(f)
    else:
        result = [-1 for i in range(n)]
        result[0] = M[0]
        result[1] = max(M[0], M[1])
        for i in range(2, n):
            result[i] = max(result[i-1], (result[i-2]+M[i]))
        f = result[-1]
        # 推算IND
        IND = []
        ind = result.index(f)
        IND.append(ind)
        while result[ind] > M[ind]:
            ind = result.index((result[ind]-M[ind]))
            IND.append(ind)
        IND = IND[::-1]             # 逆序
        IND = [i+1 for i in IND]    # 索引都+1，从1开始
    return f, IND

M = [7, 5, 8, 10, 6, 9, 11, 6, 8, 13, 15, 12, 4, 3, 9]
start = time.time()
f, IND = djjs_with_ind(M)
print('打家劫舍动态规划算法得到的最大金额为：')
print(f)
print('偷窃房子的编号分别为：')
for i in IND:
    print(i, end=' ')
end = time.time()
print()
print('打家劫舍动态规划算法用时(s)：')
print('%.8f' % (end-start))
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

结果：