最小垂直切片和

Question

我想输入一个矩阵大小的N x N并剪切一个切片，这样每个元素就在它上面的元素的正下方、左下方或右下方。而成本是切片中所有元素的总和。我如何写一个程序来做这件事呢？

例如：矩阵以列表的形式给出

[[1,2,3],
 [4,5,6],
 [7,8,9]]

它包含以下切片：

(1,4,7), (1,4,8), (1,5,7), (1,5,8), (1,5,9), 
(2,4,7), (2,4,8), (2,5,7), (2,5,8), (2,5,9), (2,6,8), (2,6,9), 
(3,5,7), (3,5,8), (3,5,9), (3,6,8), (3,6,9)

那么最低权重的切片是(1,4,7)，其和为12。

Answer 1

正如vivek提到的，你可以用一个动态程序来解决这个问题：

创建一个与输入矩阵大小相同的成本表。成本矩阵的每个元素都存储在该元素处结束的最小切片的成本。如果您还将前面的切片元素存储在此成本表中，则还可以提取末尾的实际切片(而不仅仅是其成本)。

您可以很容易地初始化成本表。只需将输入矩阵的第一行复制到表中即可。然后，我们将逐行填充表的其余部分。设C为成本矩阵，M为输入矩阵。然后：

//Initialize cost table
for col = 0 to N - 1
    C(0, col) = M(0, col)
//Run dynamic program
for row = 1 to N - 1
    for col = 0 to N - 1
        //take the minimum of the three possible predecessors:
        //make sure that the entries exist (i.e., take care of the edges, not shown here)
        C(row, col) = M(row, col) 
                       + min(C(row - 1, col - 1)), C(row - 1, col), C(row - 1, col + 1))

在此之后，您只需在C的最后一行中找到最小值，这将为您提供最小切片的成本。要获得实际的切片，请遍历您在循环期间设置的前置指针(伪代码片段中未显示)。

Answer 2

我们可以将矩阵元素看作是slices中的顶点，看作是边。然后，问题可以表示为找到从任何顶行顶点到任何底行顶点的最短路径，其中每条边的权重等于连接的元素的值(除了连接第一行的边，它还具有第一行元素的权重)。

然后，我们可以使用例如Bellman-Ford algorithm在这些条件下找到最短路径。以下是一个示例实现：

import numpy as np


