Вопрос или проблема
В своем исследовании я использую реализацию ядра ковариации Мatern в Numpy. Чтобы убедиться, что я понимаю, что происходит, я пытаюсь самостоятельно вывести формулы внутри реализации Numpy.
Например, рассмотрим анизотропное ядро Мatern с $\nu=2.5$ в двух измерениях:
$$K(x,y)=(1+d+\frac{d^2}{3})\exp(-d)$$
где
$$d=\sqrt{5}\sqrt{\sum_{i=1}^2 \left(\frac{x_i-y_i}{l_i}\right)^2}.$$
Теперь я хочу вычислить градиент по длинам масштабов для целей оптимизации гиперпараметров. Итак, я вычисляю:
$$
\partial_{l_i} K(x,y)=\frac{5}{3}\exp(-d)\left( 1 + d \right)\frac{(x_i-y_i)^2}{l_i^3}
$$
(вычислено с помощью следующего кода на Mathematica)
d = Sqrt[((x1 - y1)/l1)^2 + ((x2 - y2)/l2)^2]*Sqrt[5]
FullSimplify[
D[Exp[-d]*(1 + d + d^2/3),
l1]] /. {Sqrt[5] Sqrt[(x1 - y1)^2/l1^2 + (x2 - y2)^2/l2^2] -> tmp}
Пока все хорошо. Теперь я смотрю на реализацию Numpy (строки 1729-1752 в kernels.py):
# Нам нужно заново вычислить попарные расстояния по измерениям
если self.anisotropic:
D = (X[:, np.newaxis, :] - X[np.newaxis, :, :]) ** 2 / (
длина_масштаба**2
)
еще:
D = squareform(dists**2)[:, :, np.newaxis]
если self.nu == 0.5:
denominator = np.sqrt(D.sum(axis=2))[:, :, np.newaxis]
divide_result = np.zeros_like(D)
np.divide(
D,
denominator,
out=divide_result,
where=denominator != 0,
)
K_gradient = K[..., np.newaxis] * divide_result
elif self.nu == 1.5:
K_gradient = 3 * D * np.exp(-np.sqrt(3 * D.sum(-1)))[..., np.newaxis]
elif self.nu == 2.5:
tmp = np.sqrt(5 * D.sum(-1))[..., np.newaxis]
K_gradient = 5.0 / 3.0 * D * (tmp + 1) * np.exp(-tmp)
И мы наблюдаем множество тех же самых членов, мы заменяем $d$ на tmp. Однако, в то время как мои аналитические вычисления имеют член $\frac{(x_i-y_i)^2}{l_i^3}$, реализация numpy, похоже, имеет член $\frac{(x_i-y_i)^2}{l_i^2}$ (через определение D в верхней части приведенного кода).
Вряд ли реализация numpy ошибочная, может кто-нибудь объяснить, где я ошибся?
Большое спасибо заранее!
Ответ или решение
В вашем вопросе о реализации производной ковариационного ядра Матерна с помощью NumPy есть важные аспекты, которые стоит обсудить, чтобы прояснить возникшие неясности. Давайте детально рассмотрим отдельные компоненты как аналитического вычисления, так и реализации в NumPy.
Анализ производной ковариационного ядра Матерна
Вы рассматриваете ковариационное ядро Матерна для (\nu = 2.5) в двумерном пространстве, определяемое формулой:
[
K(x,y) = \left(1 + d + \frac{d^2}{3}\right) \exp(-d)
]
где
[
d = \sqrt{5} \sqrt{\sum_{i=1}^2 \left(\frac{x_i – y_i}{l_i}\right)^2}
]
Вы правильно отметили, что при вычислении производной по длинам масштабов (l_i) вы получили:
[
\partial_{l_i} K(x,y) = \frac{5}{3} \exp(-d) \left( 1 + d \right) \frac{(x_i – y_i)^2}{l_i^3}
]
Понимание реализации в NumPy
Обратимся к реализации, представленной в NumPy. Здесь необходимо обратить внимание на то, как вычисляется матрица расстояний (D):
D = (X[:, np.newaxis, :] - X[np.newaxis, :, :]) ** 2 / (length_scale**2)В данной строке кода расстояния вычисляются с учетом длины масштаба (l_i). Это означает, что 2D-нормировка расстояний уже учитывала масштаб, из-за чего в расчете (D) на самом деле происходит деление на (l_i^2).
Сопоставление результатов
В вашем аналитическом решении (d) учитывает (\left(\frac{x_i – y_i}{l_i}\right)^2), в то время как в NumPy вы делите на (l_i^2) при вычислении (D). Обратите внимание, что при правильной интерпретации обе реализации приводят к тому же результату, когда вы смотрите на общий вывод.
Разбор разности
Ваше вычисление приводит к выражению, которое содержит (\frac{(x_i-y_i)^2}{l_i^3}), в то время как в реализации NumPy (в частности, в строке) в первую очередь присутствует (\frac{(x_i – y_i)^2}{l_i^2}). Эта разность происходит из того, какจัดรูปแบบ (D) и связанным с ним вычислениям, что приводит к появлению разных степеней в выражении.
Заключение
Таким образом, можно заключить, что ваши вычисления по производной верны, однако для получения результата, эквивалентного реализации в NumPy, необходимо учесть, что изменения масштабов уже внедрены в выражение для (D). Поэтому, строго придерживаясь формальностей, у вас могло возникнуть впечатление о несоответствии. Однако обе реализации истинны в контексте их определения.
Если у вас остались дополнительные вопросы или требуется дальнейшая информация по данной теме, пожалуйста, дайте знать.