разделить двумерное гауссово ядро ​​на два одномерных ядра - PullRequest
Новая
       17

разделить двумерное гауссово ядро ​​на два одномерных ядра

0 голосов
/ 22 сентября 2019

Гауссовское ядро ​​вычисляется и проверяется, что его можно отделить, посмотрев на ранг ядра.

kernel = gaussian_kernel(kernel_size,sigma)
print(kernel)

[[ 0.01054991  0.02267864  0.0292689   0.02267864  0.01054991]
 [ 0.02267864  0.04875119  0.06291796  0.04875119  0.02267864]
 [ 0.0292689   0.06291796  0.0812015   0.06291796  0.0292689 ]
 [ 0.02267864  0.04875119  0.06291796  0.04875119  0.02267864]
 [ 0.01054991  0.02267864  0.0292689   0.02267864  0.01054991]]

 rank = np.linalg.matrix_rank(kernel)

if rank == 1:
    print('The Kernel is separable')
else:
    print('The kernel is not separable')

Теперь я считаю, что разделение не является правильным.Я делаю это следующим образом: 

 u,s,v = np.linalg.svd(kernel)
 k1 = (u[:,0] * np.sqrt(s[0]))[np.newaxis].T
 k2 = v[:,0] * np.sqrt(s[0])

Затем я умножил два вышеупомянутых ядра, чтобы вернуть исходное ядро.Но я не получил его. 

if not  np.all(k1 * k2 == kernel):
    print('k1 * k2 is not equal to kernel')

Я предполагаю, что разделение, которое я пытаюсь сделать с помощью SVD и далее, не является правильным.Некоторое объяснение поможет.

1 Ответ

0 голосов
/ 22 сентября 2019

матричный ранг 1 означает, что все строки равны нулю или одинаковы с точностью до масштабирования, и то же самое верно для столбцов.Они также имеют масштабирование, равное двум факторам.Поэтому вы можете восстановить их, используя что-то вроде 

I,J = np.unravel_index(np.abs(kernel).argmax(), kernel.shape)
f1 = np.nansum(kernel / (kernel[None,:,J]@kernel),1,keepdims=True)
f2 = np.nansum(kernel / (kernel@kernel[I,:,None]),0,keepdims=True)
scaling = np.sqrt(np.abs(kernel).sum()/np.abs(f1*f2).sum())
f1 *= scaling * np.sign(f1[I,0]) * np.sign(kernel[I,J])
f2 *= scaling * np.sign(f2[0,J])

Обратите внимание, что большая часть сложности связана с моей попыткой собрать как можно больше данных.Проще, но я бы предположил, что численно не совсем стабильный метод будет 

I,J = np.unravel_index(np.abs(kernel).argmax(), kernel.shape)
f1 = kernel[:,J,None]
f2 = kernel[None,I,:] / kernel[I,J]

Конечно, ваш метод также работает, когда вы правильно настроили индексирование: 

k1 = u[:,0,None] * np.sqrt(s[0])
k2 = v[None,0,:] * np.sqrt(s[0])

np.allclose(kernel, k1*k2)
# True
...