Momozowy

Rが好きだ | Momozowy

Rが好きだ

Posted: 3 minutes ago

R が好きだ。

R 言語が好きだ。

https://ja.wikipedia.org/wiki/R%E8%A8%80%E8%AA%9E

静的型付けではないし、実行速度も遅い。その上メモリも食う。統計解析に向いていると言われるが、同じように言われる Python の方が人気だ。

https://redmonk.com/sogrady/2024/09/12/language-rankings-6-24/

R の何がいいのか。

余計なことを考えず、余計なことを書かず、サッと実行できる。

基本的な構文

代入

まず、変数 x への代入はこう書ける。

x <- 1

いきなり他の言語と違って面食らったかもしれない。しかしよく見てほしい。変数 x に、左向きの矢印 <- で 1 を押し込んでいる、というふうに見えないだろうか。

他のいくつかの言語では、= を使って代入を行うと思う。

x = 1

等号である。x と 1 が等しい。つまり x は 1 ということが分かる。

x = x + 1

ちょっと待て。これはなんだ？ x と x + 1 が等しい？

Rでも等号を使った代入ができるが一部制限される。素直に <- を使った方がいい。

x <- x + 1

x に、x + 1 を押し込んでいる。なるほど。他の分野で使われている意味とも競合しない。直感的でわかりやすい。私は R が好きだ。

押し込むなら右向きだっていいじゃないか。

1 -> a
2 -> b
a -> b
a + b -> c

そうだろう。実際これも実行できる。結果として a も b も値は 1 だし、c の値は 2 だ。とはいえ、このように書きたくなるケースは稀であると思われるし、3行目なんかは C 言語のアロー演算子と混乱するかもしれない。

じゃあこれは？

a = b <- 3 -> c

繰り返し

例えば、ある数列 x を標準化するとき、数列の平均 ave と標準偏差 sd_x が与えられていれば、

y <- (x - ave) / sd_x

と書ける。注目すべきは、x, y がベクトル（数列）で、ave, sd_x はスカラーであることだ。ベクトルとスカラーの計算であるが、x - ave、これだけで数列 x の"各要素を" ave の値で減算してくれる。数学的にはスカラーをベクトルに変換してから計算する必要があるが、その操作を隠蔽してくれている。同様に、vector / sd_x は vector の各要素を sd_x で除算してくれる。これはRではリサイクリングと呼ばれている。Numpy のブロードキャストのようなものだ。

もちろん、for 文を使い x の各要素を取り出して計算することもできる。しかし、行数で言えば 3 倍だ。

y <- NULL
for (xx in x)
    y <- c(y, (xx - ave) / sd_x)

（関数c は、combine。複数の引数をとり、それらを繋げた新しいベクトル（R における一次元配列）を作る）

冒頭で、R は実行速度が遅いと言った。いいじゃないか。for 文で書けば3行のものを1行で書ける。コードを書く速度は 3 倍だ。また、他の言語を使っている時のように「for 文を回すか……待てこの言語は map があったか？いやないかでは forEach は？ない？ for 文を回すか……」などと考えなくていい。思いつくままにコードを書ける。私は R が好きだ。

R は実行速度が遅い？ちょっと待ってほしい。そのコードは for 文を乱用していないだろうか。ほとんどの場合、for 文を回すよりリサイクリングを使用した方が速い。R はそのように最適化されている。

次のコードを実行してみてほしい。

x <- iris$Sepal.Length
ave <- mean(x)
sd_x <- sd(x)

system.time(
    for (i in 1 : 1000000)
    {
        y <- (x - ave) / sd_x
    }
)

system.time(
    for (i in 1 : 1000000)
    {
        y <- NULL
        for (xx in x )
            y <- c(y, (xx - ave) / sd_x)
    }
)

私の環境では、リサイクリングの方が 150 倍速かった。簡単に書いた方が速い。私は R が好きだ。

少し実践的な話

Iris（アヤメ）という植物のデータを見てみたい。head(iris) を実行すると、iris データの先頭数行が表示される。

   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
 1          5.1         3.5          1.4         0.2  setosa
 2          4.9         3.0          1.4         0.2  setosa
 3          4.7         3.2          1.3         0.2  setosa
 4          4.6         3.1          1.5         0.2  setosa
 5          5.0         3.6          1.4         0.2  setosa
 6          5.4         3.9          1.7         0.4  setosa

そう。サンプルデータが組み込みで存在するのだ。私は R が好きだ。

Sepal と Petal、それぞれの Length と Width のデータが数値型で存在し、Species は文字型でおそらく何種類かあるのだろう。

このようなデータは、R ではデータフレームと呼ばれる。行列に見立てて、1行3列目、Petal.Length の値を取り出したい時は、いくつかの書き方ができる。

