numpy 与 valid

Question

我正在研究使用移动平均类型平滑操作对信号块进行卷积，但遇到影响下游计算的填充错误问题。

如果绘制出来，这是在块边界上引起的问题。

为了解决这个问题，我将前一个块中的 N 个样本添加到数组的 beginning 中。

问题是下面的卷积操作是否从原始数组的开头或结尾删除样本？

arr = np.array([0.9, 0.9, 0.9, 0.1, 0.1, 0.1, 1.0, 0.1, 0.1, 0.1, 0.9, 0.9, 0.9, ])
print(f"length input array {len(arr)}")
result = np.convolve(arr, [1]*3, 'valid')/3
print(f"length result array {len(result)}")
print(result)
plt.plot(result)

我设置了上面的示例来尝试确认将样本添加到数组的开头是正确的位置，但它没有帮助。输入长度为13，输出长度为11。

重叠时，卷积看起来像：

result = np.convolve(np.concatenate(previous_samples[-2:], arr), [1]*3, 'valid')/3

Answer 1

我假设（1）读者对什么是卷积有基本的了解，并且（2）使用具有奇数个元素的内核。

简短回答：

• 使用

  np.convolve(…, mode="valid")

  。

• 填充你的块

  • …从前一个块的

    (len(kernel) - 1) / 2

    最后一个元素开始，
  • ...在其末尾带有来自以下块的

    (len(kernel) - 1) / 2

    第一个元素，

  其中

  kernel

  指的是你的滑动内核。

长答案：

由于卷积可以被认为是具有滑动内核的操作，因此问题是，正如您所意识到的，边界处会发生什么。

numpy.convolve()

函数有 3 种模式来处理边界，如下所述：

• mode="valid"

  确保滑动内核永远不会离开信号。然而，这意味着结果的元素少于信号的元素；即

  len(kernel) - 1

  更少的元素。因此，如果您有一个包含 100 个元素的信号和一个包含 5 个元素的内核，那么您的结果将包含 100-(5-1)=96 个元素。操作在哪里“删除”缺失元素的问题可能有点误导；然而，考虑结果的一种方法是将其值写入滑动内核的中心，因此通过这种解释，可以对称地删除元素：开头的

  (len(kernel) - 1) / 2

  元素，结尾处的

  (len(kernel) - 1) / 2

  元素。信号结束。我们将在下面的代码示例中看到这是一个有意义的解释。

• mode="same"

  确保生成的信号与输入信号具有相同数量的值。然而，这意味着我们必须“弥补”值来抵消我们在

  mode="valid"

  中看到的元素较少的影响。因此，在这种模式下，信号在操作之前都会在开始和结束时用零填充 –

  (len(kernel) - 1) / 2

  。
最后，

• mode="full"

  确保产生信号和滑动内核之间重叠的所有潜在组合。这意味着（a）我们必须“弥补”更多的值，并且（b）结果将大于输入信号。这里的“组成”值再次为零，结果将包含

  len(signal) + len(kernel) - 1

  个元素，因此对于前面的示例来说，有 100+5-1=104 个元素。

如果您以块的形式处理信号，那么

mode="valid"

就是正确的选择：您可以确保生成的卷积块与输入块具有相同的长度，而无需“弥补”值（即不需要零） -padding），通过使用前一个块的尾随值和后一个块的前导值进行填充。对于除第一个和最后一个块之外的所有块都是如此，在这种情况下，您必须自己决定如何处理边界情况。在所有其他情况下，用前一个和后一个块的值的一半内核长度减一进行填充。

下面的代码演示了所描述的块填充产生的结果与根本不对信号进行分块完全相同。请注意，我将图中生成的卷积信号移动了

(len(kernel) - 1) / 2

值，以解释在信号开头使用

mode=valid

描述的缺失值。

import matplotlib.pyplot as plt
import numpy as np

len_kernel = 5
len_signal = 100
num_chunks = 10

rand = np.random.default_rng(seed=42)
signal = np.cumsum(rand.normal(size=len_signal))
kernel = np.ones(len_kernel) / len_kernel

signal_convolved = np.convolve(signal, kernel, mode="valid")

# Split into chunks, pad chunks: (1) at the start with `(len_kernel - 1) / 2`
# last elements from preceding chunk (except first chunk), (2) at the end with
# `(len_kernel - 1) / 2` first elements from following chunk (except last chunk)
len_arm = (len_kernel - 1) // 2  # Length of a "kernel arm"
chunks_convolved = []
for i, chunk in enumerate(chunks := np.split(signal, num_chunks)):
    start = [] if i == 0 else chunks[i - 1][-len_arm:]
    end = [] if i == len(chunks) - 1 else chunks[i + 1][:len_arm]
    chunk_conv = np.convolve(np.r_[start, chunk, end], kernel, mode="valid")
    chunks_convolved.append(chunk_conv)
chunks_convolved = np.r_[*chunks_convolved]
# Show results: just convolved (left) vs chunked, padded, and convolved (right)
assert np.allclose(signal_convolved, chunks_convolved)
x_sig, x_conv = np.arange(len_signal), np.arange(len(signal_convolved)) + len_arm
plt.subplot(121); plt.plot(x_sig, signal); plt.plot(x_conv, signal_convolved)
plt.subplot(122); plt.plot(x_sig, signal); plt.plot(x_conv, chunks_convolved)
plt.show()

最后附注：scipy 的

convolve1d()

提供了除零填充之外的更多选项来处理边界情况，因此可能是一种替代方案，具体取决于用例。