沿着每根轴应用DCT矩阵没有得到预期的结果。

Question

使用下面的代码，我希望得到的是 linalg.norm(y-z) 等于零。我遵循的理念是 在NumPy中沿三维阵列的各轴应用DFT矩阵？ 因为FFT和DCT都是可分离的单元线性变换。但不知为何，情况并非如此。

import numpy as np
from scipy.fftpack import dct 

x = np.random.rand(5,5)

D = dct(np.eye(5), norm='ortho')
y = np.dot(D,np.dot(x, D.T))
z = dct(dct(x, axis = 0 , norm = 'ortho'), axis = 1 , norm = 'ortho')

Answer 1

挑选 一个 以下解决方案中的一种：

1. 重新定义 D 通过指定 axis=0:

   D = dct(np.eye(n), axis=0, norm="ortho")
   

2. 重新定义 D 通过指定 .T 并使用默认的 axis=-1 (并不能真正推广到更高的维度...)。

   D = dct(np.eye(n), norm="ortho").T
   

3. 使用 D.T @ x 来表示沿DCT axis=0 的 x:

   y = D.T @ x @ D
   

在DFT的情况下，顺序并不重要，因为矩阵是对称的(D == D.T) 此外 变成单一的(D.conj().T @ D == 1). 但是DCT矩阵并不是对称的，所以你必须小心使用什么顺序。

考虑DCT-II的定义。

当你构造算子时 D 寓意 D @ x 取DCT-II沿行(axis=0)，那么 D 必须定义，以便正确处理变换的协方差和反方差。

---

完整的例子。

import numpy as np
from scipy.fftpack import dct 

n = 5
x = np.random.rand(n, n)
D = dct(np.eye(n), axis=0, norm="ortho")

y = D @ x @ D.T

z = x
z = dct(z, axis=0, norm="ortho")
z = dct(z, axis=1, norm="ortho")

>>> np.linalg.norm(z - y)
6.20161216470283e-16