Golden State: 2023年2月

［R］t分布におけるパーセント点を求める

> # 自由度5のt分布における上側5パーセント点
> qt(1 - 0.05, 5)
[1] 2.015048
> # 自由度5のt分布における両側5パーセント点
> qt(1 - 0.05 / 2, 5)
[1] 2.570582

2023年2月28日 (火) R（統計学） | 固定リンク | コメント (0)

［R］t分布におけるパーセント点（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.138）

> # 自由度5のt分布における上側5パーセント点
> qt(1 - 0.05, 5)
[1] 2.015048
> # 自由度5のt分布における両側5パーセント点
> qt(1 - 0.05 / 2, 5)
[1] 2.570582

2023年2月28日 (火) R（本の計算を再現） | 固定リンク | コメント (0)

［R］F分布におけるp値

> # 自由度5,7のF分布における、確率変数が2のときの下側p値
> pf(2, 5, 7)
[1] 0.8043268
> # 自由度5,7のF分布における、確率変数が2のときの上側p値
> pf(2, 5, 7, lower.tail = FALSE)
[1] 0.1956732

2023年2月27日 (月) R（統計学） | 固定リンク | コメント (0)

［R］F分布におけるp値（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.206）

> # 自由度5,7のF分布における、確率変数が2のときの下側p値
> pf(2, 5, 7)
[1] 0.8043268
> # 自由度5,7のF分布における、確率変数が2のときの上側p値
> pf(2, 5, 7, lower.tail = FALSE)
[1] 0.1956732

2023年2月27日 (月) R（本の計算を再現） | 固定リンク | コメント (0)

［R］F分布におけるパーセント点を求める

> # 自由度5,7のF分布における、下側5パーセント点
> qf(0.05, 5, 7)
[1] 0.2050915
> # 自由度5,7のF分布における、上側5パーセント点
> qf(0.05, 5, 7, lower.tail = FALSE)
[1] 3.971523
> qf(1 - 0.05, 5, 7)
[1] 3.971523

2023年2月26日 (日) R（統計学） | 固定リンク | コメント (0)

［R］F分布におけるパーセント点（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、pp.205-206）

> # 自由度5,7のF分布における、下側5パーセント点
> qf(0.05, 5, 7)
[1] 0.2050915
> # 自由度5,7のF分布における、上側5パーセント点
> qf(0.05, 5, 7, lower.tail = FALSE)
[1] 3.971523
> qf(1 - 0.05, 5, 7)
[1] 3.971523

2023年2月26日 (日) R（本の計算を再現） | 固定リンク | コメント (0)

［R］χ^2分布におけるp値を求める

> # 自由度5のχ^2分布における、確率変数が12のときの下側p値
> pchisq(12, 5)
[1] 0.9652122
> # 自由度5のχ^2分布における、確率変数が12のときの上側p値
> pchisq(12, 5, lower.tail = FALSE)
[1] 0.03478778

2023年2月24日 (金) R（統計学） | 固定リンク | コメント (0)

［R］χ^2分布におけるp値（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.153）

> # 自由度5のχ^2分布における、確率変数が12のときの下側p値
> pchisq(12, 5)
[1] 0.9652122
> # 自由度5のχ^2分布における、確率変数が12のときの上側p値
> pchisq(12, 5, lower.tail = FALSE)
[1] 0.03478778

2023年2月24日 (金) R（本の計算を再現） | 固定リンク | コメント (0)

［R］χ^2分布におけるパーセント点を求める

> # 自由度5のχ^2分布における、下側5パーセント点
> qchisq(0.05, 5)
[1] 1.145476
> # 自由度5のχ^2分布における、上側5パーセント点
> qchisq(0.05, 5, lower.tail = FALSE)
[1] 11.0705
> qchisq(1 - 0.05, 5)
[1] 11.0705

2023年2月15日 (水) R（統計学） | 固定リンク | コメント (0)

［R］χ^2分布におけるパーセント点（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.152）