m, n = 10, 10
M = np.arange(m*n).reshape(m, n) + 1
for i in range(1, m):
    M[i:] = np.roll(M[i:], 1 if i <= m // 2 else -1, axis=1)
print('Matrix:')
print(M, end='\n\n')


def edges():
    for i in range(m - 1):
        yield [(i, 0), (i + 1, 0)]
        yield [(i, 0), (i + 1, 1)]
        for j in range(1, n - 1):
            yield [(i, j), (i + 1, j - 1)]
            yield [(i, j), (i + 1, j)]
            yield [(i, j), (i + 1, j + 1)]
        yield [(i, n - 1), (i + 1, n - 1)]
        yield [(i, n - 1), (i + 1, n - 2)]


def compute_path(start):
    distance = {index: np.inf for index in np.ndindex(m, n)}
    predecessor = {index: None for index in np.ndindex(m, n)}

    distance[start] = M[start]
    for __ in range(M.size - 1):
        for u, v in edges():
            weight = M[v]
            if distance[u] + weight < distance[v]:
                distance[v] = distance[u] + weight
                predecessor[v] = u
    stop = min(filter(lambda x: x[0] == n - 1, distance), key=lambda y: distance[y])
    path = [stop]
    while predecessor[path[-1]] is not None:
        path.append(predecessor[path[-1]])
    return path[::-1], distance[stop]


paths = [compute_path((0, c)) for c in range(n)]
opt = min(paths, key=lambda x: x[1])
print('Optimal path: {}, with weight: {}'.format(*opt))
print('Vertices: ', M[list(zip(*opt[0]))])

它作为输出给出：

Matrix:
[[  1   2   3   4   5   6   7   8   9  10]
 [ 20  11  12  13  14  15  16  17  18  19]
 [ 29  30  21  22  23  24  25  26  27  28]
 [ 38  39  40  31  32  33  34  35  36  37]
 [ 47  48  49  50  41  42  43  44  45  46]
 [ 56  57  58  59  60  51  52  53  54  55]
 [ 67  68  69  70  61  62  63  64  65  66]
 [ 78  79  80  71  72  73  74  75  76  77]
 [ 89  90  81  82  83  84  85  86  87  88]
 [100  91  92  93  94  95  96  97  98  99]]

Optimal path: [(0, 0), (1, 1), (2, 2), (3, 3), (4, 4), (5, 5), (6, 4), (7, 3), (8, 2), (9, 1)], with weight: 460
Vertices:  [ 1 11 21 31 41 51 61 71 81 91]

Answer 3

这个问题可以用graph theory表示，然后使用linear programming技术解决。

通过将矩阵的所有元素视为顶点，目标是找到一组在特定约束下最小化通过矩阵的加权路径的边(例如，如何构建切片)。

我们可以使用scipy.optimize.linprog (它使用Simplex algorithm)来解决线性规划问题c.T @ x。解决方案的每个元素表示N个节点中的任何节点之间的可能连接(即解决方案向量的大小为N**2)。系数c确定连接的权重:即连接到的节点的权重(矩阵元素的值)，除了第一层外，我们还需要添加起始权重。

为了得到一个有效的解决方案，我们还需要应用几个约束：

1. 结果必须恰好具有N - 1个边。
2. 每个节点都不能连接到自己。
3. 每行只能连接到最后一行以外的一个节点。
4. 每个节点最多只能连接到一个其他节点(最后一层的节点必须没有

节点)。对于第1层到第(N-1)层中的每一层，H118必须连接到下一层。H219>G220

这些都是我们必须装配在一起的相当多的部分，因此结果代码有点冗长，一开始可能看起来令人不知所措。然而，当仔细观察时，应该可以识别出单个部件和它们的结构(我试图尽可能多地进行评论，并添加了几个指纹)。下面是示例代码：

import numpy as np
from scipy.optimize import linprog


M = np.arange(9).reshape(3, 3) + 1
print('Matrix:')
print(M)
print('\n\n')

N = len(M)

# Compute all possible connections between nodes (1: possible, 0: forbidden).
pc = np.zeros(shape=(N**2, N**2), dtype=int)

# Connect to nodes below (except the last layer).
i = np.arange(N**2 - N)
pc[i, i + N] = 1

# Connect to left nodes (except the last layer and leftmost column).
pc[i, i + N - 1] = 1
pc[i[::N], i[::N] + N - 1] = 0

# Connect to left nodes (except the last layer and rightmost column).
r = i + N + 1
mask = r < N**2
pc[i[mask], r[mask]] = 1
r = r[N-1::N]
mask = mask[N-1::N]
pc[i[N-1::N][mask], r[mask]] = 0

print('Possible connections:')
print(pc)
print('\n\n')

# Coefficients for linear programming problem represent the weight of connections.
c = np.zeros(shape=(N**2, N**2), dtype=int)
# Add weights for connections.
c = np.tile(M.ravel(), (N**2, 1))
# Add additional weights for first layer.
c[:N] += M[0, :][:, None]
print('Coefficient matrix:')
print(c)
print('\n\n')

# === Add constraints ===
A_eq_1 = np.concatenate((
    # Exactly N-1 connections.
    np.ones(N ** 4, dtype=int)[None, :],
    # No node can connect to itself.
    np.diag([1] * N**2).flatten()[None, :]
), axis=0)
b_eq_1 = np.asarray([N - 1, 0], dtype=int)
print('Exactly N-1 connections and no self-connecting nodes:')
print(A_eq_1)
print(b_eq_1)
print('\n\n')

# Each layer connects to exactly one other node (except the last layer).
A_eq_2 = np.zeros((N, N ** 4), dtype=int)
for j in range(N):
    A_eq_2[j, j * N**3:(j + 1)*N**3] = 1
b_eq_2 = np.ones(N, dtype=int)
b_eq_2[-1] = 0
print('Each layer connects to exactly one other node (except the last layer):')
print(A_eq_2)
print(b_eq_2)
print('\n\n')

# Each node connects to at most one other node (except the ones in the last layer).
N = N ** 2
A_ub_1 = np.zeros((N, N ** 2), dtype=int)
for j in range(N):
    A_ub_1[j, j * N:j * N + N] = 1
b_ub_1 = np.ones(N, dtype=int)
b_ub_1[-1] = 0
print('Each node connects to at most one other node (except the ones in the last layer):')
print(A_ub_1)
print(b_ub_1)
print('\n\n')

# Each node respects its possible succesors (i.e. forbid all other connections).
A_eq_3 = np.zeros((N, N ** 2), dtype=int)
for j in range(N):
    A_eq_3[j, j * N:j * N + N] = 1 - pc[j, :]
b_eq_3 = np.zeros(N, dtype=int)
print('Each node respects its possible succesors (i.e. forbid all other connections):')
print(A_eq_3)
print(b_eq_3)
print('\n\n')

# For the layers 1 through (N-1) each node connected to must connect to the next layer.
A_eq_4 = np.zeros((N, N ** 2), dtype=int)
for j in range(len(M), N-len(M)):
    A_eq_4[j, j::N] = 1
    A_eq_4[j, j*N:(j+1)*N] = -1
b_eq_4 = np.zeros(N, dtype=int)
print('For the layers 1 through (N-1) each node connected to must connect to the next layer:')
print(A_eq_4)
print(b_eq_4)
print('\n\n')

# Concatenate all constraints.
A_eq = np.concatenate([A_eq_1, A_eq_2, A_eq_3, A_eq_4])
b_eq = np.concatenate([b_eq_1, b_eq_2, b_eq_3, b_eq_4])

A_ub = np.concatenate([A_ub_1])
b_ub = np.concatenate([b_ub_1])

res = linprog(c.ravel(), A_ub=A_ub, b_ub=b_ub, A_eq=A_eq, b_eq=b_eq, bounds=(0, 1))
print(res.success)
print(res.x.reshape(N, N))  # Edges.

最后一个输出是结果，它的形式如下：

[[0. 0. 0. 1. 0. 0. 0. 0. 0.]
 [0. 0. 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.]
 [0. 0. 0. 0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0. 0. 0. 0.]]

这告诉我们，为了获得最小路径，我们应该将节点0(行索引)连接到节点3(列索引)，以及将节点3(行索引)连接到节点6(列索引)。这表示路径(或“切片”) (1, 4, 7)。我们可以从第一行开始推导路径，然后作为边点遍历图形：

edges = res.x.reshape(N, N)
for i, r in enumerate(edges):
    # Numerical instabilities can cause some elements to have very small values > 0.
    if r.sum() > 0.5:
        print('Connect {} -> {}'.format(i, r.argmax()))

path = [edges[:len(M)].ravel().argmax() // N]
while edges[path[-1]].max() > 0.5:
    path.append(edges[path[-1]].argmax())
print('Path: ', path)
print('Elements: ', M.ravel()[path])
print('Path weight: ', M.ravel()[path].sum())

结束语

上面的示例代码为性能改进留下了很大的空间。例如，它将节点之间所有可能的连接视为可扩展为M.size**2的解决方案。尽管我们限制了可能的连接，但通过将其包含在问题体系结构中，计算的数量仍然比我们从一开始就限制它的数量多得多。这意味着，我们可以只使用2*(M.shape[0] - 1) + 3*(M.shape[1] - 2)*(M.shape[0] - 1)，而不是M.size**2 coefficients，它只可扩展为M.size。另外，我们可以使用更小的约束矩阵，因为我们已经在问题的体系结构中构建了这些约束。以上面的代码示例为基础，应该可以对其进行相应的调整。因此，我将在这里结束，并将可能的性能改进的实现留给感兴趣的读者。

(上面的实现也只适用于方阵，但对非方阵的推广应该很简单。)