跳转到内容

R 编程/数学

来自维基教科书,开放世界中的开放书籍
?Arithmetic
?Special

线性代数

[编辑 | 编辑源代码]

内积也称为点积标量积。它是逐项乘积之和。

> u <- rep(3,3)
> v <- 1:3
> u%*%v # the inner product
     [,1]
[1,]   18

外积也称为叉积向量积。它是两个向量元素的乘积得到的矩阵。

> v <- rep(3,3)
> u <- 1:3
> u%o%v # The outer product
     [,1] [,2] [,3]
[1,]    3    3    3
[2,]    6    6    6
[3,]    9    9    9

矩阵代数

[编辑 | 编辑源代码]

如果你想创建一个新的矩阵,一种方法是使用matrix()函数。你必须输入一个数据向量,行数和/或列数,最后你可以指定是否希望 R 按行或按列(默认选项)读取你的向量,使用byrow。你也可以使用向量组合cbind()rbind()。矩阵的维数可以使用dim()函数或可替代地使用nrow()ncol().

> matrix(data = NA, nrow = 5, ncol = 5, byrow = T)
> matrix(data = 1:15, nrow = 5, ncol = 5, byrow = T)
> v1 <- 1:5
> v2 <- 5:1
> cbind(v1,v2)
> rbind(v1,v2)
> dim(X)
> nrow(X)
> ncol(X)

一些特殊的矩阵

[编辑 | 编辑源代码]

单位矩阵在对角线上为 1,对角线外为 0。

  • eye()(matlab)
  • diag(1,nrow=10,ncol=10)
  • diag(rep(1,10))

J 矩阵全部为 1

  • ones()(matlab)

一个全为 0 的矩阵

  • zeros()(matlab)
> library(matlab)
> eye(3)
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1
> ones(3)
     [,1] [,2] [,3]
[1,]    1    1    1
[2,]    1    1    1
[3,]    1    1    1
> zeros(3) 
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

对角矩阵

> diag(3)

     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1

上三角矩阵

> round(upper.tri(matrix(1, n, n))) 

for n=3
     [,1] [,2] [,3]
[1,]    0    1    1
[2,]    0    0    1
[3,]    0    0    0

If you also need the diagonal of one's 

> round(upper.tri(matrix(1, 3, 3), diag = TRUE))

      [,1] [,2] [,3]
[1,]    1    1    1
[2,]    0    1    1
[3,]    0    0    1

下三角矩阵

与上三角矩阵相同,但使用 lower.tri 代替


  • 使用以下方法创建希尔伯特矩阵hilbert()(fUtilities)。

矩阵计算

[编辑 | 编辑源代码]
> b <- matrix(nrow = 2, ncol = 2, c(1, 2, 3, 4))
> a <- matrix(nrow = 2, ncol = 2, c(1, 0, 0, -1))
> a
     [,1] [,2]
[1,]    1    0
[2,]    0   -1
> b
     [,1] [,2]
[1,]    1    3
[2,]    2    4
> a%*%b
     [,1] [,2]
[1,]    1    3
[2,]   -2   -4
> b%*%a
     [,1] [,2]
[1,]    1   -3
[2,]    2   -4
> M <- matrix(rep(2,4),nrow = 2) 
> M
     [,1] [,2]
[1,]    2    2
[2,]    2    2
> I <- eye(2) 
> I
     [,1] [,2]
[1,]    1    0
[2,]    0    1
> I %x% M 
     [,1] [,2] [,3] [,4]
[1,]    2    2    0    0
[2,]    2    2    0    0
[3,]    0    0    2    2
[4,]    0    0    2    2
> library(fUtilities)
> kron(I,M)
     [,1] [,2] [,3] [,4]
[1,]    2    2    0    0
[2,]    2    2    0    0
[3,]    0    0    2    2
[4,]    0    0    2    2

矩阵转置

[编辑 | 编辑源代码]
  • 转置矩阵
> t(M)
     [,1] [,2] [,3]
[1,]    1    0    1
[2,]    0    1    2
[3,]    0    0    1

矩阵的迹和行列式

[编辑 | 编辑源代码]
  • 使用以下方法计算矩阵的迹tr()(fUtilities)
  • 使用以下方法返回矩阵的秩rk()(fBasics:)

矩阵求逆

[编辑 | 编辑源代码]
  • 使用以下方法对矩阵求逆solve()inv()(fUtilities)。我们也可以使用广义逆计算ginv()MASS 包中。
> M <- cbind(c(1,0,1),c(0,1,2),c(0,0,1))
> solve(M)
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]   -1   -2    1
> solve(M)%*%M
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1

求解线性方程

[编辑 | 编辑源代码]
> m=matrix(nrow=2,ncol=2,c(1,-.8,1,.2))
> m
     [,1] [,2]
[1,]  1.0  1.0
[2,] -0.8  0.2
> 
> l=matrix(c(1.0+25.0/18,25.0/18.0))
> l
         [,1]
[1,] 2.388889
[2,] 1.388889
> 
> k=solve(m,l)
> k
           [,1]
[1,] -0.9111111
[2,]  3.3000000
> 
> m%*%k          #checking the answer
         [,1]
[1,] 2.388889
[2,] 1.388889
>


特征值、特征向量和特征空间

[编辑 | 编辑源代码]
  • 特征值和特征向量
> eigen(M)
$values
[1] 1 1 1

$vectors
     [,1]          [,2]          [,3]