> # 自由度5のχ^2分布における、下側5パーセント点
> qchisq(0.05, 5)
[1] 1.145476
> # 自由度5のχ^2分布における、上側5パーセント点
> qchisq(0.05, 5, lower.tail = FALSE)
[1] 11.0705
> qchisq(1 - 0.05, 5)
[1] 11.0705

2023年2月15日 (水) R（本の計算を再現） | 固定リンク | コメント (0)

［R］連立一次方程式を解く

solve関数を使う。以下の連立一次方程式を解く（解はx=3,y=2）。

 x + 2y = 7
3x - 4y = 1

計算する。

> mxaa <- matrix(c(1, 2, 3, -4), 2, 2, byrow = TRUE)
> mxy <- matrix(c(7, 1), 2, 1)
> print(mxaa)
     [,1] [,2]
[1,]    1    2
[2,]    3   -4
> print(mxy)
     [,1]
[1,]    7
[2,]    1
> solve(mxaa, mxy)
     [,1]
[1,]    3
[2,]    2

同じく以下の連立一次方程式を解く（解はx=2,y=-3,z=-14）。

 3x - 3y + z =  1
 3x + 2y     =  0
-1x - 5y + z = -1

計算する。

> mxaa <- matrix(c(3, 3, -1, -3, 2, -5, 1, 0, 1), 3, 3)
> mxy <- matrix(c(1, 0, -1), 3, 1)
> print(mxaa)
     [,1] [,2] [,3]
[1,]    3   -3    1
[2,]    3    2    0
[3,]   -1   -5    1
> print(mxy)
     [,1]
[1,]    1
[2,]    0
[3,]   -1
> solve(mxaa, mxy)
     [,1]
[1,]    2
[2,]   -3
[3,]  -14

2023年2月13日 (月) R（行列） | 固定リンク | コメント (0)

［R］一般逆行列を求める

MASSパッケージのginv関数を使う。

> library(MASS)
> mx <- matrix(1:16, 4, 4, byrow = TRUE)
> print(mx)
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    5    6    7    8
[3,]    9   10   11   12
[4,]   13   14   15   16
> det(mx)
[1] 4.733165e-30
> solve(mx)
 solve.default(mx) でエラー: 
   システムは数値的に特異です: 条件数の逆数 = 2.95756e-18 
> ginv(mx)
        [,1]    [,2]    [,3]    [,4]
[1,] -0.2850 -0.1450 -0.0050  0.1350
[2,] -0.1075 -0.0525  0.0025  0.0575
[3,]  0.0700  0.0400  0.0100 -0.0200
[4,]  0.2475  0.1325  0.0175 -0.0975

一般逆行列の定義を満たしているか否か確認。

> mxi <- ginv(mx)
> mx %*% mxi %*% mx
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    5    6    7    8
[3,]    9   10   11   12
[4,]   13   14   15   16

当然、正方行列でなくとも求めることができる。

> mx <- matrix(1:8, 2, 4, byrow = TRUE)
> print(mx)
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    5    6    7    8
> det(mx)
 determinant.matrix(x, logarithm = TRUE, ...) でエラー: 
   'x' は正方行列でなければなりません 
> solve(mx)
 solve.default(mx) でエラー:  'a' (2 x 4) は正方行列である必要があります 
> ginv(mx)
       [,1]          [,2]
[1,] -0.550  2.500000e-01
[2,] -0.225  1.250000e-01
[3,]  0.100 -2.081668e-17
[4,]  0.425 -1.250000e-01
> mxi <- ginv(mx)
> mx %*% mxi %*% mx
     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    5    6    7    8

ginv関数はムーア・ペンローズ一般逆行列を数値的に求める関数のため、求まった一般逆行列は近似値であることに注意。ムーア・ペンローズ一般逆行列の定義の4つの式について、比較をした例は以下のとおり。すべての要素がTRUEにならないことがわかる。

> mx %*% mxi %*% mx == mx
      [,1]  [,2]  [,3] [,4]
[1,] FALSE FALSE  TRUE TRUE
[2,] FALSE FALSE FALSE TRUE
> mxi %*% mx %*% mxi == mxi
      [,1]  [,2]