iris[1, 3]
iris[1, "Petal.Length"]
iris$Petal.Length[1]

1つ目は、単純に何行何列目かを指定する書き方。 2つ目は、何列目かの代わりに列名 "Petal.Length" を指定する書き方。 3つ目は、データフレームから "Petal.Length" の列を取り出し、先頭から何番目のデータか、を指定する書き方。

列名で指定できるのは嬉しい。決して、3列目を取り出したいのではなく、Petal.Length の値を取り出したいのである。目的のデータが何列目か数えなくて済むし、可読性が上がる。私は R が好きだ。

言い忘れていたが、R では配列などの先頭要素のインデックスは 1 である。ポインタのアドレスからどのくらい進むか、ではない。何個目のデータか、である。長さ N の配列の最後は N 番目である。その方が直感的だし、今はメモリのことは考えずデータに集中したい。 R は低レイヤを隠蔽してくれる。私は R が好きだ。

Iris に話を戻して、どんな Species があるか見てみよう。 levels 関数は、データの因子を抽出してくれる。

levels(iris[, 5])
levels(iris[, "Species"])
levels(iris$Species)

先ほどと同じように、上から列番号を指定する書き方、列名を指定する書き方、列名で列ごと取り出す書き方である。上 2 つでは、行番号を指定していない。そうすると全ての行、すなわち列全体を取り出せる。その結果、"setosa", "versicolor", "virginica" の3種類のアヤメに関するデータであることが分かる。

ここで思う。種ごとに Petal.Length の平均を知りたい。こんな時どうするか。まず、種ごとの Petal.Length を取り出す必要がある。その後 mean 関数を使えば、for 文を使うことなく計算できる。どうやって種ごとのデータを取り出すか？行方向に for を回して、if で条件分岐しないといけないのか？

そんなことはない。ここでもリサイクリングを使える。まずは "setosa" だけ見てみよう。その行が "setosa" かどうかは、比較演算子 == を使えばいいだろう。

iris$Species == "setosa"

これもリサイクリングされるのだろうか？もちろんである。

結果として、次のような bool 値のベクトルが手に入る。

TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE
TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE
TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE FALSE FALSE FALSE FALSE
...

これをどう使うのか。そのまま添字として渡せば、TRUE の位置にあるものだけを取り出せる。データフレームから値を取り出す4つ目の方法である（データフレームだけでなく、行列やベクトルでも同様である）。

まず、長さ 3 のベクトルで簡単な例で見てみる。

c(11, 12, 13)[c(FALSE, TRUE, TRUE)]

この例では、cによって3つの数値をつなげてベクトルにし、その添字として bool 型のベクトルを指定している。この時、TRUE の位置にある値、12 と 13 が取り出される。

そうすると、次のコードで "setosa" の行のみ取り出すことができると分かる。

iris[iris$Species == "setosa", ]

今考えているのは Petal.Length の平均である。列を次のように指定し、これを mean 関数に渡せばいい。

mean(iris[iris$Species == "setosa", "Petal.Length"])

時間を測ってみよう。

system.time(
    for (i in 1:10000)
        mean(iris[iris$Species == "setosa", "Petal.Length"])
)

system.time(
    for (i in 1:10000)
    {
        setosa.petal.length <- NULL
        for (index in 1 : nrow(iris))
            if (iris[index, "Species"] == "setosa")
                setosa.petal.length <- c(setosa.petal.length, iris[index, "Petal.Length"])
        mean(setosa.petal.length)
    }
)

やはり 150 倍ほど速いし、コードは 1/5 だ。iris の行数が 150 行だから 150 倍なのか？

"setosa" 以外はどうすればいいだろうか。今度こそ for 文で回すしかないのか？いや、lapply 関数が使える。これはリスト（一次元配列のようなもの）用の apply 関数であり、第一引数にリストを取り、そのリストの各要素に第二引数の関数を"適用"する。

lapply(levels(iris$Species), function(sp) {
    mean(iris[iris$Species == sp, "Petal.Length"])
})

[[1]]
[1] 1.462

[[2]]
[1] 4.26 

[[3]]
[1] 5.552

念のため時間を測ってみよう。

system.time(
    for (i in 1:10000)
    {
        for (sp in levels(iris$Species))
            mean(iris[iris$Species == sp, "Petal.Length"])
    }
)
system.time(
    for (i in 1:10000)
    {
        lapply(levels(iris$Species), function(sp) {
            mean(iris[iris$Species == sp, "Petal.Length"])
        })
    }
)

実行速度の差は、なんと誤差範囲内だ。種類が 3 つだけだからか？いや違う。levels(iris$Species) を 50 個つなげても、両者に速度差は出ない。なぜなら、やっていることが同じだからだ。

