Я использую дополнение Шура (d
в моем коде ниже) для обновления обратной матрицы. Предположим, что у меня есть квадратная обратимая матрица A
размерности n-1
, а обратная ей - A.inverse
, затем я добавляю одно измерение в матрицу A
, получая new.A
, теперь с помощью дополнения Шура я получаю new.inverse
, используя A.inverse
и new.A
. Код представлен ниже:
new.inverse = function(A.inverse, new.A){
n = dim(A.inverse)[1]
sigma12 = new.A[1:n, n+1, drop=FALSE]
sigma21 = new.A[n+1, 1:n, drop=FALSE]
sigma22 = new.A[n+1, n+1, drop=FALSE]
d = (sigma22 - sigma21 %*% A.inverse %*% sigma12)^(-1)
a = A.inverse + A.inverse %*% sigma12 %*% d %*% sigma21 %*% A.inverse
b = - A.inverse %*% sigma12 %*% d
c = - d %*% sigma21 %*% A.inverse
return(rbind(cbind(a, b), c(c, d)))
}
Код Шура enet - это точная форма расчета (не приближение), например solve(new.A) = new.inverse(A.inverse, new.A)
, но когда я использую эту формулу много раз в al oop (Приблизительно 1000), new.inverse(A.inverse, new.A)
расходится с solve(new.A)
. Подсказка, как я могу решить эту проблему с числами?
Пример проблемы:
xx=5*runif(2)+0 #x coordinates
yy=5*runif(2)+0 #y coordinates
npoints = cbind(xx,yy)
distmatN = as.matrix(dist(npoints, method = "euclidean", upper=TRUE))
cov0 = exp(-distmatN)
inv.cov0 = solve(cov0)
res=vector()
for(i in 1:700){
xx_new = 5*runif(1)+0; yy_new = 5*runif(1)+0; new_point = cbind(xx_new,yy_new) #creating the new point
distmatNN = as.matrix(dist(rbind(npoints,new_point), method = "euclidean", upper=F)) #creating the new matrix distance with the new point
cov_new = exp(-distmatNN)
npoints = rbind(npoints,new_point, deparse.level = 0)
cov0 = cov_new
inv.cov0 = inverse.update.in(inv.cov0, cov0)
res[i] = max(abs(inv.cov0 - solve(cov0)))
}#;res
summary(res)
Значение res
должно быть очень маленьким, но оно увеличивается.