R での処理は一連の関数を用いて行なわれます。 既に多くの関数が定義されており、それだけでも一通りの解析ができるものの、R の真の実力は "自分で関数を書く" ことで発揮されます。 R には、最初から関数を簡便に定義するための多くの便利な機構が組み込まれており、R の組み込み関数の多くも、R で書かれた関数です。関数はエディタで定義ファイルに書き込み source("ファイル名") で読み込むか、R 起動中にインタラクティブに定義します。 いずれにしても、一旦定義されれば、自前の関数も、R 固有の関数も外見からは全く同等に使えます。
> foo <- function(引数リスト) 関数本体 # 関数本体が一つの式なら {} は不要
> foo <- function(引数リスト) {関数本体} # 関数本体が複数(もしくは長い)式なら {} で囲む
> "foo" <- function(引数リスト) 関数本体 # この形式もあり(二項演算子等はこの形式を使う)
警告:R の組み込み関数は自由に再定義可能である。うっかり既存の関数と同じ名前の関数を定義すると、(その関数を内部で使っている他の関数を含み)本来の意味で機能しなくなるのでくれぐれも注意! 二項演算子 + すら再定義可能であるが、もしそうしたら悲惨なことになる。そうした再定義はセーブ付きで終了すれば、.Rdata ファイルに定義が残るので、再実行の際復活する。おかしいと思ったら ls() でオブジェクトの一覧を点検し、rm() 関数で消去する。組み込み関数は rm() 関数で消去されないので、本来の関数定義が復活する。
関数本体は R の式を必要なだけ並べたものです。一つの行に複数の式を並べる時は、セミコロンで区切ります。行末に特別な記号(C における ; など)を置く必要はありません。一つの式が、複数の行に渡っても構いません。
# 例
> foo <- function(x, y) { # foo という名前の関数を二つの引数 x, y で定義
z <- x + y; w <- x * y # セミコロンで区切って一行に
if (z > w) t <- z # 好みで if (z > w) {t <- z} としても可
else {t <- w; w <- z; z <-t} # 波括弧でくくるとブロック化される
for (i in 1:10) { w <- w + i; z <- z - i}
cat("z = ", z, " w = ", w, "?n")
return (list(w, z)) # 二つの返り値
}
> foo(2, 3) # すぐ実行できる
> foo <- function () cat("R is fantastic! ?n") # 引数なし
> foo1 <- function (x) x^2 # 引数一つ(自乗値を計算)
> foo2 <- function (x,y) x*y # 引数二つ(積を計算)
> foo3 <- function (x,y,z) { # 引数3つ(本体が長いので括弧で括る)
w <- exp(x)*(y+z)
cat("x=",x,"y=",y,"z=",z, "result=",w)
}
> foo <- function(x=1, y=2) x+y # 既定値つき(名前つき)引数 > foo() # 既定値 x=1, y=2 と解釈 [1] 3 > foo(0.1, 20) # x=0.1, y=20 と解釈 [1] 20.1 > foo(x=0) # 既定値 y=2 が使われる [1] 2 > foo(y=0) # 既定値 x=2 が使われる [1] 1 > foo(x=5, y=6) [1] 11 > foo(6) # 既定値 y=2 が使われる [1] 8 > foo(6,) # 上と同じ [1] 8 > foo(,6) # 既定値 x=1 が使われる [1] 7
次のような(場合により値が異なる)既定値つきの引数を与えることもできる。
> test <- function (x, y=ifelse(x>=0, 1, -1)) {cat(x,y,"?n")}
> test(1, 2) # y 引数を与えればそれが第二引数の値になる
1 2
> test(1) # y 引数を省略すると、x 引数に応じた既定値が第二引数の値になる
1 1
> test(-1)
-1 -1
digitsBase(x, base = 2, ndigits = 1 + floor(log(max(x), base))) という関数は次の様にさまざまな引数の与え方で実行できる (1) digitsBase(10) base と ndigits は既定値が使われる (2) digitsBase(10, base = 8) base と ndigits は既定値が使われる (3) digitsBase(10, ndigits = 5) base に既定値 2 が使われる (4) digitsBase(10, base = 8, ndigits = 5) 既定値なし
関数 foo <- function (x, mab = 3, mbc = 0.1, bar = 10) {関数本体} は次のような引数の指定で実行できる
foo(3)
foo(3, ma=4, mb=1)
foo(3, b=4)
> test <- function(x,y) print(x+y) > test(1,2) [1] 3 > test(x=1,2) [1] 3 > test(1,y=2) [1] 3 > test(x=1,y=2) [1] 3 > test(y=2, x=1) [1] 3
引数として関数オブジェクトを与えることができる。 普通関数オブジェクト名を引数として与える:
test <- function(x, fun=sin) {plot(x, fun(x))}
test(1:100/100) # 既定の関数 sin が使われる
test(x, fun=cos) # 関数として cos が使われる
test(x, fun=function(z) 1/(1+exp(z))) # 関数定義を直接与える
関数自身がオプション引数をもち、それを特定の値に指定したいときは少し注意がいる。
util <- function(x , opta=1, optb=0) # オプション引数付の関数 test(1:100/100, fun=test(x, opta=10, optb=-1) ) # エラーになる my.util <- function(x) util(x, opta=10, optb=-1) # オプションを固定した関数を定義 test(1:100/100, fun=my.util) # 正しい指定方法 # もしくは引数として関数定義を直接書いても良い test(1:100/100, fun=function(x) util(x, opta=10, optb=-1))
> foo <- function(x, ...) {
if (x >= 0) return(x*cos(...[1]) + x^2*sin(...[2]))
else return(exp(x))}
> foo(1,c(2,3)) # ... にベクトル値
[1] -0.2750268
> foo(-1) # ... 部分の引数無し。実際この場合必要無い
[1] 0.3678794
> test1 <- function(x, y=3, ...) (x+...[1])^y > test2 <- function(x, ..., y=3) (x+...[1])^y # ... 引数の後に y 引数 > test1(1); test2(1) [1] 8 # = (1+1)^3 [1] 8 # = (1+1)^3 # ...[1] = 1 とされている (?) > test1(1,2); test2(1,2) [1] 4 # (1+1)^2 [1] 27 # (1+2)^3 # ...[1] = 2 とされている (問題無し) > test1(1,2,3) [1] 16 # (1+3)^3 > test2(1,2,3) # ... 引数の後の引数は名前付きで与える必要 Error in ...[1] : incorrect number of dimensions > test2(1,2,y=3) [1] 27 # (1+2)^3 # y 引数を名前付きで与えた (問題無し)
参考: ... 引数は関数引数欄に明示的に与えなくても不思議な既定値を持つ(ほとんどバグ?) > test3 <- function(x, y=3, ...) print(...[1]) > test3(1) [1] 1 # ...[1] = 1 とされている > test3(3,2) [1] 1 # ...[1] = 1 とされている > test3(3,2,1) [1] 1 # ...[1] = 1 とされている > test4 <- function(x, y=3, ...) print(c(...[1],...[2],...[3])) > test4(3) [1] 1 2 3 # ... = c(1,2,3) とされている
その他大勢引数が実際は複数の引数になるとき、k 番目の項目を取り出すには list(...)[[ ]] のようにすれば良い(2006/6/22 補足)。
> test <- function(i,...) print(list(...)[[i]])
> test(1,1,"abc",list(1:4),matrix(1:4,2,2))
[1] 1
> test(2,1,"abc",list(1:4),matrix(1:4,2,2))
[1] "abc"
> test(3,1,"abc",list(1:4),matrix(1:4,2,2))
[[1]]
[1] 1 2 3 4
> test(4,1,"abc",list(1:4),matrix(1:4,2,2))
[,1] [,2]
[1,] 1 3
[2,] 2 4
従って ... 引数が実際は幾つの引数かを知るには length(list(...)) とすれば良い。
関数は返り値を持っても良いし、持たなくてもよい。返り値のない関数は、その本体で実行される実行文の 副作用 のために使われる
もし return(), invisible() 関数により返り値を陽に指定しなければ,(関数本体の最後に実行される式が値を返す限り)それが関数の暗黙の返り値とされる
> foo1 <- function (x) x^2 # x^2 の値が暗黙の返り値
> foo <- function (x,y) {return(x+y)}
> foo(1,2)
[1] 3 # x+y が返り値
> foo <- function (x, y) {return(x); cat("sum =", x+y)}
> foo(1,2) # cat は実行されない
[1] 1
> foo <- function (x,y) {if (x>=y) return(x) else return(y)} # min(x,y) を計算
> foo(1,2)
[1] 2
> foo(2,1)
[1] 2
> foo <- function (x,y) {return(x,y,x+y)}
> res <- foo(1,2)
> res[[1]]
[1] 1
> res[[2]]
[1] 2
> res[[3]]
[1] 3
> foo <- function (x,y) {return(x=x,y=y,sum=x+y)} # ラベル付き返り値
> res <- foo(1,2)
> res$x # res[[1]], res[["x"]] でも良い
[1] 1
> res$y # res[[2]], res[["y"]] でも良い
[1] 2
> res$sum # res[[3]], res[["sum"]] でも良い
[1] 3
例えば、非常に大量で、可読性の低い値を返す関数を考えると、返り値を変数に付値する場合以外コンソールに出力したくない事が起きる。invisible() はコンソールに表示されない返り値を返す。
> foo <- functionn (x) return(x) > foo(1) [1] 1 > foo <- functionn (x) invisible(x) > foo(1) # 何も表示されない > x <- foo(1) > x # 結果は変数 x に付値されている [1] 1
注意:invisible() 関数は関数の返り値の指定に使われるだけでなく、次のような使い方もできる。
> test <- function(x) {print(x); x <- x+1; x}
> test(3)
[1] 3 # print(x) による出力
[1] 4 # 関数 test の返り値
> temp <- invisible(test(3))
[1] 3
> temp
[1] 4
普通最後の結果だけが返り値になるという特徴は、途中結果を返したい時は困ることがある。 解決策は必要な途中結果を print(), cat() 関数でコンソールに出力させるか、途中結果をベクトル、リストに蓄積しておき、最終的な返り値にすることである。apply 関数群は後者を暗黙の内に実行するので便利である。
# 途中結果を print() 関数でコンソールに出力する
> test <- function(x) {for (i in 1:3) print(x+i)}
> test(3)
[1] 4
[1] 5
[1] 6
> temp <- test(3)
[1] 4
[1] 5
[1] 6
> temp # しかし返り値は依然として最後のものだけ
[1] 6
途中結果をベクトルに貯めておき、最後に返り値とするやりかた
> test2 <- function(x) {res <- numeric(0); for (i in 1:3) res <- c(res, x+i)}
> test2(3)
> test2 <- function(x) {res <- numeric(0); for (i in 1:3) res <- c(res, x+i); res}
> test2(3)
[1] 4 5 6 # 途中結果のベクトル返り値
結果として返す値の個数(ベクトルの要素数)がわかっている場合は,numeric(0) のような,メモリー確保ではなく,正確な個数のメモリーを確保するのがよいのはいうまでもない。
> test9 <- function(n) {res <- numeric(n); for(i in 1:n) res[i] <- factorial(i); return(res)}
> test9(10)
関数中で永続付値演算子を使うと、関数終了後も残る値を変数に付値できる。関数実行前には存在しなかった変数を作ることもできる。永続付値演算子 <<-, ->> は現在の実行環境の親環境(しばしば R の大局的環境(.GlobalEnv という名前の隠しオブジェクト)中のオブジェクトに対する操作*1であり、特に関数中で用いるときは注意がいる。一方で、親環境のオブジェクトへの付値はその関数が終了後も残り、特に大局的環境中のオブジェクトは任意の関数から直ちに参照でき、特に副関数に複雑な引数を引き渡すとき便利である。永続付値演算子を関数中で用いるとき陥りやす間違いは次の例でみられる:
> v <- 1
> test <- function () {v <- 2; v <<- 3; cat("v=", v, "\n")}
> test()
v= 2
> v # いつのまにか値が変わっている
[1] 3
この例では関数内部での変数 v と永続付値で生成された大局的環境(つまり関数 test を実行した環境 .GlobalEnv )中の変数 v という、同じ名前だが実体が異なる二つの変数が登場する。関数中では単なる変数 v は局所的な変数を指し、関数実行後は残らない。一方大局的変数 v は関数実行後も残る。
もう少し微妙な例:
> test <- function () {v <- 1:3; v[2] <<- 0; cat("v=", v, "?n")}
> test()
v= 1 2 3 # 局所的な v の値
> v
[1] 1 0 3 # 大局的な変数 v がつくり出され、いつの間にかベクトル値を持っている)
大局的変数への永続付値の便利な使い方は、他の関数で使うべき値を暗黙のうちに引き渡すことである(次の例を参照)
> foo <- function (x, y, z) {.w <<- list(x, y, z)} # 大局的変数 .w に付値
> foo2 <- function () {.w[[1]] + .w[[2]] + .w[[3]]} # 大局的変数 .w を暗黙のうちに参照
> foo(1, 2, 3) # これで大局的変数 .w (=list(1 , 2, 3)) が作られる
> foo2()
[1] 6
教訓: 関数中の永続付値は便利だが使用には注意がいる。さもないと、理解不能なバグに悩まされ、睡眠不足に陥る可能性がある。使うときはできるだけ個性的な名前を用いる。決して x, y 等のありふれた変数名、もしくは関数中の局所的変数と紛らわしい変数名は使わない!
R のすべての実行式は対応する環境を持つ。普通重要になる環境は大局的環境(R 起動時の環境).GlobalEnv と関数呼び出しの際の関数内から見た親環境(その関数を呼び出した環境)parent.frame() である。そうした環境に変数名とその値を登録すれば、結果として関数中から値を返せることになる。以下の例は関数中で作られたオブジェクト mat を他の環境に mat という名前で登録する方法を示す:
# 関数の実行環境の親環境(その関数を呼び出した環境)に登録
# mat <<- mat と同じことになる
assign("mat", mat, parent.frame())
assign("mat", mat, .GlobalEnv) # 大局的環境に登録
stop() 関数は、現在の式の実行を中断し、その引数として与えられたメッセージを表示*2し、それから options(error=...) で指定された、エラー時動作を実行する。既定では トップレベルのプロンプトに戻る
# x が負またはゼロなら中断、正なら log(x) を返す関数 > foo <- function(x) ifelse( x<=0, stop(message="non-positive arg!"), log(x)) > foo(3) [1] 1.098612 > foo(-3) Error in rep(yes, length = length(ans)) : non-positive arg! # エラー
警告メッセージを生成する。 実際の結果は options("warn") の値に依存する
> foo <- function(x) {
if (x<=0) {return(NA); warning(message="NA produced")}
else return(log(x)
}
> foo(3)
[1] 1.098612
> foo(-3)
Error in rep(yes, length = length(ans)) : non-positive arg! # エラーメッセージを表示し中断
引数としてカンマで区切って並べられた条件が満たされないときにプログラムをストップする。
(stopif という関数仕様の方が素直で便利だと思うのだが)
> a <-3 > b <- 4 : > stopifnot(a == 3, b==5) Error: b == 5 is not TRUE
デバッグは困難な問題であるが、最大の困難ではない。デバッグ中は少なくとも関数にエラーがあることだけは確信できるが、デバッグにより一見エラーがなくなった後に、真の困難が登場する。つまり この関数は本当に求めるものを計算しているのだろうか。
最も安易な関数のデバッグ法は、関数中に cat() を沢山埋め込み、途中の変数値を出力させることである。
> foo <- function (x) {
cat("arg x=", x, "?n")
log(x) }
> foo(-1)
arg x= -1 # デバッグメッセージ (なるほど log 関数は正の値しか受け付けない!)
[1] NaN
Warning message:
NaNs produced in: log(x)
デバッグ用に cat() を沢山埋め込んだ後は、すぐに消さないで、それらをコメントアウトして残しておくのが賢明。複数の行を一度にコメントアウトする方法は以下を参照。
> foo <- function (x) {
#cat("arg x=", x, "?n") # 一行コメントアウト
log(x) }
> foo <- function (x, y) {
if (0) { cat("1st arg =", x, "?n")
cat("2nd arg =", y, "?n") } # 二行まとめてコメントアウト
log(x)*log(y) }
注意:if (0) { と } の中には,R の文法上正しいものしか書けない(syntax check がかかる)
> bar <- function(x) {
+ if (0) {
+ この行はエラーにならないけど,次の行はエラーになる(アッタリ・メーダ)
> x <- b±root(b2-4ac)
> }
エラー: 予想外の '}' です ( " }" の)
> return(x+1)
エラー: オブジェクト 'x' がありません
> }
エラー: 予想外の '}' です ( "}" の)
# browser 関数埋め込み例
> foo <- function(x, y) {z <- x*y; browser(); w <- x + y; return(z, w)}
> a <- foo(2,3) # 関数実行
Called from: foo(2, 3) # 中断メッセージ
Browse[1]> ls() # 現在の環境内のオブジェクトを一覧する
[1] "x" "y" "z"
Browse[1]> x # 変数 x の値は?
[1] 2
Browse[1]> y
[1] 3
Browse[1]> z
[1] 6
Browse[1]> w # 変数 w はこの時点では存在しない
Error: Object "w" not found
Browse[1]> c # 関数実行継続
> a$z
[1] 6
> a$w
[1] 5
> foo <- function(x, y) {z <- x*y; w <- x + y; return(z, w)}
> debug(foo) # デバッグ用フラグを立てる
> foo(2,3) # 関数実行
debugging in: foo(2, 3)
debug: {
z <- x * y
w <- x + y
return(z, w)
}
Browse[1]> ls() # 現在の関数実行環境中のオブジェクト(引数変数だけ)
[1] "x" "y"
Browse[1]> x
[1] 2
Browse[1]> n
debug: z <- x * y
Browse[1]> ls()
[1] "x" "y"
Browse[1]> n
debug: w <- x + y
Browse[1]> ls()
[1] "x" "y" "z" # 変数 z が作られている
Browse[1]> z
[1] 6
Browse[1]> n
debug: return(z, w)
Browse[1]> ls()
[1] "w" "x" "y" "z" # 変数 w も作られている
Browse[1]> w # 変数 w の値
[1] 5
Browse[1]> n
exiting from: foo(2, 3) # 関数実行終了
名前の無い関数があり得る。例えば関数中である特別な値を一時的に作りたい時などに使う
> (function(z){dnorm(z)})((0:8)/2) # N(0,1) の密度関数の値
[1] 0.3989422804 0.3520653268 0.2419707245 0.1295175957 0.0539909665
[6] 0.0175283005 0.0044318484 0.0008726827 0.0001338302
> foo <- function (x) {y <- (function(z){dnorm(z)})((0:8)/2); x*y} # N(0, x^2) の密度関数
> foo(2)
[1] 0.7978845608 0.7041306535 0.4839414490 0.2590351913 0.1079819330
[6] 0.0350566010 0.0088636968 0.0017453654 0.0002676605
> sapply(airquality, mean) # NA を含むデータファイルの場合には,平均値が計算されない
Ozone Solar.R Wind Temp Month Day
NA NA 9.957516 77.882353 6.993464 15.803922 # Ozone, Solar.R には欠損値 NA があるので,平均値が NA になる
> sapply(airquality, function(x) mean(x, na.rm=TRUE)) # 無名関数を使えばよい
Ozone Solar.R Wind Temp Month Day
42.129310 185.931507 9.957516 77.882353 6.993464 15.803922
> sapply(airquality, mean, na.rm=TRUE) # この場合には,以下のようにすればよいのだけど
Ozone Solar.R Wind Temp Month Day
42.129310 185.931507 9.957516 77.882353 6.993464 15.803922
関数の定義中に、その関数自身が登場しても良い。当然、登場した関数の値が合理的に定義可能である必要がある。
> foo <- function(x) {ifelse(x==1, 1, x*foo(x-1))} # 即席階乗関数
> foo(1) # 1!
[1] 1
> foo(2)
[1] 2 # 2! (foo(2) を計算する際、foo(1) は既に定義済み)
> foo(10)
[1] 3628800 # 10!
R の関数内部には別の関数定義を置き、その場で実行できる。但し、内部関数は単独では存在しない。(単独では存在しない,というのではなく,foo の外では,foo2 は存在しないということ。つまり,foo の外で foo2() とやっても,「そんな関数はない」といわれるということ。)
> foo <- function() {
foo2 <- function() { # 関数内部での関数定義
x <- 1:3; y <- 23:34; z <- cbind(x, y); return(x, y, z)}
c <- foo2() # その場で実行
return(c$z, c$x, c$y) }
> x <- foo()
> x[[2]]
[1] 1 2 3
> x[[3]]
[1] 23 24 25 26 27 28 29 30 31 32 33 34
> foo2() # 関数 foo2 は foo 終了後は存在しない
Error: couldn't find function "foo2"
> foo <- function (x) {
y <- x
temp <- function (a) a*y # temp 定義中の変数 y は直前の行で定義された変数
temp(2) }
> foo(4)
[1] 8
こんなのもあり
> foo <- function(x) {
+ foo2 <- function(x) {foo3 <- function(x) x+1; foo3(x)}
+ foo2(x)}
> foo(3) # foo は foo2 を呼出し、foo2 は foo3 を呼び出す
[1] 4
R における +, -, *, / 等の組み込み二項演算子も特殊な形の関数である。R には独自に二項演算子を定義する機構がある。詳しくは 二項演算子定義の例 を参照。
> "%+=%" <- function(x, y) {x <<- x + y} # C 風のインクリメント演算子
> x <- 10; x %+=% 5
> x
[1] 15 # x は 5 だけ増えている
注意:C だと,以下の演算と同じことをする x+= 5*7 の後,xは 45 になる,Rでは表示されるものと,xの内容が異なるし,Cと同じ結果にもならない
> x <- 10 > x %+=% 5*7 [1] 105 > x [1] 15
R の関数は R の任意のオブジェクトを返すことができる。当然関数オブジェクトを返すことも可能。次の例を参照。関数生成関数 foo への引数が、返り値関数の内部パラメータになっていることを注意(つまり返り値である関数 Foo の附属環境に登録されている)。
> foo <- function(a, b) {function (x, y) a*x + b*y}
> Foo <- foo(2, 3) # Foo(x, y) は関数 2*x + 3*y
> Foo(1, 2)
[1] 8
> Foo # 定義を見る (a=1, b=2 の情報は附属環境に含まれており表に見えない)
function(x,y) a*x+b*y
<environment: 0x846e354>
> Foo <- foo(-2, -3) # Foo(x, y) は関数 (-2)*x + (-3)*y
> Foo(1, 2)
[1] -8
R の スコープ規則 (ある名前のオブジェクトの居場所を決定する規則) について簡単に知っておくことは、関数の振舞いを理解するのに役立つ。
> x <- 3 # 関数 foo の実行環境中の変数 x
> foo <- function (x) { # 関数 foo はその実行環境中のオブジェクト、引数 x は附属環境中のオブジェクト(?)
+ y <- 5 # 関数 foo の附属環境中の変数 y
# 関数 foo2 は foo の附属環境中のオブジェクト、変数 z は関数 foo2 の附属環境中の変数
# foo2 を計算する際、x と y (とその値)はその親環境のものとされる
+ foo2 <- function(z) {x+y+z} + foo2(4)}
> foo(1) # x+y+z=1+5+4 を計算
[1] 10
on.exit(EXPR) は関数が正常・異常終了した際に EXPR を実行する。関数中で変更したR の基本的設定を元に戻すのに便利。on.exit(EXPR, add=TRUE) は先の on.exit 関数の実行内容に付け加える。
> foo <- function() {
old.par <- par(no.readonly = TRUE) # 現在の作図パラメータ記録
on.exit(par(oldpar) # 関数終了時に元に戻す
on.exit(dev.off(), add=TRUE) # 更に、関数終了時にデバイスを閉じる
png("foo.png") # png デバイスを開く
--- 幾つかの作図命令 ---
}
R 言語で頻繁に使われる plot() 等の関数は総称的(generic)な関数と呼ばれる。 これらの関数は引数である R オブジェクト x に応じて、実際に使われる関数が決まる。例えば x が線形回帰関数 lm() の返り値であれば、x はクラス属性 "lm" を持ち、plot() はこれを確認すると、線形回帰オブジェクトのプロット用に設計されたプロット関数 plot.lm() を起動する。関数 plot.lm() はプロット関数 plot() の一つのメソッド(method)と呼ぶ。もちろん、クラス "lm" であることがあらかじめわかるなら plot.lm() 関数を直接使っても良い。plot() 関数はその他にも既定で多くのメソッドを持つ。また、ユーザーが独自のクラスに対する独自のメソッドを定義することも可能である。この機構は、多くの同種の関数を同じ関数で代表させることで、コードの大幅な簡略化を可能にする。
# 例:総称的関数とメソッド > methods(plot) # 総称的関数 plot() に対するメソッドの一覧を得る [1] plot.HoltWinters* plot.POSIXct plot.POSIXlt [4] plot.TukeyHSD plot.acf* plot.data.frame (途中略) [28] plot.table plot.ts plot.tskernel* Non-visible functions are asterisked
> methods(class = "lm") # クラス "lm" のオブジェクトに適用可能なメソッドの一覧 [1] add1.lm alias.lm anova.lm case.names.lm [5] coef.lm confint.lm cooks.distance.lm deviance.lm (途中略) [33] variable.names.lm vcov.lm
例えば Unix と MS Windows では仕様が異なる関数を統一的に定義したいとする。
> .Platform$OS.type # 使用OS情報
[1] "unix"
> foo <- if (.Platform$OS.type = "unix" ) {function() cat("using *unix.\n"}
else {function() cat("Maybe using MS Windows?\n")}
> system.time(for (i in 1:10000) mean(runif(500))) [1] 2.00 0.06 4.15 0.00 0.00
時間関数の返すベクトルは順に (最後の二つは普通ゼロ)
Rは所詮インタプリター言語であり、コンパイラー型言語に比べれば速度が劣るのは致し方ない。それでも、R言語の特性をうまく使えば信じられないほどの速度を達成することも夢では無い。いくつかのヒントをあげる。
論理判断は如何なる言語でも相対的に時間を食う。避けられないことも多いが、工夫で避けられることもある。R の組み込み関数はできるだけ速度が早くなるように工夫(しばしば C や Fortran サブルーティンを使用)されているので、それを使うことが好ましい。
> foo <- function(x) {if(x<=0) return(-x) else return(x)} # 手製の絶対値関数
> system.time(for (i in 1:100000) {x=runif(1)-0.5; foo(x)} )
[1] 2.44 0.06 2.80 0.00 0.00
> system.time(for (i in 1:100000) {x=runif(1)-0.5; abs(x)} ) # 餅は餅屋の例え
[1] 1.90 0.01 2.16 0.00 0.00
> system.time(for (i in 1:100000) {x=runif(1)-0.5; sqrt(x^2)} ) # これでも少しはまし
[1] 2.31 0.04 2.55 0.00 0.00
> system.time(for (i in 1:100000) {x=runif(1)-0.5; max(x,-x)} )
[1] 2.29 0.08 2.61 0.00 0.00
> foo <- function (x) {ifelse (x==0, 1, 0)} # 引数 x は 0 か 1 とする
> x <- floor(2*runif(100000)) # 長さ 10000 の 0, 1 ベクトル
> system.time(foo(x))
[1] 0.54 0.02 0.59 0.00 0.00
> foo2 <- function (x) 1-x # 実は foo と同値
> system.time(foo2(x))
[1] 0.00 0.01 0.01 0.00 0.00
R はベクトル・配列演算を内部 C ルーティンで処理し、高速になるように最適化されている。for ループはしばしばベクトル・配列演算で代用できる。しかもコードが短くなり、可読性が高まると言うおまけもつく。
> foo <- function (X, Y) { # 二つの 100x100 行列の愚直な掛け算
Z <- matrix(0,ncol=100, nrow=100)
for (i in 1:100)
for (j in 1:100)
for (k in 1:100) Z[i,j] <- Z[i,j] + X[i,k]*Y[k,j]
Z}
> X <- matrix(runif(10000), ncol=100) # 100x100 行列を定義
> Y <- matrix(runif(10000), ncol=100) # 同上
> system.time(Z<-foo(X, Y))
[1] 23.05 0.00 23.54 0.00 0.00
> system.time(Z <- X%*%Y) # 専用演算子を使用
[1] 0 0 0 0 0 # 一瞬!
上の例に出てくる 0/1 ベクトルを作るのには,最適な関数がある。sample 関数である。
> system.time(x <- floor(2*runif(10000000)))
ユーザ システム 経過
0.804 0.216 1.015
> system.time(y <- sample(0:1, 10000000, replace=TRUE))
ユーザ システム 経過
0.432 0.102 0.531
関数の実行速度をあげるためには、どの部分で時間を喰っているか確認する必要がある。普通、関数の実行時間の大部分は極少数のコードの実行に当てられていることが多い。その部分を改良すれば、大幅な速度の向上が期待できる。Rprof() *3はコード中の各実行単位の使用回数を、一定のサンプリング期間毎にチェックし、レポートファイル(既定では Rporf.out)を作成する。結果はレポートファイルを直接みても良いし、命令 summaryRprof() で簡潔に一覧することもできる。
> Rprof() # 既定のログファイル Rprof.out とサンプリング間隔 0.02秒 でプロファイル開始
> x <- y <- matrix(runif(10^6),ncol=1000) # 大きな行列
> z <- x%*%y # その掛け算
> Rprof(NULL) # プロファイル終了
> summaryRprof() # 既定のログファイル Rprof.out からプロファイル結果要約
$by.self
self.time self.pct total.time total.pct
"%*%" 2.26 71.5 2.26 71.5 # 行列掛け算がやはり一番時間を喰う
"runif" 0.42 13.3 0.42 13.3 # 一様乱数生成がその次
"as.vector" 0.24 7.6 0.66 20.9 # (恐らく) ベクトルを行列化するオーバーヘッド
"matrix" 0.24 7.6 0.90 28.5 # (恐らく) ベクトルを行列化するオーバーヘッド
$by.total
total.time total.pct self.time self.pct
"%*%" 2.26 71.5 2.26 71.5
"matrix" 0.90 28.5 0.24 7.6
"as.vector" 0.66 20.9 0.24 7.6
"runif" 0.42 13.3 0.42 13.3
$sampling.time
[1] 3.16
使い勝手の向上はしばしば速度と反比例するが、より好ましいことが多い
R の多くの組み込み関数は意味がある限りできるだけベクトル化されている。引数にベクトルを与えることができ、対応する返り値もベクトルになる。こうすることにより、むき出しの for ループを隠し、コードの短縮化ができ、可読性も高まる。
> x=1:10 # 等差数列
> x^2
[1] 1 4 9 16 25 36 49 64 81 100 # 巾乗関数もベクトル化されている
> x <- matrix(1:16, ncol=4)
> x^2 # それどころか行列化されている!
[,1] [,2] [,3] [,4] # 実は R では配列は特殊なベクトルとして表現されているので意外ではない
[1,] 1 25 81 169
[2,] 4 36 100 196
[3,] 9 49 121 225
[4,] 16 64 144 256
> choose(10, 0:10) # 二項係数もベクトル化されている!
[1] 1 10 45 120 210 252 210 120 45 10 1
> choose(10:13, 10)
[1] 1 11 66 286
> choose(10:13, 0:10) # はてこれは何を計算してくれたのだろう?
[1] 1 11 66 286 210 462 924 1716 45 55 66
> runif(10)>0.5 # 論理判断もベクトル化されている!
[1] TRUE TRUE TRUE FALSE FALSE FALSE TRUE FALSE FALSE TRUE
ベクトル化された式だけを使えば結果は自然にベクトル化される。
> foo <- function (x) sin(x>0.5) # R では論理値 TRUE, FALSE は数値 1, 0 とも解釈される > foo(runif(10)) [1] 0.000000 0.841471 0.000000 0.000000 0.841471 0.000000 0.841471 0.841471 [9] 0.841471 0.841471
無理矢理ベクトル化する方法の例
> foo <- function(x) if (x<0) return(2) else return(3) # ベクトル化されていない関数の例
> foo(runif(10) - 0.5) # x の最初の要素だけ受け付ける
[1] 2
Warning message:
the condition has length > 1 and only the first element will be used in:
if (x < 0) return(2) else return(3)
> foo <- function (x) { # そのベクトル化版
temp <- function(x) if (x<0) return(2) else return(3)
sapply(x, temp) } # sapply 関数は x の各要素に temp を施した結果のベクトルを返す関数
> foo(runif(10)-0.5)
[1] 2 2 2 2 3 3 3 3 2 2
> ifelse(runif(10) - 0.5 < 0, 2, 3) # 実はベクトル化された関数 ifelse を使えばすんだこと
[1] 3 2 2 3 3 3 3 3 2 2
> (2:3)[(runif(10) < 0.5)+1] # トリッキーな方法を追加
[1] 2 3 2 3 3 3 2 2 3 3
引数には自然な制約条件があることもある。そうした場合、範囲外の引数が与えられたらエラーにするか、警告を出すようにしておくと、デバッグがし易い。
演算途中の変数値にも自然な制約範囲があることがある。もし制約範囲外なら、残りの演算がエラーになるか、無意味な結果を与えそうなら、エラーもしくは警告を出すようにすると、デバッグがし易くなる。
幾つかの変数には容易に意味が分かる名前を付け、合理的な既定値を与えておく
他人が使う可能性のあるコードには積極的にコメントを要所要所に入れておくと、意味が掴み易い。自分しか使わない場合も、後から何をしているのか分からなくなることはよくあること。
必要な機能を細分、パーツ化し小関数化する(Unix 精神)。 短い関数はチェックがし易い。また既に既成のものが存在する可能性も高い。部分部分を関数化することにより全体の流れが見易くなる。
できるだけ短くしたエレガントなコード、速度を最適化したアクロバティックなコードは、後から読むと、何が何だか自分でも理解不能になりがち。
苦労して書いた関数は、しばしばスパゲッティコードになり、醜くくなる。最初から新たに書くとこつも分かっており、よりきれいなコードが書けることが多い。