[1,]    0  2.220446e-16  0.000000e+00
[2,]    0  0.000000e+00  1.110223e-16
[3,]    1 -1.000000e+00 -1.000000e+00
  • 使用以下方法计算矩阵的范数norm()(fUtilities)。
  • 检查矩阵是否为正定矩阵isPositiveDefinite()(fUtilities)。
  • 使矩阵为正定矩阵makePositiveDefinite()(fUtilities)。
  • 计算行统计量和列统计量 (fUtilities)。
  • 提取矩阵的上部和下部triang()Triang()(fUtilities)。
  • 另请参阅 matrixmatlabmatrixcalcmatrixStats 包。

对数和指数

[编辑 | 编辑源代码]

我们有幂函数10^3 10**3 ,对数和指数log(2.71), log10(10),exp(1).

> 10^3 # exponent
[1] 1000
> 10**3 # exponent
[1] 1000
> exp(1) # exponential
[1] 2.718282
> log(2.71) # natural logarithm
[1] 0.9969486
> log10(1000) # base 10 logarithm
[1] 3
> log(1000,base = 10) # base 10 logarithm
[1] 3


多项式方程

[编辑 | 编辑源代码]

求解 ,其中 为给定数字,使用命令

> polyroot(c(n,...,b,a))

例如,要计算方程 的根,可以执行以下操作

> polyroot(c(-3,-5,2))
 [1] -0.5+0i  3.0-0i

解可以读为 .

另请参阅 polynommultipol

符号计算

[编辑 | 编辑源代码]

R 可以给出表达式的导数。你需要使用expression()函数将你的函数转换为表达式。否则你会收到错误消息。

以下是一些示例

> D(expression(x^n),"x")
x^(n - 1) * n
> D(expression(exp(a*x)),"x")
exp(a * x) * a
> D(expression(1/x),"x")
-(1/x^2)
> D(expression(x^3),"x")
3 * x^2
> D(expression(pnorm(x)),"x")
dnorm(x)
> D(expression(dnorm(x)),"x")
-(x * dnorm(x))

数值逼近

[编辑 | 编辑源代码]
  • numDeriv

R 可以执行一维积分。例如,我们可以在正态分布的密度之间积分

> integrate(dnorm,-Inf,Inf)
1 with absolute error < 9.4e-05
> integrate(dnorm,-1.96,1.96)
0.9500042 with absolute error < 1.0e-11
> integrate(dnorm,-1.64,1.64)
0.8989948 with absolute error < 6.8e-14
# we can also store the result in an object
> ci90 <- integrate(dnorm,-1.64,1.64)
> ci90$value
[1] 0.8989948
> integrate(dnorm,-1.64,1.64)$value
[1] 0.8989948

查看 adapt 包以了解多元积分。

> library(adapt)
> ?adapt
> ir2pi <- 1/sqrt(2*pi)
> fred <- function(z) { ir2pi^length(z) * exp(-0.5 * sum(z * z))}
> 
> adapt(2, lo = c(-5,-5), up = c(5,5), functn = fred)
       value       relerr       minpts       lenwrk        ifail 
    1.039222 0.0007911264          231           73            0 
> adapt(2, lo = c(-5,-5), up = c(5,5), functn = fred, eps = 1e-4)
       value       relerr       minpts       lenwrk        ifail 
    1.000237 1.653498e-05          655          143            0 
> adapt(2, lo = c(-5,-5), up = c(5,5), functn = fred, eps = 1e-6)
      value      relerr      minpts      lenwrk       ifail 
   1.000039 3.22439e-07        1719         283           0
  • 另请参阅integrate.gh()ecoreg 包中。
  • n 个数字中长度为 k 的组合数:
> choose(100, 5)
[1] 75287520
  • 并集和交集
> union(1:10, 5:7)
[1]  1  2  3  4  5  6  7  8  9 10
> intersect(1:10, 5:7)
[1] 5 6 7

阶乘函数

[编辑 | 编辑源代码]

factorial返回整数的阶乘。这也可以使用以下方法计算prod()(乘积) 应用于 1 到目标数字之间的整数向量。

> factorial(3)
[1] 6
> prod(1:3)
[1] 6

请注意,根据惯例 .factorial()在 0 返回 1。情况并非如此prod()函数。

> factorial(0)
[1] 1
> prod(0)
[1] 0

阶乘数可能非常大,无法计算高值。

> factorial(170)
[1] 7.257416e+306
> factorial(171)
[1] Inf
Message d'avis :
In factorial(171) : value out of range in 'gammafn'

模运算和欧几里德除法

[编辑 | 编辑源代码]
  • 模运算和整数除法(即欧几里德除法)
> 5%%2
[1] 1
>5%/%2
[1] 2

注意:R 受 问题 的影响,该问题与非整数和欧几里德除法有关。

> .5%/%.1 # we get 4 instead of 5
[1] 4
> .5%%.1 # we get .1 instead of 0
[1] 0.1
  • pi常数
  • cos(), sin(), tan()三角函数。

符号计算

[编辑 | 编辑源代码]

rSymPy (rsympy) 在 R 中提供 sympy (链接) 函数。

如果您想做更多符号计算,请参阅 Maxima[1]、SAGE[2]、Mathematica[3]

另请参阅

[编辑 | 编辑源代码]

以下命令提供有关与贝塔函数和伽马函数相关的特殊数学函数的帮助。

?Special

参考资料

[编辑 | 编辑源代码]
  1. Maxima 是开源软件 http://maxima.sourceforge.net/
  2. SAGE 是一个包含 R 和 Maxima 的开源包:http://www.sagemath.org/
  3. Mathematica 不是开源软件 http://www.wolfram.com/products/mathematica/index.html
华夏公益教科书