[1,] FALSE FALSE
[2,] FALSE FALSE
[3,] FALSE FALSE
[4,] FALSE FALSE
> t(mx %*% mxi) == mx %*% mxi
      [,1]  [,2]
[1,]  TRUE FALSE
[2,] FALSE  TRUE
> t(mxi %*% mx) == mxi %*% mx
      [,1]  [,2]  [,3]  [,4]
[1,]  TRUE FALSE FALSE FALSE
[2,] FALSE  TRUE FALSE FALSE
[3,] FALSE FALSE  TRUE FALSE
[4,] FALSE FALSE FALSE  TRUE

2023年2月 9日 (木) R（行列） | 固定リンク | コメント (0)

［R］逆行列を求める

solve関数を使う。

> mx <- matrix(3, 1, 1)
> print(mx)
     [,1]
[1,]    3
> solve(mx)
          [,1]
[1,] 0.3333333
> mx <- matrix(1:4, 2, 2, byrow = TRUE)
> print(mx)
     [,1] [,2]
[1,]    1    2
[2,]    3    4
> solve(mx)
     [,1] [,2]
[1,] -2.0  1.0
[2,]  1.5 -0.5
> mx <- matrix(c(3, 3, -1, -3, 2, -5, 1, 0, 1), 3, 3)
> print(mx)
     [,1] [,2] [,3]
[1,]    3   -3    1
[2,]    3    2    0
[3,]   -1   -5    1
> solve(mx)
     [,1] [,2] [,3]
[1,]  1.0   -1 -1.0
[2,] -1.5    2  1.5
[3,] -6.5    9  7.5

正則ではない行列（行列式の値が0）の場合は、solve関数はエラーを返す。

> mx <- matrix(1:9, 3, 3)
> print(mx)
     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6    9
> det(mx)
[1] 0
> solve(mx)
 solve.default(mx) でエラー: 
  Lapack routine dgesv: システムは正確に特異です: U[3,3] = 0

2023年2月 7日 (火) R（行列） | 固定リンク | コメント (0)

［R］転置行列を求める

t関数を使う。

> mx <- matrix(1:9, 3, 3, byrow = TRUE)
> print(mx)
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6
[3,]    7    8    9
> t(mx)
     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6    9
> mx <- matrix(1:6, 2, 3, byrow = TRUE)
> print(mx)
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6
> t(mx)
     [,1] [,2]
[1,]    1    4
[2,]    2    5
[3,]    3    6

2023年2月 6日 (月) R（行列） | 固定リンク | コメント (0)

［Python］行列の掛け算

dotメソッドを使う。最後の例のとおり、掛けられる行列の列数と掛ける行列の行数が同じでなければ計算することはできない。

>>> mx1 = np.array([[-1, 0, 3], [8, 1, -5]])
>>> mx2 = np.array([[3, 4, 0], [-2, 1, 2]])
>>> mx3 = np.array([[4, -3, 7], [2, 0, -1], [1, 5, 0]])
>>> print(mx1)
[[-1  0  3]
 [ 8  1 -5]]
>>> print(mx2)
[[ 3  4  0]
 [-2  1  2]]
>>> print(mx3)
[[ 4 -3  7]
 [ 2  0 -1]
 [ 1  5  0]]
>>> np.dot(mx1, mx3)
array([[ -1,  18,  -7],
       [ 29, -49,  55]])
>>> np.dot(mx2, mx3)
array([[ 20,  -9,  17],
       [ -4,  16, -15]])
>>> np.dot(mx1, mx2)
Traceback (most recent call last):
  File "", line 1, in 
  File "<__array_function__ internals>", line 5, in dot
ValueError: shapes (2,3) and (2,3) not aligned: 3 (dim 1) != 2 (dim 0)

2023年2月 4日 (土) Python | 固定リンク | コメント (0)

［R］行列式を求める

det関数を使う。

> mx <- matrix(3, 1, 1)
> print(mx)
     [,1]
[1,]    3
> det(mx)
[1] 3
> mx <- matrix(1:4, 2, 2, byrow = TRUE)
> print(mx)
     [,1] [,2]
