跳轉到內容

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
華夏公益教科書