具有容差的两个数据集的最长公共子序列答案

【问题标题】：Longest common subsequence of two datasets with tolerance具有容差的两个数据集的最长公共子序列
【发布时间】：2020-07-31 13:02:57
【问题描述】：

我有四个一维np.arrays：x1, y1, x2, y2，其中x1和y2具有相同的长度，x2和y2也具有相同的长度，因为它们是对应的x和y值对于一个数据集。 len(x1) 和 len(x2) 总是不同的。我们现在假设len(x1) > len(x2)。这两个数组总是有共同的值，但是以一种特殊的方式：这些值并不完全相同，只是在一个公差范围内（因为数值错误等）。容差 = 0.01 的示例：

x1 = np.array([0, 1.01, 1.09, 1.53, -9.001, 1.2, -52, 1.011])
x2 = np.array([1, 1.1, 1.2, 1.5, -9, 82])

我只想保留公共值（以容差方式）。使用较短的数组作为参考，在本例中为 x2。 x2中的第一个值是1，在x1中有一个对应的值，即1.01。 Next：1.2在x2、1.2中也有对应的值。 1.5的值没有对应的值，因为1.53超差，所以过滤掉等等。完整的结果应该是：

x1 = np.array([1.01, 1.09, -9.001, 1.2])
x2 = np.array([1, 1.1, -9, 1.2])

为了更进一步，基于以这种方式过滤 x 值，我想过滤两个数据集的相同索引的 y 值，换句话说，我想找到两个数据集的最长公共子序列。请注意，这里的排序很重要，因为它与 y 值相关（如果我们先对 x 进行 argsort 排序并重新索引 x 和 y 并不重要）。

我根据this的答案尝试了什么：

def longest_common_subseq(x1, x2, y1, y2, tol=0.02):
    # sort them first to keep x and y connected
    idx1 = np.argsort(x1)
    x1, y1 = x1[idx1], y1[idx1]
    idx2 = np.argsort(x2)
    x2, y2 = x2[idx2], y2[idx2]
    
    # here I assumed that len(x2) < len(x1)
    idx = (np.abs(x1[:,None] - x2) <= tol).any(axis=1)
    
    return x1[idx], x2[idx], y1[idx], y2[idx]

在这种情况下，y 值可以是任意的，只有形状必须与x1 和x2 匹配。例如：

y1 = np.array([0, 1, 2, 3, 4, 5, 6, 7])
y2 = np.array([-1, 0, 3, 7, 11, -2])

尝试运行上面的函数引发

IndexError: boolean index did not match indexed array along dimension 0.

我明白：索引数组的长度是错误的，因为x1 和x2 的长度不同，到目前为止我做不到。有没有很好的方法来实现这一点？

编辑：

如果多个值在公差范围内，则应选择最接近的值。

【问题讨论】：

任务定义不够好。如果x1=[1, 1.02, 1.03] 和x2=[1.021, 1.017] 会发生什么？多个答案在tol=0.02 之内，x1 中的相同答案最接近x2 中的不同元素。什么是正确的配对？
@Aguy 这几乎肯定不会发生，因为x1 和x2 原本是相同的数组，其中一个只是经过切片和某些操作。您提到的正确答案是x1=[1.02]，x2=[1.021]。我们使用较短的数组作为参考，遍历它的元素，查看x2 中是否有任何对应的值对应所需的tol，然后选择最接近的。一旦我们选择了一个值，我们需要在之后排除它以避免重复。我知道这听起来很复杂，我同意它的定义不够好，也许我稍后会改写这个问题。

标签： python numpy

【解决方案1】：

一种简单的方法是找出所有元素之间的距离：

dist = np.abs(x1 - x2[:, None])

既然你说通常你不会有多个元素在任何其他元素的公差范围内，你可以这样做

i2, i1 = np.nonzero(dist < tol)

如果你有多个匹配，你可以先修剪匹配：

i1 = np.argmin(dist, axis=1)
i2 = np.flatnonzero(dist[np.arange(x2.size), i1] < tol)
i1 = i1[i2]

如果原始数据已排序，则索引也会排序（它们将是对角线）。这意味着您可以通过检查索引之间的间距来检查子序列长度。匹配序列的两个索引都将递增 1。

mask = (np.diff(i1) == 1) & (np.diff(i2) == 1)
# smear the mask to include both endpoints
mask = np.r_[False, mask] | np.r_[mask, False]
# pad the mask to ensure proper indexing and find the changeover points
locs = np.diff(np.r_[False, mask, False])
inds = np.flatnonzero(locs)
lengths = inds[1::2] - inds[::2]

您可以从上面的数量中找到最长运行的指数：

k = np.argmax(lengths)
start = inds[2 * k]
stop = inds[2 * k + 1]
longest_x1 = x1[i1[start:stop]]
longest_y1 = y1[i1[start:stop]]
longest_x2 = x2[i2[start:stop]]
longest_y2 = y2[i2[start:stop]]

【讨论】：

它（几乎）有效。行i2 = np.flatnonzero(dist[i1, np.arange(x2.size)] < tol) 引发IndexError，但是如果我把那部分放在外面，它会返回正确的结果。知道为什么会这样吗？
@PéterLeéh。我翻转了那个索引。现在修好了。如果您需要在多个相邻元素之间进行选择，则只需要这三行来代替第一个 sn-p。

【解决方案2】：

我认为这应该可以解决问题：

def longest_common_subseq(x1, x2, y1, y2, tol=0.02):
    # sort them first to keep x and y connected
    idx1 = np.argsort(x1)
    x1, y1 = x1[idx1], y1[idx1]
    idx2 = np.argsort(x2)
    x2, y2 = x2[idx2], y2[idx2]
    
    # here I assumed that len(x2) < len(x1)
    difference = np.abs(x1[:,None] - x2) <= tol
    no_multiples = difference.cumsum(axis=0).cumsum(axis=0) == 1
    out_idx1 = no_multiples.any(axis=1)
    out_idx2 = no_multiples.any(axis=0)
    return x1[out_idx1], x2[out_idx2], y1[out_idx1], y2[out_idx2]

分解，这段代码

difference = np.abs(x1[:,None] - x2) <= tol
no_multiples = difference.cumsum(axis=0).cumsum(axis=0) == 1
out_idx1 = no_multiples.any(axis=1)

与您上面的函数做同样的事情，但我使用了this post 中的 cumsum 技巧来消除容差内的多个值。

然后你需要另一个轴的第二组索引来避免IndexError。这就是这条线的作用

out_idx2 = no_multiples.any(axis=0)

【讨论】：

但是您的解决方案为tol=0.02 返回了正确的结果，请尝试使用tol=0.2：返回的数组具有不同的长度。现在我试图理解为什么会这样。