R语言中矩阵运算.docx
- 文档编号:6619495
- 上传时间:2023-05-10
- 格式:DOCX
- 页数:17
- 大小:46.59KB
R语言中矩阵运算.docx
《R语言中矩阵运算.docx》由会员分享,可在线阅读,更多相关《R语言中矩阵运算.docx(17页珍藏版)》请在冰点文库上搜索。
R语言中矩阵运算
R语言中矩阵运算
目录:
矩阵的生成,矩阵的四则运算,矩阵的矩阵运算,矩阵的分解。
1.矩阵的生成
1_1将向量定义成数组
向量只有定义了维数向量(dim属性)后才能被看作是数组.比如:
>z=1:
12;
>dim(z)=c(3,4);
AA>z;
[,1][,2][,3][,4]
[1,]
1
4
7
10
[2,]
2
5
8
11
[3,]
3
6
9
12
注意:
生成矩阵是按列排列的。
1_2用array()函数构造多维数组
用法为:
array(data=NA,dim=length(data),dimnames=NULL)
参数描述:
data:
是一个向量数据。
dim:
是数组各维的长度,缺省时为原向量的长度。
dimname:
是数组维的名字,缺省时为空。
例子:
>x=array(1:
20,dim=c(4,5))
>x
[,1][,2][,3][,4][,5]
[1,]
1
5
9
13
17
[2,]
2
6
10
14
18
[3,]
3
7
11
15
19
[4,]
4
8
12
16
20
1_3用matrix()函数构造矩阵
函数matrix)是构造矩阵(二维数组)的函数,其构造形式为matrix(data=NA,nrow=1,ncol=1,byrow=FALSE,dimnames=NULL)
其中data是一个向量数据,nrow是矩阵的行数,ncol是矩阵的列数.当byrow=TRUE时,生成矩阵的数据按行放置,缺省时相当于byrow=t,数据按列放置.dimname。
是数组维的名字,缺省时为空.A
如构造一个3x5阶的矩阵
>A=matrix(1:
15,nrow=3,byrow=TRUE)
>A
[,1][,2][,3][,4][,5]
[1,]
1
2
3
4
5
[2,]
6
7
8
9
10
[3,]
11
12
13
14
15
2.矩阵的四则运算
可以对数组之间进行四则运算(+、一、*、/),这时进行的是数组对应元素的四则运算。
一般情况下参加运算的矩阵或者数组的维数是相同的,但也可以计算不同维的,这是要将对应的元素补足。
3.矩阵的矩阵运算
3_1运算
对于矩阵A,函数t(A)表示矩阵A的转置,如:
>A=matrix(1:
6,nrow=2);
>A;
[,1][,2][,3]
[1,]
1
3
5
[2,]
>t(A);
2
[,1][,2]
4
6
[1,]
1
2
[2,]
3
4
[3,]
5
6
3_2求方阵的行列式
函数det()是求矩阵行列式的值,如
>det(matrix(1:
4,ncol=2));
[1]-2
3_3向量的内积
对于n维向量x,可以看成nxl阶矩阵或lxn阶矩阵。
若x与y是相同维数的向量,则x%*%Y表示x与y作内积.例如,
>x=1:
5;Y=2*1:
5Z
>x%*%y
[,1]
[1,]110
函数crossprod()是内积运算函数(表示交叉乘积),crossprod(x,y)计算向量x与y的内积,即t(x)%*%y'。
crossprod(x)表示x与x的内积.
类似地,tcrossprod(x,y)表示’x%*%t(Y),’即x与y的外积,也称为叉积。
tcrossprod(x)
表示x与x作外积.如:
>x=1:
5;y=2*1:
5;
>crossprod(x);[,1]
[1,]55
>crossprod(x,y);
[,1]
[1,]110
>tcrossprod(x);
[,1][,2][,3][,4][,5]
[1,]
1
2
3
4
5
[2,]
2
4
6
8
10
[3,]
3
6
9
12
15
[4,]
4
8
12
16
20
[5,]
5
10
15
20
25
>tcrossprod(x,y);
[,1][,2][,3][,4][,5]
[1,]
2
4
6
8
10
[2,]
4
8
12
16
20
[3,]
6
12
18
24
30
[4,]
8
16
24
32
40
[5,]
10
20
30
40
50
3_4向量的外积(叉积)
设x和y是n维向量,则x%o%y表示x与y作外积.例如
>x%o%y;
[,1][,2][,3][,4][,5]
[1,]246810
[2,]
4
8
12
16
20
[3,]
6
12
18
24
30
[4,]
8
16
24
32
40
[5,]
10
20
30
40
50
outer()是更为强大的外积运算函数,outer(x,y)计算向量二与y的外积,它等价于
x%o%y
函数。
outer()的一般调用格式为
outer(x,y,fun=”*”)
其中x,y矩阵(或向量),fun是作外积运算函数,缺省值为乘法运算。
函数outer()在绘制三维曲面时非常有用,它可生成一个x和y的网格。
3_5矩阵的乘法
设A和B为两个矩阵,通常意义下的矩阵乘法是通过A%*%B来完成,crossprod(A,B)
表示的是
t(A)%*%B,而tcrossprod(A,B)表示的是A%*%t(B)。
最后我们通过运算知道x%*%A%*%x
为二次型。
例子:
>A=array(1:
9,dim=(c(3,3)))
>B=array(9:
1,dim=(c(3,3)))
>A%*%B;
[,1][,2][,3]
[1,]905418
[2,]
114
69
24
[3,]
138
84
30
>crossprod(A,B)==t(A)%*%B;[,1][,2][,3]
[1,]TRUETRUETRUE
[2,]TRUETRUETRUE[3,]TRUETRUETRUE
>tcrossprod(A,B)==A%*%t(B);[,1][,2][,3]
[1,]TRUETRUETRUE
[2,]TRUETRUETRUE[3,]TRUETRUETRUE
3_6生成对角阵和矩阵取对角运算
函数diag()依赖于它的变量,当v是一个向量时,diag(v)表示以v的元素为对角线元素的对角阵.当M是一个矩阵时,则diag(M)表示的是取M对角线上的元素的向量.如
>v=c(1,4,5);
>diag(v);
[,1][,2][,3]
[1,]
1
0
0
[2,]
0
4
0
[3,]
0
0
5
>M=array(1:
9,dim=c(3,3));
>diag(M);[1]159
3_7解线性方程组和求矩阵的逆矩阵
若求解线性方程组Ax=b,其命令形式为solve(A,b),求矩阵A的逆,其命令形式为solve(A).设矩阵A=t(array(c(1:
8,10),dim=c(3,3))),b<-c(1,1,1),则解方程组Ax=b的解x和求矩阵A的逆矩阵的命令如下:
>A=t(array(c(1:
8,10),dim=c(3,3)));
>b=c(1,1,1);
>x=solve(A,b);
>x;
[1]-1.000000e+001.000000e+003.806634e-16
>solve(A);
[,1][,2][,3]
[1,]-0.6666667-1.3333331
[2,]-0.66666673.666667-2
[3,]1.0000000-2.0000001
3_8求矩阵的特征值与特征向量
函数eigen(Sm)是求对称矩阵Sm的特征值与特征向量,其命令形式为:
ev=eigen(Sm),
则ev存放着对称矩阵Sm特征值和特征向量,是由列表形式给出的,其中ev$values是Sm
的特征值构成的向量,ev$vectors是Sm的特征向量构成的矩阵.如
>Sm=crossprod(A,A);
>ev=eigen(Sm);
>ev;
$values
[1]303.195336180.765907390.03875643
$vectors
[,1][,2][,3]
[1,]-0.46466750.8332863550.2995295
[2,]-0.5537546-0.009499485-0.8326258
[3,]-0.6909703-0.5527599940.4658502
4.矩阵的分解
4_1特征值分解
(1).定义:
对N阶方阵A,x为标量,v是非零的N维列向量,且满足Ax=xv,则称x为矩阵A
的特征值,v是相对应于x的特征向量。
特征值的全体成为A的谱。
(2).在r中的实现:
在r中利用函数eigen(A)来求矩阵的特征值和特征向量,具体的调用格式为:
以矩阵A为例说明此问题
>A=array(c(1,1,1,4,2,1,9,3,1),dim=c(3,3));
>D=eigen(A);
>D;
$values
[1]5.8284271-2.0000000
$vectors
0.1715729
[,1]
[,2]
[,3]
[1,]-0.8597736-9.486833e-01
0.5384820
[2,]-0.43464986.474883e-17-0.7872938
[3,]-0.26808393.162278e-010.3003425
(3).特征值分解的性质:
我们知道当所求的的特征向量构成的矩阵可逆时会满足
solve(vectors)%*%A%*%vectors=diag(values),下面进行验证。
>solve(vectors)%*%A%*%vectors;
[,1][,2][,3][1,]5.828427e+008.339683e-16-1.285213e-15[2,]1.211325e-15-2.000000e+002.704000e-16[3,]-3.471971e-16-1.607126e-161.715729e-01
结果的精度还是比较高的。
4_2矩阵的奇异值分解
函数svd(A)是对矩阵A作奇异值分解,即A=U%*%D%*%t(V),其中U,V是正交阵,D为对角阵,也就是矩阵A的奇异值.svd(A)的返回值也是列表,svd(A)$d表示矩阵A的奇异值,即矩阵D的对角线上的元素.svd(A)$u对应的是正交阵U,svd(A)$v对应的是正交阵
V.例如,
>A<-t(array(c(1:
8,10),dim=c(3,3)))
>SVD=svd(A);
>SVD;
$d
[1]17.4125052
0.8751614
0.1968665
$u
[,1]
[,2]
[,3]
[1,]-0.2093373
0.96438514
0.1616762
[2,]-0.50384850.03532145-0.8630696
[3,]-0.8380421-0.26213299
0.4785099
$v
[,1]
[,2]
[,3]
[1,]-0.4646675-0.8332863550.2995295
[2,]-0.55375460.009499485-0.8326258
[3,]-0.69097030.5527599940.4658502
>attach(SVD);
Thefollowingobject(s)aremaskedfrom'SVD(position3)':
d,u,v
>u%*%diag(d)%*%t(v);
[,1][,2][,3]
[1,]123
[2,]
4
5
6
[3,]
>A;
7
8
10
[,1][,2][,3]
[1,]123
[2,]456
[3,]7810
4_3qr分解
设A为m*n矩阵,如果存在m*m酉矩阵Q(即Q(H)Q=QQ(H)=I)和m*n阶梯形矩阵R,使得A=QR,那么此分解称为QR分解。
QR分解在解决最小二乘问题、特征值计算等方面有着十分重要的作用。
#建立矩阵
>A=(array(c(1:
12),dim=c(4,3)));
>A;
[,1][,2][,3]
[1,]
1
5
9
[2,]
2
6
10
[3,]
3
7
11
[4,]
4
8
12
#进行矩阵分解
>QR=qr(A);QR
$qr
[,1][,2][,3]
[1,]-5.4772256-12.7801930-2.008316e+01
[2,]0.3651484-3.2659863-6.531973e+00
[3,]
0.5477226
-0.3781696
7.880925e-16
[4,]
$rank
0.7302967
-0.9124744
9.277920e-01
[1]2
$qraux
[1]1.1825741.1561351.373098
$pivot[1]123
attr(,"class")
[1]"qr"
#提取Q,R并验证分解的正确性。
>Q=qr.Q(QR);
>R=qr.R(QR);
>Q%*%R;
[,1][,2][,3]
[1,]
1
5
9
[2,]
2
6
10
[3,]
3
7
11
[4,]4812
4_4Schur分解
引言:
从特征值的分解中可以看出,特征值的分解是有条件的,如果特征向量不是线性无关的,那么对于一个矩阵来说便不能采用特征值分解的方法对矩阵进行分解。
例如对于矩阵A=t(array(c(6,12,19,-9,-20,-33,4,9,15),dim=c(3,3))进行特征值分解有:
>A=t(array(c(6,12,19,-9,-20,-33,4,9,15),dim=c(3,3)));
>A;
[,1][,2][,3]
[1,]61219
[2,]-9-20-33
[3,]4915
>det(A);
[1]-1
>W=eigen(A);
>W;
$values
[1]1+0i1-0i-1+0i
$vectors
[,1][,2][,3][1,]-0.4082483-0i-0.4082483+0i-0.4740998+0i
[2,]0.8164966+0i0.8164966+0i0.8127426+0i[3,]-0.4082483+0i-0.4082483-0i-0.3386427+0i
>attach(W);
Thefollowingobject(s)aremaskedfrom'W(position3)':
values,vectors
>det(vectors);
错误于determinant.matrix(x,logarithm=TRUE,...):
目前还不能算复数矩阵的行列式
>det(Re(vectors));[1]-7.599489e-19
>solve(vectors)
[,1][,2][,3][1,]0.000000+78209959i0.00000+78209959i-9.26965+78209959i[2,]0.000000-78209959i0.00000-78209959i-9.10153-78209959i[3,]3.691206+0i11.07362+0i18.45603+0i
很明显vectors不是一个可逆矩阵此时进行特征值分辨这种方法便不可行,对于这种情况我们可以作Schur分解。
描述:
对于任意的方针A,其Schur分解的形式为:
A=USU(H),其中U是标准的正交矩阵(即满足UU(H)=I),S为上三角矩阵,并且对角线上的元素为A的特征值。
由于此函数在包Matrix中,所以使用之前必须调入。
并且注意matrix和Matrix的区别。
例子:
>A=Matrix(c(6,12,19,-9,-20,-33,4,9,15),ncol=3,byrow=TRUE);
>A;
3x3Matrixofclass"dgeMatrix"[,1][,2][,3]
[1,]61219
[2,]-9-20-33
[3,]4915
>library(Matrix);
>Sch=Schur(A,vectors=TRUE);
>Q=Sch@Q;
>Q=as.matrix(Q)
>attach(Sch);
错误于attach(Sch):
'attach'只适用于串列,数据框和环境
>Q%*%T%*%t(Q)
3x3Matrixofclass"dgeMatrix"[,1][,2][,3]
[1,]61219
[2,]-9-20-33
[3,]4915
4_5Cholesky分解(柯利分解)
描述:
正定矩阵:
设A是n阶实系数矩阵,如果对任何非零向量X=(x1,...xn)都有t(X)AX>0,就称A正定(PositiveDefinite)。
正定矩阵在相合变换下可化为标准型,即单位矩阵。
Cholesky分解:
对任意的正定矩阵A,存在上三角矩阵R,使A=t(R)%*%R,则称为A的Cholesky分解(柯利分解)。
例子:
>#输入矩阵
>m=matrix(c(5,1,1,3),ncol=2);
>m;
[,1][,2]
[1,]
5
1
[2,]
1
3
>#矩阵分解
>CH=chol(m);
>#验证结果
>t(CH)%*%CH;
[,1][,2]
[1,]
5
1
[2,]
1
3
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 语言 矩阵 运算