Я пишу функцию, чтобы взять расстояние Махаланобиса между двумя векторами. Я понимаю, что это достигается с помощью уравнения a ‘* C ^ -1 * b, где a и b – векторы, а C – ковариационная matrix. Мой вопрос в том, есть ли эффективный способ найти обратную матрицу без использования исключения Гаусса-Йордана, или нет никакого способа обойти это? Я ищу способ сделать это сам, а не с любыми предопределенными функциями.
Я знаю, что C – эрмитова положительно определенная matrix, так можно ли каким-либо образом алгоритмически воспользоваться этим фактом? Или есть какой-то умный способ вычислить расстояние Махаланобиса, не вычисляя обратную ковариацию? Любая помощь будет оценена по достоинству.
*** Изменить: неверное уравнение расстояния Махаланобиса указано неверно. Это должно быть x ‘* C ^ -1 * x, где x = (ba), а b и a – два вектора, расстояние которых мы пытаемся найти (спасибо LRPurser). Таким образом, решение, выбранное в выбранном ответе, выглядит следующим образом:
d = x ‘* b, где b = C ^ -1 * x C * b = x, поэтому разрешите для b с использованием факторизации LU или факторизации LDL.
Вы можете (и должны!) Использовать декомпозицию LU вместо явной инвертирования матрицы: решение C x = b
с использованием разложения имеет лучшие числовые свойства, чем вычисление C^-1
и умножение вектора b
.
Поскольку ваша matrix симметрична, декомпозиция LU эффективно эквивалентна декомпозиции LDL * , которую вы действительно должны использовать в вашем случае. Поскольку ваша matrix также положительно определена, вы должны иметь возможность выполнять это разложение без поворота.
Изменить: обратите внимание, что для этого приложения вам не нужно решать всю проблему C x = b
.
Вместо этого, учитывая C = LDL*
и разностный вектор v = ab
, разрешите L* y = v
для y
(что в два раза больше, чем полный решатель LU).
Тогда y^t D^-1 y = v^t C^-1 v
можно вычислить в линейном времени.
Первое расстояние Махаланобиса (MD) является нормированным расстоянием относительно неопределенности в измерении двух векторов. Когда C=Indentity matrix
, MD сводится к евклидову расстоянию, и, следовательно, произведение сводится к векторной норме. Также MD всегда положительно определен или больше нуля для всех ненулевых векторов. По вашей формулировке с соответствующим выбором векторов a
и b
a*C^-1*b
может быть меньше нуля. Надеемся, что разность векторов, которые вы ищете, равна x=(ab)
что делает вычисление x^t*C^-1*x
где x^t
– транспонирование вектора x
. Также обратите внимание, что MD=sqrt(x^t*C^-1*x)
Так как ваша matrix симметрична и положительно определена, вы можете использовать разложение Холески (MatLab-chol)
которое использует половину операций как LU и численно больше стабильный. chol(C)=L
где C=L*L^t
где L
– нижняя треугольная matrix, а L^t
– транспонирование L
которая делает ее верхней треугольной. Ваш алгоритм должен идти примерно так:
(Matlab)
x=ab; L=chol(C); z=L\x; MD=z'*z; MD=sqrt(MD);