如何对两个数组进行排序以获得平滑度
Wol*_*ger 6 python arrays sorting algorithm
说我有同样大小的两个数组A和B.为了明确,我们假设它们是二维的.两个数组中的值表示一些数据,这些数据应平滑地取决于数组中的位置.但是,其中的一些值A已与其对应的值交换B.我想扭转这种交换,但我很难找到一个标准来告诉我何时应该交换两个值.
一个例子应该有帮助.下面是一些python代码,它创建了两个这样的数组,并随机交换了一些元素
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
import random
### create meshgrid ###
x = np.linspace(-10,10,15);
y = np.linspace(-10,10,11);
[X,Y] = np.meshgrid(x,y);
### two sufficiently smooth functions on the meshgrid ###
A = -X**2-Y**2;
B = X**2-Y**2-100;
### plot ###
ax=plt.subplot(2,2,1)
im1=ax.imshow(A,extent=[-10, 10, -10, 10])
ax.set_title('A')
ax2=plt.subplot(2,2,2)
im2=ax2.imshow(B,extent=[-10, 10, -10, 10])
ax2.set_title('B')
### randomly exchange a few of the elements of A and B ###
for i in np.arange(0,15):
for j in np.arange(0,11):
randNumb = random.random();
if randNumb>0.8:
mem=A[j,i];
A[j,i] = B[j,i];
B[j,i] = mem;
### plot for comparison ###
ax=plt.subplot(2,2,3)
im1=ax.imshow(A,extent=[-10, 10, -10, 10])
ax2=plt.subplot(2,2,4)
im2=ax2.imshow(B,extent=[-10, 10, -10, 10])
plt.show()
这导致以下情节:
上面两个图是原始阵列A和B,下面两个是混洗版本.现在,练习是为了逆转这个过程,即从洗牌版本创建原始阵列A和B.
关于'平滑'的意思.当然,算法仅在原始数据实际上足够平滑时才起作用,这意味着一个阵列中的相邻数据点对于所有点的值都足够接近.解决方案可以假设是这种情况.
另请注意,在上面的示例中,此练习非常容易通过眼睛进行.我在实现这个问题时遇到的问题是找到一个很好的标准来告诉我何时交换到单元格.
注意,允许反向变换重新标记A和B.
一种稳健的方法是将 4 个相邻像素值的平均值与每个像素中实际存在的值进行比较。A也就是说，对于每个像素，计算和中 4 个相邻像素的平均值，并将它们与和B中该像素的实际值进行比较。以下条件效果很好，并且实际上是一种最小二乘法： AB
if ( (A[i, j] - A_mean)**2 + (B[i, j] - B_mean)**2
> (A[i, j] - B_mean)**2 + (B[i, j] - A_mean)**2
):
# Do swap
这里，A_mean和B_mean是 4 个相邻像素的平均值。
另一件需要考虑的重要事情是，一次扫描所有像素不一定足够：可能会发生这样的情况：在一次扫描之后，校正交换使得上述条件能够识别更多应该交换的像素。因此，我们必须扫描数组，直到找到“稳定状态”。
完整代码
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
import random
### create meshgrid ###
x = np.linspace(-10,10,15);
y = np.linspace(-10,10,11);
[X,Y] = np.meshgrid(x,y);
### two sufficiently smooth functions on the meshgrid ###
A = -X**2-Y**2;
B = X**2-Y**2-100;
### plot ###
ax=plt.subplot(3,2,1)
im1=ax.imshow(A,extent=[-10, 10, -10, 10])
ax.set_title('A')
ax2=plt.subplot(3,2,2)
im2=ax2.imshow(B,extent=[-10, 10, -10, 10])
ax2.set_title('B')
### randomly exchange a few of the elements of A and B ###
for i in np.arange(0,15):
for j in np.arange(0,11):
randNumb = random.random();
if randNumb>0.8:
mem=A[j,i];
A[j,i] = B[j,i];
B[j,i] = mem;
### plot for comparison ###
ax=plt.subplot(3,2,3)
im1=ax.imshow(A,extent=[-10, 10, -10, 10])
ax2=plt.subplot(3,2,4)
im2=ax2.imshow(B,extent=[-10, 10, -10, 10])
### swap back ###
N, M = A.shape
swapped = True
while swapped:
swapped = False
for i in range(N):
for j in range(M):
A_mean = np.mean([A[i - 1 , j - 1 ],
A[i - 1 , (j + 1)%M],
A[(i + 1)%N, j - 1 ],
A[(i + 1)%N, (j + 1)%M],
])
B_mean = np.mean([B[i - 1 , j - 1 ],
B[i - 1 , (j + 1)%M],
B[(i + 1)%N, j - 1 ],
B[(i + 1)%N, (j + 1)%M],
])
if ( (A[i, j] - A_mean)**2 + (B[i, j] - B_mean)**2
> (A[i, j] - B_mean)**2 + (B[i, j] - A_mean)**2
):
# Do swap
A[i, j], B[i, j] = B[i, j], A[i, j]
swapped = True
### plot result ###
ax=plt.subplot(3,2,5)
im1=ax.imshow(A,extent=[-10, 10, -10, 10])
ax2=plt.subplot(3,2,6)
im2=ax2.imshow(B,extent=[-10, 10, -10, 10])
plt.show()
请注意，上面的代码认为数组是周期性的，因为边界处的像素的相邻像素被选择为相对边界上的像素（您在示例中提供的数组就是这种情况）。当索引变为负数时，这种索引包装会自动发生，但当索引大于或等于数组的给定维度时，不会自动发生，这就是%使用模运算符的原因。
作为额外的技巧，请注意我如何交换A[i, j]并且B[i, j]不需要临时mem变量。另外，我的外循环位于第一维度（长度为 11 的维度），而内循环位于第二维度（长度为 15 的维度）。这比执行反向循环顺序更快，因为现在每个元素都是按连续顺序访问的（它们在内存中实际存在的顺序）。
最后请注意，此方法不能保证始终有效。可能会发生交换太多附近像素而无法找到“正确”解决方案的情况。然而，无论您选择什么标准来确定两个像素是否应该交换，情况都会如此。
编辑（非周期数组）
对于非周期性阵列，边界像素的邻居数少于 4 个（边缘像素 3 个，角像素 2 个）。沿着这些思路：
A_neighbors = []
B_neighbors = []
if i > 0 and j > 0:
A_neighbors.append(A[i - 1, j - 1])
B_neighbors.append(B[i - 1, j - 1])
if i > 0 and j < M - 1:
A_neighbors.append(A[i - 1, j + 1])
B_neighbors.append(B[i - 1, j + 1])
if i < N - 1 and j > 0:
A_neighbors.append(A[i + 1, j - 1])
B_neighbors.append(B[i + 1, j - 1])
if i < N - 1 and j < M - 1:
A_neighbors.append(A[i + 1, j + 1])
B_neighbors.append(B[i + 1, j + 1])
A_mean = np.mean(A_neighbors)
B_mean = np.mean(B_neighbors)
请注意，邻居越少，该方法的鲁棒性就越差。您还可以尝试使用最近的 8 个像素（即，包括对角邻居），而不仅仅是 4 个。