使用MASS库 ginv() 计算任意矩阵的广义逆 (Moore-Penrose)

前面讲了矩阵的逆, 例如矩阵A的逆用 A ˉ 1 表示, 使用solve(A)可与计算A的逆.
并且满足 A %*%  ==  1  %*% A   ==   diag(min(ncol(A), nrow(A)))  ==   Identity matrix
如下 : 

并非所有的矩阵都有逆, 但是所有的矩阵都可有广义逆.
如果 广义逆满足下列条件,  n×m矩阵A+是矩阵A的Moore-Penrose逆:
用A+表示A的广义逆. 那么必须满足(后两个是对称性) : 
AA+A=A
A+AA+=A+
(AA+)T=AA+
(A+A)T=A+A


R中MASS包中的ginv()函数可以计算矩阵的Moore-Penrose逆. 例如:
> A=matrix(1:12,nrow=3,ncol=4)
> A
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12

使用solve只能计算正方矩阵的逆.
> solve(A)
Error in solve.default(A) : 'a' (3 x 4) must be square

需要用到MASS的ginv函数.
> library(MASS)
> ginv(A)
             [,1]        [,2]        [,3]
[1,] -0.483333333 -0.03333333  0.41666667
[2,] -0.244444444 -0.01111111  0.22222222
[3,] -0.005555556  0.01111111  0.02777778
[4,]  0.233333333  0.03333333 -0.16666667

验证 :
1. a %*% A %*% a =  a

> a %*% A %*% a
             [,1]        [,2]        [,3]
[1,] -0.483333333 -0.03333333  0.41666667
[2,] -0.244444444 -0.01111111  0.22222222
[3,] -0.005555556  0.01111111  0.02777778
[4,]  0.233333333  0.03333333 -0.16666667
> a
             [,1]        [,2]        [,3]
[1,] -0.483333333 -0.03333333  0.41666667
[2,] -0.244444444 -0.01111111  0.22222222
[3,] -0.005555556  0.01111111  0.02777778
[4,]  0.233333333  0.03333333 -0.16666667
精度问题, 不能直接划等号
> (a %*% A %*% a) == a
      [,1]  [,2]  [,3]
[1,] FALSE FALSE FALSE
[2,] FALSE FALSE FALSE
[3,] FALSE FALSE FALSE
[4,] FALSE FALSE FALSE


2.  A %*% a %*% A =  A
> A %*% a %*% A
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12
> A
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12
精度问题, 不能直接划等号
> (A %*% a %*% A)  ==  A
      [,1]  [,2]  [,3]  [,4]
[1,] FALSE FALSE FALSE FALSE
[2,] FALSE FALSE FALSE FALSE
[3,] FALSE FALSE FALSE FALSE


3. 
验证对称性
A %*% a ==  t(A %*% a)
a %*% A ==  t(a %*% A)

> A %*% a
           [,1]      [,2]       [,3]
[1,]  0.8333333 0.3333333 -0.1666667
[2,]  0.3333333 0.3333333  0.3333333
[3,] -0.1666667 0.3333333  0.8333333
> t(A %*% a)
           [,1]      [,2]       [,3]
[1,]  0.8333333 0.3333333 -0.1666667
[2,]  0.3333333 0.3333333  0.3333333
[3,] -0.1666667 0.3333333  0.8333333
精度问题, 不能直接划等号
> t(A %*% a) == A %*% a
      [,1]  [,2]  [,3]
[1,]  TRUE FALSE FALSE
[2,] FALSE  TRUE  TRUE
[3,] FALSE  TRUE  TRUE

> t(a %*% A)
     [,1] [,2] [,3] [,4]
[1,]  0.7  0.4  0.1 -0.2
[2,]  0.4  0.3  0.2  0.1
[3,]  0.1  0.2  0.3  0.4
[4,] -0.2  0.1  0.4  0.7
> a %*% A
     [,1] [,2] [,3] [,4]
[1,]  0.7  0.4  0.1 -0.2
[2,]  0.4  0.3  0.2  0.1
[3,]  0.1  0.2  0.3  0.4
[4,] -0.2  0.1  0.4  0.7
精度问题, 不能直接划等号
> t(a %*% A) == a %*% A
      [,1]  [,2]  [,3]  [,4]
[1,]  TRUE FALSE FALSE FALSE
[2,] FALSE  TRUE  TRUE FALSE
[3,] FALSE  TRUE  TRUE FALSE
[4,] FALSE FALSE FALSE  TRUE

[参考]
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值