「for 文を使うのは R 初心者」と言われたりするが、apply 関数が速いと思っているようではまだまだだ。

じゃあ for 文でもいいのだろうか。このような場合はいいかもしれない。だが私は apply を使う。「この場合は for を使ってもいいだろうか？あの場合は？」そんなことは考えたくない。どんな時でも applyだ。私は R が好きだ。

関数型

R は関数型プログラミングをサポートする。ここでモナドとか圏論とか、そういう難しい話をしたいわけではない。

まず、関数の定義は function 関数で行い、変数に代入することで名前をつける。そしてその呼び出しは、 関数名(引数1, ...) とする。特に違和感はないだろう。

add <- function(a, b, c)
{
    a + b + c
}

add(1, 2, 3)

また、関数名をシングルクォートで囲んでもいい。

'add'(1, 2, 3)

演算子

代入演算子は関数だ。シングルクォートで囲むことで、関数として呼び出すことができる。

'<-'(a, 2)
'='(b, 3)

結果として、変数 a に 2、b に 3 が代入される。

比較演算子は関数だ。シングルクォートで囲むことで、関数として呼び出すことができる。

'>='(2, 2) ## 結果は TRUE
'>'(2, 2) ## 結果は FALSE

ベクトルや行列などから添字を指定して値を取り出す [ は関数だ。シングルクォートで囲むことで、関数として呼び出すことができる。

a <- c(11, 12, 13)
a[1] ## 結果は 11
'['(a, 2) ## 結果は 12

これ何に使うの？よく分からないけど面白い。私はRが好きだ。

算術演算子は関数だ。シングルクォートで囲むことで、関数として呼び出すことができる。

'+'(1, 1) ## 結果は 2
'-'(7, 2) ## 結果は 5
'*'(2, 0) ## 結果は 0
'/'(2, 2) ## 結果は 1

ポーランド記法である。なんとも関数型言語のようだ。しかし、伝えたいのは記法ではなく、関数として呼び出せるということだ。その説明のために、まずはリサイクリングを思い出したい。

行列とスカラーの計算は、行列の各要素とスカラーで計算される。次のコードでは要素が全て $0$ の3行2列の行列に $1$ を足している。数学的にも正しい。

matrix(numeric(6), nrow = 3) + 1

\begin{bmatrix} 0 & 0\\ 0 & 0\\ 0 & 0\\ \end{bmatrix} + 1 = \begin{bmatrix} 1 & 1\\ 1 & 1\\ 1 & 1\\ \end{bmatrix}

行列とベクトルの加算、減算は、数学的には定義されない。しかし R においては、ベクトルがリサイクリングにより "行方向に" 繰り返し拡張され、被加算行列と同じサイズの行列として加算される。例えば4行3列の行列（要素は全て 0）とベクトル (1, 2, 3, 4, 5, 6) を加算しようとすると、下記のようになる。

matrix(numeric(12), nrow = 4) + c(1, 2, 3, 4, 5, 6)

\begin{bmatrix} 0 & 0 & 0\\ 0 & 0 & 0\\ 0 & 0 & 0\\ 0 & 0 & 0\\ \end{bmatrix} + Recycling \begin{pmatrix} \begin{bmatrix} 1\\ 2\\ 3\\ 4\\ 5\\ 6\\ \end{bmatrix} \end{pmatrix} = \begin{bmatrix} 0 & 0 & 0\\ 0 & 0 & 0\\ 0 & 0 & 0\\ 0 & 0 & 0\\ \end{bmatrix} + \begin{bmatrix} 1 & 5 & 3\\ 2 & 6 & 4\\ 3 & 1 & 5\\ 4 & 2 & 6\\ \end{bmatrix} = \begin{bmatrix} 1 & 5 & 3\\ 2 & 6 & 4\\ 3 & 1 & 5\\ 4 & 2 & 6\\ \end{bmatrix}

重要なのは、リサイクリングでは "行方向に" 拡張されるということだ。

次の式を考えたい。そう。ニューラルネットワークの層だ。例えば、 $x$ は長さ3のベクトル、 $W$ は3行2列の行列、そして $b$ は長さ3のベクトルとする。

x^TW + b

R では、行列の積は %*% 演算子（もちろん関数）で計算できる（* だと要素どうしの乗算になる）。この時、 $x$ を転置するのは当然なので明示的に行う必要がない。それでも転置するなら t(x) だ。計算結果は、x %*% W、x %*% W + b 、いずれも1行2列の行列となる。そして drop 関数を使って長さ2のベクトルにすることもできる。特に疑問はないだろう。

x %*% W + b

ここからが本題である。バッチ処理を考えたい。 3つの $x$ を取り出し、それぞれを行とした3行3列の行列 $X$ を考える。 $XW$ が3行2列になるので、 $b$ も3行2列に拡張する必要がある。さてどうするか。

\begin{bmatrix} x_{11} & x_{12} & x_{13}\\ x_{21} & x_{22} & x_{23}\\ x_{31} & x_{32} & x_{33}\\ \end{bmatrix} \begin{bmatrix} w_{11} & w_{12}\\ w_{21} & w_{22}\\ w_{31} & w_{32}\\ \end{bmatrix} + \begin{bmatrix} b_1 & b_2\\ b_1 & b_2\\ b_1 & b_2\\ \end{bmatrix}

自動でリサイクリングしてくれると嬉しいのだが。実際に計算して見てみよう。 $b$ の影響だけを見たいので、 $X$ 、 $W$ の要素は全て $0$ とした。

X <- matrix(numeric(9), nrow = 3)
W <- matrix(numeric(6), nrow = 3)
b <- c(1, 2)
X %*% W + b

これを実行すると警告も出ずエラーになることなく、次の結果が表示される。

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

b をリサイクリングしてくれたようだが、期待通りの結果ではない。下記左のように拡張して欲しかったのだが、行方向 に拡張され、右のようになってしまっている。

\begin{bmatrix} b_1 & b_2\\ b_1 & b_2\\ b_1 & b_2\\ \end{bmatrix} , \begin{bmatrix} b_1 & b_2\\ b_2 & b_1\\ b_1 & b_2\\ \end{bmatrix}

ではどうするか。X %*% W のサイズを取得し、行方向に重ねてくれる rbind 関数を行の数だけ繰り返すか。

XW <- X %*% W

B <- b
for (i in 2 : nrow(XW))
    B <- rbind(B, b)

XW + B

行列を作成する matrix 関数において、第一引数のベクトルを行として 列方向 に拡張する byrow オプションを使って拡張するか。

XW <- X %*% W
B <- matrix(b, ncol = length(b), nrow = nrow(XW), byrow = TRUE)
XW + B

いいや、どちらも私の好きな R の書き方ではない。

apply 関数はどうだろうか。第一引数に行列を取り、第二引数でその行列から行(1)、または列（2）どちらを取り出すかを指定し、第三引数で指定した関数に取り出した行（または列）と、第四引数以降で指定した値を渡す。そして第三引数の関数の返り値を 列方向 に結合して行列とする。

次のコードで説明すると、X %*% W の結果の行列から、1行(1)ずつ取り出しそれを o、b を p として無名関数に渡し、加算して返している。

apply(X %*% W, 1, function(o, p) { o + p }, b)

結果を見てみよう。

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

うーん。惜しい。惜しいが、t 関数で転置すればいいだけだ。それよりも function(o, p) { o + p } だ。加算するだけなら '+' 関数が使えるじゃないか。

t(apply(X %*% W, 1, '+', b))

素晴らしい。ようやく辿り着いた。長い旅路の末、算術演算子が関数であることの嬉しさを伝えることができた。ここまで読んでくれた人に感謝を伝えたい。私はあなたが好きだ。いるといいのだが。

実を言うと、転置しているところが気になる。もっと言うと、転置させられているところが気に食わない。この計算において、私は転置したいわけじゃあない。こんなコードは全然クリーンじゃあない。そんな時は sweep 関数だ。第三引数を、第二引数の方向に拡張して第一引数と同じサイズにし、第一引数と合わせて第四引数の関数に渡してくれる。

sweep(X %*% W, 2, b, '+')

私が求めていたのはこの簡潔さ、そして分かりやすさだ。美しい。私はRが好きだ。

環境構築

R をインストールする。以上だ。必要なものは揃っている。なにも include することはない。なにも import しなくとも、複雑な関数が使える。

私にとっての R は、よく手に馴染んだ関数電卓のようなものなのだ。

いや、もう少しだけ。Emacs と、Emacs Speaks Statistics パッケージだけは欲しい。 https://ess.r-project.org

これで Emacs で R のコードを書き、C-c C-c を押して、カーソルのある行を実行できる。Emacs の中に引きこもっていられる。私は Emacs も好きだ。

終了

私は R が好きだ。

私は R が好きだが、精密な機械の制御やエッジコンピューティングに使用しないくらいの分別は持ち合わせているつもりだ。しかし、R 以外の言語で書いていると、「こんな時、 R だったらもっと簡潔に書けるのに」と思ってしまうのだ。

私が愛を語っているだけの文章を、ここまで読んでくれて本当に感謝している。 R をインストールし、インタプリタを起動していくつかの関数を実行してみただろうか。少しでも R を好きになってくれていたら嬉しい。

インタプリタを終了するには quit()、または短く q() と入力して終了できる。そう。関数である。私は R が好きだ。あなたはどうだろうか。