R语言——矩阵运算

R语言的矩阵运算

创建矩阵向量；矩阵加减，乘积；矩阵的逆；行列式的值；特征值与特征向量；QR分解；奇异值分解；广义逆；backsolve与fowardsolve函数；取矩阵的上下三角元素；向量化算子等。

1、创建向量

2、创建矩阵

在R中可以用函数matrix()来创建一个矩阵。

3、矩阵转置

A为m×n矩阵，求A’的转置矩阵在R中可用函数t()，例如：

4、矩阵相加减

5、数与矩阵相乘

A为m×n矩阵，c>0，在R中求cA可用符号：“*”，例如：

6、矩阵相乘

A为m×n矩阵，B为n×k矩阵，在R中求AB可用符号：“％*％”，例如：

若A为n×m矩阵，要得到A’B，可用函数crossprod()，该函数计算结果与t(A)%*%B相同，但是效率更高。例如：

矩阵Hadamard积：若A={aij}m×n, B={bij}m×n, 则矩阵的Hadamard积定义为：
A⊙B={aij bij }m×n,R中Hadamard积可以直接运用运算符“*”例如：

7、矩阵对角元素相关运算

例如要取一个方阵的对角元素，

8、矩阵求逆

矩阵求逆可用函数solve()，应用solve(a, b)运算结果是解线性方程组ax = b，若b缺省，则系统默认为单位矩阵，因此可用其进行矩阵求逆，例如：

9、矩阵的特征值和特征向量

矩阵A的谱分解为A=UΛU’,其中Λ是由A的特征值组成的对角矩阵，U的列为A的特征值对应的特征向量，在R中可以用函数eigen()函数得到U和Λ，

10、矩阵的Choleskey分解

对于正定矩阵A，可对其进行Choleskey分解，即：A=P’P，其中P为上三角矩阵，在R中可以用函数chol()进行Choleskey分解，例如：

11、矩阵的奇异值分解

A为m×n矩阵，rank(A)= r, 可以分解为：A=UDV’,其中U’U=V’V=I。在R中可以用函数scd()进行奇异值分解，例如：

12、矩阵QR值分解

A为m×n矩阵可以进行QR分解，A=QR，其中：Q’Q＝I，在R中可以用函数qr()进行QR分解，例如：

rank项返回矩阵的秩，qr项包含了矩阵Q和R的信息，要得到矩阵Q和R，可以用函数qr.Q()和qr.R()作用qr()的返回结果，例如：

13、矩阵的广义逆

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逆，例如：

14 矩阵Kronecker积

n×m矩阵A与h×k矩阵B的kronecker积为一个nh×mk维矩阵，
在R中kronecker积可以用函数kronecker()来计算，例如：

15 矩阵的维数

在R中很容易得到一个矩阵的维数，函数dim()将返回一个矩阵的维数，nrow()返回行数，ncol()返回列数，例如：

16 矩阵的行和、列和、行平均与列平均

在R中很容易求得一个矩阵的各行的和、平均数与列的和、平均数，例如：

上述关于矩阵行和列的操作，还可以使用apply()函数实现。

args(apply)
function (X, MARGIN, FUN, …)
其中：x为矩阵，MARGIN用来指定是对行运算还是对列运算，MARGIN＝1表示对行运算，MARGIN＝2表示对列运算，FUN用来指定运算函数, …用来给定FUN中需要的其它的参数，例如：
apply(A,1,sum)
[1] 22 26 30
apply(A,1,mean)
[1] 5.5 6.5 7.5
apply(A,2,sum)
[1] 6 15 24 33
apply(A,2,mean)
[1] 2 5 8 11
apply()函数功能强大，我们可以对矩阵的行或者列进行其它运算，例如：
计算每一列的方差
A=matrix(rnorm(100),20,5)
apply(A,2,var)
[1] 0.4641787 1.4331070 0.3186012 1.3042711 0.5238485
apply(A,2,function(x,a)xa,a=2)
[,1] [,2] [,3] [,4]
[1,] 2 8 14 20
[2,] 4 10 16 22
[3,] 6 12 18 24
注意：apply(A,2,function(x,a)xa,a=2)与A*2效果相同，此处旨在说明如何应用alpply函数。

