提问人:JoeBass 提问时间:1/27/2016 最后编辑:CommunityJoeBass 更新时间:1/27/2016 访问量:3683
R chol 和正半定矩阵
R chol and positive semi-definite matrix
问:
我有以下矩阵:
j <- matrix(c(1,1,.5,1,1,.5,.5,.5,1), nrow=3, ncol=3)
这是正半定的,因为所有特征值都是 >= 0。
> eigen(j, symmetric = TRUE)
$values
[1] 2.3660254 0.6339746 0.0000000
$vectors
[,1] [,2] [,3]
[1,] -0.6279630 -0.3250576 7.071068e-01
[2,] -0.6279630 -0.3250576 -7.071068e-01
[3,] -0.4597008 0.8880738 -1.942890e-15
然而,cholesky 分解失败了......
> chol(j)
Error in chol.default(j) :
the leading minor of order 2 is not positive definite
我还从互联网上改编了一些代码......
cholesky_matrix <- function(A){
# http://rosettacode.org/wiki/Cholesky_decomposition#C
L <- matrix(0,nrow=nrow(A),ncol=ncol(A))
colnames(L) <- colnames(A)
rownames(L) <- rownames(A)
m <- ncol(L)
for(i in 1:m){
for(j in 1:i){
s <- 0
if(j > 1){
for(k in 1:(j-1)){
s <- s + L[i,k]*L[j,k]
}
}
if(i == j){
L[i,j] <- sqrt(A[i,i] - s)
} else {
L[i,j] <- (1 / L[j,j])*(A[i,j] - s)
}
}
}
return(L)
}
NaN 也“失败”。
> cholesky_matrix(j)
[,1] [,2] [,3]
[1,] 1.0 0 0
[2,] 1.0 0 0
[3,] 0.5 NaN NaN
有谁知道发生了什么?为什么我的分解失败?
答:
3赞
Hong Ooi
1/27/2016
#1
矩阵的特征值为
> eigen(j)
$values
[1] 2.366025e+00 6.339746e-01 4.440892e-16
最后一个实际上是零,在数值精度的范围内。每:?chol
计算实对称正定方阵的 Choleski 因式分解。
(强调我的)
也就是说,您仍然可以通过设置 来获得分解,它能够处理半确定性:pivot=TRUE
> chol(j, pivot=TRUE)
[,1] [,2] [,3]
[1,] 1 0.5000000 1
[2,] 0 0.8660254 0
[3,] 0 0.0000000 0
attr(,"pivot")
[1] 1 3 2
attr(,"rank")
[1] 2
Warning message:
In chol.default(j, pivot = TRUE) :
the matrix is either rank-deficient or indefinite
评论
0赞
JoeBass
1/27/2016
Cholesky 仍然适用于正半定正确吗?R 对警告消息如此具体有什么原因吗?
0赞
1/27/2017
我不知道,但这是一个巨大的拖累
1赞
Stéphane Laurent
11/25/2017
请注意,在这种情况下,平等不成立。看。t(C)%*%C=M
?chol
下一个:R !按因子水平复制
评论