R语言中与矩阵相关的所有操作(上)

您所在的位置:网站首页 r语言生成矩阵(1,2,3,4;2,3,4,1) R语言中与矩阵相关的所有操作(上)

R语言中与矩阵相关的所有操作(上)

2024-07-17 06:28| 来源: 网络整理| 查看: 265

主要包括以下内容:创建矩阵向量;矩阵加减,乘积;矩阵的逆;行列式的值;特征值与特征向量;QR分解;奇异值分解;广义逆;backsolve与fowardsolve函数;取矩阵的上下三角元素;向量化算子等.

1   创建一个向量在R中可以用函数c()来创建一个向量,例如:> x=c(1,2,3,4)> x[1] 1 2 3 4 2   创建一个矩阵在R中可以用函数matrix()来创建一个矩阵,应用该函数时需要输入必要的参数值。> args(matrix)function (data = NA, nrow = 1, ncol = 1, byrow = FALSE, dimnames = NULL) data项为必要的矩阵元素,nrow为行数,ncol为列数,注意nrow与ncol的乘积应为矩阵元素个数,byrow项控制排列元素时是否按行进行,dimnames给定行和列的名称。例如:> matrix(1:12,nrow=3,ncol=4)    [,1] [,2] [,3] [,4][1,]   1   4   7   10[2,]   2   5   8   11[3,]   3   6   9   12> rowname[1] "r1" "r2" "r3"> colname=c("c1","c2","c3","c4")> colname[1] "c1" "c2" "c3" "c4"> matrix(1:12,nrow=3,ncol=4,dimnames=list(rowname,colname))  c1 c2 c3 c4r1 1 4 7 10r2 2 5 8 113   矩阵转置A为m×n矩阵,求A'在R中可用函数t(),例如:> 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> t(A)   [,1] [,2] [,3][1,]   1   2   3[2,]   4   5   6[3,]   7   8   9[4,]   10   11   12若将函数t()作用于一个向量x,则R默认x为列向量,返回结果为一个行向量,例如:> x[1] 1 2 3 4 5 6 7 8 9 10> t(x)  [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10][1,]   1   2   3   4   5   6   7   8   9   10> class(x)[1] "integer"> class(t(x))[1] "matrix"若想得到一个列向量,可用t(t(x)),例如:> x[1] 1 2 3 4 5 6 7 8 9 10> t(t(x))    [,1][1,]   1[2,]   2[3,]   3[4,]   4[5,]   5[6,]   6[7,]   7[8,]   8[9,]   9[10,]  104   矩阵相加减在R中对同行同列矩阵相加减,可用符号:“+”、“-”,例如:> A=B=matrix(1:12,nrow=3,ncol=4)> A+B    [,1] [,2] [,3] [,4][1,]   2   8   14   20[2,]   4   10   16   22[3,]   6   12   18   24> A-B   [,1] [,2] [,3] [,4][1,]   0   0   0   0[2,]   0   0   0   0[3,]   0   0   0   05   数与矩阵相乘A为m×n矩阵,c>0,在R中求cA可用符号:“*”,例如:> c=2> c*A    [,1] [,2] [,3] [,4][1,]   2   8   14   20[2,]   4   10  16   22[3,]   6   12  18   246   矩阵相乘A为m×n矩阵,B为n×k矩阵,在R中求AB可用符号:“%*%”,例如:> A=matrix(1:12,nrow=3,ncol=4)> B=matrix(1:12,nrow=4,ncol=3)> A%*%B    [,1] [,2] [,3][1,]   70  158 246[2,]   80  184 288[3,]   90  210 330若A为n×m矩阵,要得到A'B,可用函数crossprod(),该函数计算结果与t(A)%*%B相同,但是效率更高。例如:> A=matrix(1:12,nrow=4,ncol=3)> B=matrix(1:12,nrow=4,ncol=3)> t(A)%*%B    [,1] [,2] [,3][1,]  30   70 110[2,]  70  174 278[3,] 110  278 446> crossprod(A,B)    [,1] [,2] [,3][1,]  30  70 110[2,]  70 174 278[3,] 110 278 446矩阵Hadamard积:若A={aij}m×n, B={bij}m×n, 则矩阵的Hadamard积定义为:A⊙B={aij bij }m×n,R中Hadamard积可以直接运用运算符“*”例如:> A=matrix(1:16,4,4)> A    [,1] [,2] [,3] [,4][1,]   1   5   9   13[2,]   2   6   10   14[3,]   3   7   11   15[4,]   4   8   12   16> B=A> A*B    [,1] [,2] [,3] [,4][1,]   1   25   81 169[2,]   4   36 100 196[3,]   9   49 121 225[4,]   16   64 144 256R中这两个运算符的区别区加以注意。7   矩阵对角元素相关运算例如要取一个方阵的对角元素,> A=matrix(1:16,nrow=4,ncol=4)> A    [,1] [,2] [,3] [,4][1,]   1   5   9   13[2,]   2   6   10   14[3,]   3   7   11   15[4,]   4   8   12   16> diag(A)[1] 1 6 11 16对一个向量应用diag()函数将产生以这个向量为对角元素的对角矩阵,例如:> diag(diag(A))    [,1] [,2] [,3] [,4][1,]   1   0   0   0[2,]   0   6   0   0[3,]   0   0   11   0[4,]   0   0   0   16对一个正整数z应用diag()函数将产生以z维单位矩阵,例如:> diag(3)    [,1] [,2] [,3][1,]   1   0   0[2,]   0   1   0[3,]   0   0   18   矩阵求逆矩阵求逆可用函数solve(),应用solve(a, b)运算结果是解线性方程组ax = b,若b缺省,则系统默认为单位矩阵,因此可用其进行矩阵求逆,例如:> a=matrix(rnorm(16),4,4)> a            [,1]     [,2]     [,3]     [,4][1,] 1.6986019   0.5239738 0.2332094 0.3174184[2,] -0.2010667 1.0913013 -1.2093734   0.8096514[3,] -0.1797628 -0.7573283 0.2864535 1.3679963[4,] -0.2217916 -0.3754700 0.1696771 -1.2424030> solve(a)              [,1]     [,2]     [,3]     [,4][1,] 0.9096360 0.54057479 0.7234861 1.3813059[2,] -0.6464172 -0.91849017 -1.7546836 -2.6957775[3,] -0.7841661 -1.78780083 -1.5795262 -3.1046207[4,] -0.0741260 -0.06308603 0.1854137 -0.66078519   矩阵的特征值与特征向量矩阵A的谱分解为A=UΛU',其中Λ是由A的特征值组成的对角矩阵,U的列为A的特征值对应的特征向量,在R中可以用函数eigen()函数得到U和Λ,> args(eigen)function (x, symmetric, only.values = FALSE, EISPACK = FALSE)其中:x为矩阵,symmetric项指定矩阵x是否为对称矩阵,若不指定,系统将自动检测x是否为对称矩阵。例如:> A=diag(4)+1> A  [,1] [,2] [,3] [,4][1,]   2   1   1   1[2,]   1   2   1   1[3,]   1   1   2   1[4,]   1   1   1   2> A.eigen=eigen(A,symmetric=T)> A.eigen$values[1] 5 1 1 1$vectors        [,1]     [,2]       [,3]     [,4][1,] 0.5 0.8660254 0.000000e+00 0.0000000[2,] 0.5 -0.2886751 -6.408849e-17 0.8164966[3,] 0.5 -0.2886751 -7.071068e-01 -0.4082483[4,] 0.5 -0.2886751 7.071068e-01 -0.4082483> A.eigen$vectors%*%diag(A.eigen$values)%*%t(A.eigen$vectors)  [,1] [,2] [,3] [,4][1,]   2   1   1   1[2,]   1   2   1   1[3,]   1   1   2   1[4,]   1   1   1   2> t(A.eigen$vectors)%*%A.eigen$vectors            [,1]       [,2]         [,3]         [,4][1,] 1.000000e+00 4.377466e-17 1.626303e-17 -5.095750e-18[2,] 4.377466e-17 1.000000e+00 -1.694066e-18 6.349359e-18[3,] 1.626303e-17 -1.694066e-18 1.000000e+00 -1.088268e-16[4,] -5.095750e-18 6.349359e-18 -1.088268e-16 1.000000e+0010   矩阵的Choleskey分解  对于正定矩阵A,可对其进行Choleskey分解,即:A=P'P,其中P为上三角矩阵,在R中可以用函数chol()进行Choleskey分解,例如:> A  [,1] [,2] [,3] [,4][1,]   2   1   1   1[2,]   1   2   1   1[3,]   1   1   2   1[4,]   1   1   1   2> chol(A)        [,1]     [,2]     [,3]     [,4][1,] 1.414214 0.7071068 0.7071068 0.7071068[2,] 0.000000 1.2247449 0.4082483 0.4082483[3,] 0.000000 0.0000000 1.1547005 0.2886751[4,] 0.000000 0.0000000 0.0000000 1.1180340> t(chol(A))%*%chol(A)  [,1] [,2] [,3] [,4][1,]   2   1   1   1[2,]   1   2   1   1[3,]   1   1   2   1[4,]   1   1   1   2> crossprod(chol(A),chol(A))  [,1] [,2] [,3] [,4][1,]   2   1   1   1[2,]   1   2   1   1[3,]   1   1   2   1[4,]   1   1   1   2若矩阵为对称正定矩阵,可以利用Choleskey分解求行列式的值,如:> prod(diag(chol(A))^2)[1] 5> det(A)[1] 5若矩阵为对称正定矩阵,可以利用Choleskey分解求矩阵的逆,这时用函数chol2inv(),这种用法更有效。如:> chol2inv(chol(A))      [,1] [,2] [,3] [,4][1,] 0.8 -0.2 -0.2 -0.2[2,] -0.2 0.8 -0.2 -0.2[3,] -0.2 -0.2 0.8 -0.2[4,] -0.2 -0.2 -0.2 0.8> solve(A)  [,1] [,2] [,3] [,4][1,] 0.8 -0.2 -0.2 -0.2[2,] -0.2 0.8 -0.2 -0.2[3,] -0.2 -0.2 0.8 -0.2[4,] -0.2 -0.2 -0.2 0.811   矩阵奇异值分解  A为m×n矩阵,rank(A)= r, 可以分解为:A=UDV',其中U'U=V'V=I。在R中可以用函数scd()进行奇异值分解,例如:> A=matrix(1:18,3,6)> A  [,1] [,2] [,3] [,4] [,5] [,6][1,]   1   4   7   10   13   16[2,]   2   5   8   11   14   17[3,]   3   6   9   12   15   18> svd(A)$d[1] 4.589453e+01 1.640705e+00 3.627301e-16  $u          [,1]     [,2]     [,3][1,] -0.5290354 0.74394551 0.4082483[2,] -0.5760715 0.03840487 -0.8164966[3,] -0.6231077 -0.66713577 0.4082483$v          [,1]     [,2]     [,3][1,] -0.07736219 -0.7196003 -0.18918124[2,] -0.19033085 -0.5089325 0.42405898[3,] -0.30329950 -0.2982646 -0.45330031[4,] -0.41626816 -0.0875968 -0.01637004[5,] -0.52923682 0.1230711 0.64231130[6,] -0.64220548 0.3337389 -0.40751869> A.svd=svd(A)> A.svd$u%*%diag(A.svd$d)%*%t(A.svd$v)  [,1] [,2] [,3] [,4] [,5] [,6][1,]   1   4   7   10   13   16[2,]   2   5   8   11   14   17[3,]   3   6   9   12   15   18> t(A.svd$u)%*%A.svd$u            [,1]       [,2]       [,3][1,] 1.000000e+00 -1.169312e-16 -3.016793e-17[2,] -1.169312e-16 1.000000e+00 -3.678156e-17[3,] -3.016793e-17 -3.678156e-17 1.000000e+00> t(A.svd$v)%*%A.svd$v        [,1]       [,2]       [,3][1,] 1.000000e+00 8.248068e-17 -3.903128e-18[2,] 8.248068e-17 1.000000e+00 -2.103352e-17[3,] -3.903128e-18 -2.103352e-17 1.000000e+0012   矩阵QR分解A为m×n矩阵可以进行QR分解,A=QR,其中:Q'Q=I,在R中可以用函数qr()进行QR分解,例如:> A=matrix(1:16,4,4)> qr(A)$qr      [,1]     [,2]       [,3]       [,4][1,] -5.4772256 -12.7801930 -2.008316e+01 -2.738613e+01[2,] 0.3651484 -3.2659863 -6.531973e+00 -9.797959e+00[3,] 0.5477226 -0.3781696 2.641083e-15 2.056562e-15[4,] 0.7302967 -0.9124744 8.583032e-01 -2.111449e-16$rank[1] 2$qraux[1] 1.182574e+00 1.156135e+00 1.513143e+00 2.111449e-16$pivot[1] 1 2 3 4attr(,"class")[1] "qr"rank项返回矩阵的秩,qr项包含了矩阵Q和R的信息,要得到矩阵Q和R,可以用函数qr.Q()和qr.R()作用qr()的返回结果,例如:> qr.R(qr(A))      [,1]     [,2]       [,3]       [,4][1,] -5.477226 -12.780193 -2.008316e+01 -2.738613e+01[2,] 0.000000 -3.265986 -6.531973e+00 -9.797959e+00[3,] 0.000000   0.000000 2.641083e-15 2.056562e-15[4,] 0.000000   0.000000 0.000000e+00 -2.111449e-16> qr.Q(qr(A))      [,1]       [,2]     [,3]     [,4][1,] -0.1825742 -8.164966e-01 -0.4000874 -0.37407225[2,] -0.3651484 -4.082483e-01 0.2546329 0.79697056[3,] -0.5477226 -8.131516e-19 0.6909965 -0.47172438[4,] -0.7302967 4.082483e-01 -0.5455419 0.04882607> qr.Q(qr(A))%*%qr.R(qr(A))  [,1] [,2] [,3] [,4][1,]   1   5   9   13[2,]   2   6   10   14[3,]   3   7   11   15[4,]   4   8   12   16> t(qr.Q(qr(A)))%*%qr.Q(qr(A))        [,1]       [,2]       [,3]       [,4][1,] 1.000000e+00 -1.457168e-16 -6.760001e-17 -7.659550e-17[2,] -1.457168e-16 1.000000e+00 -4.269046e-17 7.011739e-17[3,] -6.760001e-17 -4.269046e-17 1.000000e+00 -1.596437e-16[4,] -7.659550e-17 7.011739e-17 -1.596437e-16 1.000000e+00> qr.X(qr(A))  [,1] [,2] [,3] [,4][1,]   1   5   9   13[2,]   2   6   10   14[3,]   3   7   11   15[4,]   4   8   12   1613   矩阵广义逆(Moore-Penrose)  n×m矩阵A+称为m×n矩阵A的Moore-Penrose逆,如果它满足下列条件:①   A A+A=A;②A+A A+= A+;③(A A+)H=A A+;④(A+A)H= A+A在R的MASS包中的函数ginv()可计算矩阵A的Moore-Penrose逆,例如:library(“MASS”)> A  [,1] [,2] [,3] [,4][1,]   1   5   9   13[2,]   2   6   10   14[3,]   3   7   11   15[4,]   4   8   12   16> ginv(A)    [,1]   [,2] [,3]   [,4][1,] -0.285 -0.1075 0.07 0.2475[2,] -0.145 -0.0525 0.04 0.1325[3,] -0.005 0.0025 0.01 0.0175[4,] 0.135 0.0575 -0.02 -0.0975验证性质1:> A%*%ginv(A)%*%A  [,1] [,2] [,3] [,4][1,]   1   5   9   13[2,]   2   6   10   14[3,]   3   7   11   15[4,]   4   8   12   16验证性质2:> ginv(A)%*%A%*%ginv(A)    [,1]   [,2] [,3]   [,4][1,] -0.285 -0.1075 0.07 0.2475[2,] -0.145 -0.0525 0.04 0.1325[3,] -0.005 0.0025 0.01 0.0175[4,] 0.135 0.0575 -0.02 -0.0975验证性质3:> t(A%*%ginv(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%*%ginv(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验证性质4:> t(ginv(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> ginv(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.714   矩阵Kronecker积  n×m矩阵A与h×k矩阵B的kronecker积为一个nh×mk维矩阵,在R中kronecker积可以用函数kronecker()来计算,例如:> A=matrix(1:4,2,2)> B=matrix(rep(1,4),2,2)> A  [,1] [,2][1,]   1   3[2,]   2   4> B  [,1] [,2][1,]   1   1[2,]   1   1> kronecker(A,B)  [,1] [,2] [,3] [,4][1,]   1   1   3   3[2,]   1   1   3   3[3,]   2   2   4   4[4,]   2   2   4   415   矩阵的维数  在R中很容易得到一个矩阵的维数,函数dim()将返回一个矩阵的维数,nrow()返回行数,ncol()返回列数,例如:  > A=matrix(1:12,3,4)> A  [,1] [,2] [,3] [,4][1,]   1   4   7   10[2,]   2   5   8   11[3,]   3   6   9   12> nrow(A)[1] 3> ncol(A)[1] 4

——————————————

往期精彩:

我造的假我自己打,Adobe推出“反PS”

微软删除人脸识别,除了隐私,更重要的可能是性别歧视与种族主义

亚马逊在中国失败,而中国却在亚马逊成功

640?wx_fmt=png



【本文地址】


今日新闻


推荐新闻


    CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3