如何使用Theano将“任意”操作放入滑动窗口？

Question

我想在矩阵X上定义一些函数。例如，mean(pow(X - X0, 2))，其中X0是另一个矩阵(X0是固定/常数)。为了更具体地说明，让我们假设X和X0都是10 x 10矩阵。操作的结果是一个实数。

现在我有了一个大矩阵(比方说500 x 500)。我想将上面定义的操作应用于“大”矩阵的所有10 x 10子矩阵。换句话说，我想将10 x 10窗口滑过“大”矩阵。每个窗口的位置，我应该得到一个真实的数字。因此，作为最后的结果，我需要得到一个实值矩阵(或2D张量)(它的形状应该是491 x 491)。

我想要的是接近卷积层，但不是完全相同的，因为我想使用均方偏差，而不是用神经元表示的线性函数。

Answer 1

这只是一个解决办法，希望它就足够了。我假设你的函数是由矩阵元素上的一个运算和一个平均值，即一个缩放和组成的。因此，只需将Y看作

Y = np.power(X-X0, 2)

所以我们只需要确定一个加窗的平均值。请注意，对于一维情形，可以确定具有适当向量的矩阵乘积来计算平均值。

h = np.array([0, 1, 1, 0])  # same dimension as y
m1 = np.dot(h, y) / 2 
m2 = (y[1] + y[2]) / 2
print(m1 == m2)  # True

2D情况类似，但有两个矩阵乘法，一个用于行，另一个用于列。例如。

m_2 = np.dot(np.dot(h, Y), h) / 2**2

要构造一个滑动窗口，我们需要建立一个移动窗口的矩阵。

H = [[1, 1, 1, 0, 0, ..., 0],
     [0, 1, 1, 1, 0, ..., 0],
            .
            .
            .
     [0, ..., 0, 0, 1, 1, 1]] 

计算所有的总和

S = np.dot(np.dot(H, Y), H.T)

具有(n, n)窗口的(m, m)矩阵的完整示例如下

import numpy as np

n, m = 500, 10
X0 = np.ones((n, n))
X = np.random.rand(n, n)
Y = np.power(X-X0, 2)

h = np.concatenate((np.ones(m), np.zeros(n-m)))  # window at position 0
H = np.vstack((np.roll(h, k) for k in range(n+1-m)))  # slide the window 
M = np.dot(np.dot(H,Y), H.T) / m**2  # calculate the mean
print(M.shape)  # (491, 491)

构建H的另一种但可能稍微不那么有效的方法是

H = np.sum(np.diag(np.ones(n-k), k)[:-m+1, :] for k in range(m))

更新

用这种方法也可以计算平均平方偏差。为此，我们将向量恒等式|x-x0|^2 = (x-x0).T (x-x0) = x.T x - 2 x0.T x + x0.T x0 (空间表示标量或矩阵乘法，.T表示转置向量)推广到矩阵情形：

我们假设W是一个包含块(m.m)恒等矩阵的(m,n)矩阵，它能够通过Y = W Z W.T提取(k0,k1)-th (m,m)子矩阵，其中Z是包含数据的(n,n)矩阵。计算差

 D = Y - X0 = Y = W Z W.T - X0 

非常简单，其中X0和D是(m,m)矩阵。元素平方和的平方根称为Frobenius范数。基于这些身份，我们可以将平方和写成

s = sum_{i,j} D_{i,j}^2 = trace(D.T D) = trace((W Z W.T - X0).T (H Z H.T - X0))
  = trace(W Z.T W.T W Z W.T) - 2 trace(X0.T W Z W.T) + trace(X0.T X0)
  =: Y0 + Y1 + Y2

Y0一词可以解释为H Z H.T，而above.The项Y1的方法可以解释为Z上的加权均值，Y2是一个常数，只需确定一次。因此，一项可能的执行是：

import numpy as np

n, m = 500, 10
x0 = np.ones(m)
Z = np.random.rand(n, n)

Y0 = Z**2
h0 = np.concatenate((np.ones(m), np.zeros(n-m)))
H0 = np.vstack((np.roll(h0, k) for k in range(n+1-m)))
M0 = np.dot(np.dot(H0, Y0), H0.T)

h1 = np.concatenate((-2*x0, np.zeros(n-m)))
H1 = np.vstack((np.roll(h1, k) for k in range(n+1-m)))
M1 = np.dot(np.dot(H1, Z), H0.T)

Y2 = np.dot(x0, x0)
M = (M0 + M1) / m**2 + Y2