[1,]    1    2
[2,]    3    4
> det(mx)
[1] -2
> mx <- matrix(c(0, -4, 0, -1, 3, 6, -2, 0, 2, 1, 3, 0, -1, 5, -1, 2), 4, 4)
> print(mx)
     [,1] [,2] [,3] [,4]
[1,]    0    3    2   -1
[2,]   -4    6    1    5
[3,]    0   -2    3   -1
[4,]   -1    0    0    2
> det(mx)
[1] 28

2023年2月 4日 (土) R（行列） | 固定リンク | コメント (0)

［R］等差数列を作成する

seqコマンドを使う。引数は1～3つ与えることができ、それぞれ以下のようになる。最後の例のとおり、最後に公差を加えて末項に一致しない場合は、その1つ前の項までが作成される。

seq(初項)
seq(初項, 末項) ※公差は1
seq(初項, 末項, 公差) ※末項に届かないときはその1つ前の項まで

> seq(3)
[1] 1 2 3
> seq(2, 5)
[1] 2 3 4 5
> seq(3, 12, 4)
[1]  3  7 11

負数や実数の指定も可能。

> seq(-1, -11, -2)
[1]  -1  -3  -5  -7  -9 -11
> seq(1.2, 2.8, 0.4)
[1] 1.2 1.6 2.0 2.4 2.8

2023年2月 3日 (金) R（数学） | 固定リンク | コメント (0)

［R］行列の掛け算

%*% 演算子を使う。最後の例のとおり、掛けられる行列の列数と掛ける行列の行数が同じでなければ計算することはできない。

> mx1 <- matrix(c(-1, 8, 0, 1, 3, -5), 2, 3)
> mx2 <- matrix(c(3, -2, 4, 1, 0, 2), 2, 3)
> mx3 <- matrix(c(4, 2, 1, -3, 0, 5, 7, -1, 0), 3, 3)
> print(mx1)
     [,1] [,2] [,3]
[1,]   -1    0    3
[2,]    8    1   -5
> print(mx2)
     [,1] [,2] [,3]
[1,]    3    4    0
[2,]   -2    1    2
> print(mx3)
     [,1] [,2] [,3]
[1,]    4   -3    7
[2,]    2    0   -1
[3,]    1    5    0
> mx1 %*% mx3
     [,1] [,2] [,3]
[1,]   -1   18   -7
[2,]   29  -49   55
> mx2 %*% mx3
     [,1] [,2] [,3]
[1,]   20   -9   17
[2,]   -4   16  -15
> mx1 %*% mx2
 mx1 %*% mx2 でエラー:  適切な引数ではありません

2023年2月 2日 (木) R（行列） | 固定リンク | コメント (0)

Golden State

パソコンで数値計算を行うためのクックブック。R言語などの便利帳、備忘録、よく使う機能のメモ。

2023年2月28日 (火)

［R］t分布におけるパーセント点を求める

［R］t分布におけるパーセント点（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.138）

2023年2月27日 (月)

［R］F分布におけるp値

［R］F分布におけるp値（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.206）

2023年2月26日 (日)

［R］F分布におけるパーセント点を求める

［R］F分布におけるパーセント点（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、pp.205-206）

2023年2月24日 (金)

［R］χ^2分布におけるp値を求める

［R］χ^2分布におけるp値（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.153）

2023年2月15日 (水)

［R］χ^2分布におけるパーセント点を求める

［R］χ^2分布におけるパーセント点（「まずはこの一冊から意味が分かる統計解析」（ペレ出版）、p.152）

2023年2月13日 (月)

［R］連立一次方程式を解く

2023年2月 9日 (木)

［R］一般逆行列を求める

2023年2月 7日 (火)

［R］逆行列を求める

2023年2月 6日 (月)

［R］転置行列を求める

2023年2月 4日 (土)

［Python］行列の掛け算

［R］行列式を求める

2023年2月 3日 (金)

［R］等差数列を作成する

2023年2月 2日 (木)

［R］行列の掛け算

■

■■

■■■

カテゴリー