17 矩阵X’X的逆

在统计计算中，我们常常需要计算这样矩阵的逆，如OLS估计中求系数矩阵。R中的包“strucchange”提供了有效的计算方法。

args(solveCrossprod)
function (X, method = c(“qr”, “chol”, “solve”))
其中：method指定求逆方法，选用“qr”效率最高，选用“chol”精度最高，选用“slove”与slove(crossprod(x,x))效果相同，例如：
A=matrix(rnorm(16),4,4)
solveCrossprod(A,method=”qr”)
[,1] [,2] [,3] [,4]
[1,] 0.6132102 -0.1543924 -0.2900796 0.2054730
[2,] -0.1543924 0.4779277 0.1859490 -0.2097302
[3,] -0.2900796 0.1859490 0.6931232 -0.3162961
[4,] 0.2054730 -0.2097302 -0.3162961 0.3447627
solveCrossprod(A,method=”chol”)
[,1] [,2] [,3] [,4]
[1,] 0.6132102 -0.1543924 -0.2900796 0.2054730
[2,] -0.1543924 0.4779277 0.1859490 -0.2097302
[3,] -0.2900796 0.1859490 0.6931232 -0.3162961
[4,] 0.2054730 -0.2097302 -0.3162961 0.3447627
solveCrossprod(A,method=”solve”)
[,1] [,2] [,3] [,4]
[1,] 0.6132102 -0.1543924 -0.2900796 0.2054730
[2,] -0.1543924 0.4779277 0.1859490 -0.2097302
[3,] -0.2900796 0.1859490 0.6931232 -0.3162961
[4,] 0.2054730 -0.2097302 -0.3162961 0.3447627
solve(crossprod(A,A))
[,1] [,2] [,3] [,4]
[1,] 0.6132102 -0.1543924 -0.2900796 0.2054730
[2,] -0.1543924 0.4779277 0.1859490 -0.2097302
[3,] -0.2900796 0.1859490 0.6931232 -0.3162961
[4,] 0.2054730 -0.2097302 -0.3162961 0.3447627

18 取矩阵的上、下三角部分

在R中，我们可以很方便的取到一个矩阵的上、下三角部分的元素，函数lower.tri()和函数upper.tri()提供了有效的方法。

args(lower.tri)
function (x, diag = FALSE)
函数将返回一个逻辑值矩阵，其中下三角部分为真，上三角部分为假，选项diag为真时包含对角元素，为假时不包含对角元素。upper.tri()的效果与之孑然相反。例如：
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
lower.tri(A)
[,1] [,2] [,3] [,4]
[1,] FALSE FALSE FALSE FALSE
[2,] TRUE FALSE FALSE FALSE
[3,] TRUE TRUE FALSE FALSE
[4,] TRUE TRUE TRUE FALSE
lower.tri(A,diag=T)
[,1] [,2] [,3] [,4]
[1,] TRUE FALSE FALSE FALSE
[2,] TRUE TRUE FALSE FALSE
[3,] TRUE TRUE TRUE FALSE
[4,] TRUE TRUE TRUE TRUE
upper.tri(A)
[,1] [,2] [,3] [,4]
[1,] FALSE TRUE TRUE TRUE
[2,] FALSE FALSE TRUE TRUE
[3,] FALSE FALSE FALSE TRUE
[4,] FALSE FALSE FALSE FALSE
upper.tri(A,diag=T)
[,1] [,2] [,3] [,4]
[1,] TRUE TRUE TRUE TRUE
[2,] FALSE TRUE TRUE TRUE
[3,] FALSE FALSE TRUE TRUE
[4,] FALSE FALSE FALSE TRUE
A[lower.tri(A)]=0
A
[,1] [,2] [,3] [,4]
[1,] 1 5 9 13
[2,] 0 6 10 14
[3,] 0 0 11 15
[4,] 0 0 0 16
A[upper.tri(A)]=0
A
[,1] [,2] [,3] [,4]
[1,] 1 0 0 0
[2,] 2 6 0 0
[3,] 3 7 11 0
[4,] 4 8 12 16

19 backsolve&fowardsolve函数

20 row()与col()函数

21 行列式的值

22 向量化算子

23 时间序列的滞后值