R バージョン3に寄せて
をテンプレートにして作成
[
トップ
] [
新規
|
一覧
|
検索
|
最終更新
|
ヘルプ
]
開始行:
東京工業大学の[[間瀬茂]]です．拙著「Rプログラミングマニュ...
拙著の基本方針は「Rのプログラミング機能を網羅的に紹介する...
最後の仕事になります)．出版時期はまだ未定です．
改訂作業中に特に「!!!」と思った点が幾つかありました．これ...
もし，私の誤解，理解不足がありましたら，コメントいただけ...
注意：頁タイトルが誤解を与えるかもしれませんが，ここに書...
(独り言：Rバージョン4が出る頃，私はまだ生きているかな......
----
#contents
----
*「apply関数族は簡潔な既述ができるが，実体は for ループを...
lapply, vapply, sapply, rapply 関数は内部関数を用いて計算...
以下は sapply 関数と for ループを使った同値な処理をそれぞ...
foo <- function(x) x # つまらない関数
X <- numeric(1e6) # 念のため，予め変数を用意
# sapply関数を用いた処理を100回繰り返す
Y <- sapply(1:100, function(i) system.time(X <- sapply(1...
# そのユーザ時間の要約
> c(fivenum(Y[1,]),mean(Y[1,]),sd(Y[1,]))
[1] 2.048 2.240 2.280 2.388 2.604 2.303 0.09771
# forループを用いた処理を100回繰り返す
Z <- sapply(1:100, function(x) system.time(for(i in 1:1e...
# そのユーザ時間要約
> c(fivenum(Z[1,]),mean(Z[1,]),sd(Z[1,]))
[1] 2.548 2.592 2.606 2.624 2.684 2.609 0.02449
- OS やメモリーなどの環境によるのかもしれませんが，以下の...
示されたプログラムでは system.time(X <- sapply(1:1e6,foo)...
foo <- function(x) x # つまらない関数
sim <- function(n, m) {
X <- numeric(m) # 念のため，予め変数を用意
Y <- replicate(n, {gc(); gc(); system.time(X <- sapply(...
Z <- replicate(n, {gc(); gc(); system.time(for(i in 1:m...
return(list(s=Y, f=Z))
}
得られた測定データの最小値，最大値，平均値，標準偏差を以...
n=1e4 n=1e6
sapply for sapply for
min 0.00700 0.01600 2.11900 1.73000
max 0.01900 0.02300 2.91500 2.76100
mean 0.00868 < 0.01838 2.45032 > 1.84305
sd 0.00140 0.00183 0.12012 0.22043
n=1e4 のとき（左図）は for は sapply の倍くらい時間がかか...
&ref(fig1.png); &ref(fig2.png);
- 早速のコメントありがとうございます．実行環境により様子...
n=1e5 sapply
0.132 0.140 0.144 0.148 0.180 0.145 0.00770
n=1e5 for ループ
0.248 0.256 0.260 0.264 0.292 0.260 0.00780
n=1e4 sapply
0.0120 0.0120 0.0120 0.0140 0.0160 0.0130 0.00174
n=1e4 forループ
0.0240 0.0240 0.0240 0.0280 0.0320 0.0253 0.00205
- あ, sapplyはlapplyのただのラッパーです. 性能のバラツキ...
#comment
*「Rシステム」はあるが「R言語」は無い． [#x5862e1a]
現在「R言語」という呼び方がしばしばみられますが，私は二重...
例えば，公式マニュアル An Introduction to R には
”a well developed, simple and effective programming lang...
which includes conditionals, loops, user defined recursi...
and output facilities. (Indeed most of the system suppli...
themselves written in the S language.)”
公式FAQ集の記述では
"We can regard S as a language with three current implem...
or“engines”, the “old S engine” (S version 3; S-PLUS 3.x...
the “new S engine”(S version 4; S-PLUS 5.x and above), a...
Given this understanding, asking for “the differences be...
really amounts to asking for the specifics of the R impl...
the S language, i.e., the difference between the R and S...
#comment
* Tierny氏の compiler パッケージが推奨パッケージになった ...
R開発メンバーのTierny氏の compiler パッケージのおかげで，...
関数 enableJIT は実行時(Just-In-Time)コンパイルを指示しま...
関数 compilePKGS はパッケージを読み込む際に，含まれる関数...
# lapply関数の初期の定義
> la1 <- function(X,FUN,...) {
FUN <- match.fun(FUN)
if (!is.list(X))
X <- as.list(X)
rval <- vector("list",length(X))
for(i in seq(along=X))
rval[i] <- list(FUN(X[[i]],...))
names(rval) <- names(X)
return(rval) }
# 時間のかかるまずいコード
> tmp <- function(n) {
x <- numeric(0)
for(i in 1:n)
x <- c(x,runif(1))
return(x) }
# コンパイルする
> la1c <- cmpfun(la1)
> tmpc <- cmpfun(tmp)
> la1c
function(X,FUN, ...) {
FUN <- match.fun(FUN)
if (!is.list(X)) X <- as.list(X)
rval <- vector("list",length(X))
for(i in seq(along=X)) rval[i] <- list(FUN(X[[i]],...))
names(rval) <- names(X)
return(rval) }
<bytecode: 0xa9a3388> # バイトコンパイルされている証拠
> tmpc
function(n) {x <- numeric(0); for(i in 1:n) x <- c(x,run...
<bytecode: 0xa64e254>
# 以下は la1(1:100,tmp) を1万回実行したユーザ時間の要約
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.016 0.020 0.020 0.021 0.020 0.024 # JITコンパ...
0.004 0.012 0.012 0.011 0.012 0.024 # enableJIT(1)
0.004 0.012 0.012 0.011 0.012 0.016 # enableJIT(2)
0.008 0.012 0.012 0.011 0.012 0.016 # enableJIT(3)
#comment
* バイト型データと関連関数の導入 [#o50edeef]
Rの原子オブジェクトに新しく(何時から？)バイト(raw)型デー...
# 255は16進法ではff
> x <- as.raw(255); x
[1] ff
> str(x)
raw ff
> x <- "A test string"
# xのバイト列表現
> y <- charToRaw(x); y
[1] 41 20 74 65 73 74 20 73 74 72 69 6e 67
# 文字列に戻す
> rawToChar(y)
[1] "A test string"
# 整数8のビット列表現
> intToBits(8L)
[1] 00 00 00 01 00 00 00 00 00 00 00 00 00 00 00 00 00 ...
[26] 00 00 00 00 00 00 00
#comment
* 隠し変数 .Last.value [#oa840bee]
Rのトップレベル(つまりRコンソールのプロンプト > に対して...
続けて二回実行することが勧められています．
> x <- log(10); x
[1] 2.302585
> .Last.value
[1] 2.302585
# library関数の返り値は不可視
> library(splines)
# .Last.valueには記録
> .Last.value
[1] "splines" "stats"
--以下省略--
# 巨大なオブジェクトx
> x <- rnorm(1e6)
# .Last.valueにも記録
> identical(x,.Last.value)
[1] TRUE
# この時点で.Last.valueはすでに更新されている
> str(.Last.value)
TRUE
ユーザレベルでも .Last.value 変数の便利な使い方があります...
> sin(1)
[1] 0.841471
> res <- .Last.value; res
[1] 0.841471
#comment
* マルチコア並列処理パッケージ parallel が推奨パッケージ...
Rの推奨パッケージに新しく付け加わった並列処理パッケージ p...
# クラスタを作り並列処理する例
> require(parallel)
# 4コアのクラスタを生成
> cl <- makeCluster(4)
# 変数xxと関数fooを定義
> xx <- 1; foo <- function(y) xx + y
# 小プロセスにそれを移出，オブジェクト名は文字列で与える
> clusterExport(cl,c("xx","foo"))
# sapply関数の並列版parSapplyを用いfoo(1:10)をクラスタに...
> parSapply(cl,1:10,foo)
[1] 2 3 4 5 6 7 8 9 10 11
# 2変数関数を定義し，小プロセスにそれを移出
> bar <- function(x,y) x + y
> clusterExport(cl,"bar")
# 多変数apply関数であるmapply関数の並列版clusterMap．以...
> clusterMap(cl,bar,c(a=1,b=2,c=3),
c(10,0,-10))
$a # つまりbar(a=1,10)
[1] 11
$b
[1] 2
$c
[1] -7
> stopCluster(cl) # 最後にクラスタを閉じる
関数 pvec はクラスタの生成，apply 関数風の並列計算と結果...
次の結果は日付文字列一万個を as.POSIXct 関数で日付オブジ...
# 日付文字列一万個のベクトルdatesをas.POSIXctで日付オブ...
# すべて経過時間を示す
> dates <- sprintf('%04d-%02d-%02d',as.integer(2000+rnor...
as.integer(runif(1e4,1,12)),as.intege...
> N <- 1:1000
# 非並列化版
> x <- sapply(N, function(i) system.time(as.POSIXct(date...
> summary(x)
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.1750 0.1770 0.1780 0.1795 0.1810 0.2530
# 以下子プロセス数n=1,2,...,10,100で並列計算の経過時間を...
> x <- sapply(N, function(i) system.time(pvec(dates,as.P...
# summary(x)の結果
最小値 第1四分偏差 中央値 平均 第3四分偏差 ...
0.175 0.177 0.178 0.180 0.181 0....
n=1 0.176 0.177 0.178 0.180 0.182 0....
n=2 0.094 0.099 0.100 0.101 0.102 0....
n=3 0.066 0.071 0.072 0.072 0.073 0....
n=4 0.052 0.057 0.059 0.064 0.068 0....
n=5 0.059 0.070 0.072 0.071 0.073 0....
n=6 0.051 0.060 0.061 0.062 0.063 0....
n=7 0.051 0.053 0.054 0.056 0.056 0....
n=8 0.048 0.050 0.052 0.055 0.059 0....
n=9 0.053 0.058 0.061 0.061 0.063 0....
n=10 0.052 0.057 0.059 0.060 0.061 0....
n=100 0.148 0.181 0.192 0.190 0.201 0....
- 推奨パッケージ内のbootはsnow及びmulticoreを使えるように...
#comment
* long vectors 機構の登場 [#w5a197ef]
大規模データ処理の必要性から，バージョン3.0のRよりlong ve...
しかし，これは十分なメモリ(long vectors オブジェクトは結...
一次元ベクトルの長さは依然として 2^{32}-1 以下である．
そうしたベクトルに対する数値演算には，対応関数が long vec...
必要がある．詳細は help(”long vectors”) を参照せよ．
long vectors に対する添字は(整数を表す)倍精度実数であり，...
実際のところ，私はそんなにハイスペックの計算機を持ってい...
#comment
* 国際化の進展 [#qd936bd4]
中間さんのRの日本語化に端を発し，Rの国際化が始まったこと...
現在のRでは特に UTF エンコーディングのサポートが充実して...
\unnnn, \u{nnnn}, \Unnnnnnnn, \U{nnnnnnnn} は
4桁もしくは8桁の16進数を用いたUnicodeによる文字表現を与え...
# キリル(ロシア)文字の例
> c("\u0414","\u{0411}")
[1] "Д" "Б"
Rの推奨パッケージ tools 中には，非アスキー文字や多バイト...
# intToUtf8(946)でもよい
> intToUtf8(0x03B2L)
[1] "β"
# 元に戻る
> utf8ToInt(intToUtf8(0x03B2L))
[1] 946
> intToUtf8(945:950)
[1] "αβγδεζ"
> intToUtf8(945:950,TRUE)
[1] "α" "β" "γ" "δ" "ε" "ζ"
# エキゾチックな半角文字がいっぱい
> x <- intToUtf8(as.integer(
c(160:383,0x0192,0x02C6,0x02C7,0x02CA,
0x02D8,0x02D9,0x02DD,0x200C,0x2018,
0x2019,0x201C,0x201D,0x2020,0x2022,
0x2026,0x20AC)),multiple=TRUE)
# 見やすくする
> matrix(x,ncol=16,byrow=TRUE)
---出力省略---
#comment
* 群盲象を撫でる [#ac32d51b]
R3.0の基本パッケージ中には総計1,264(内部関数は488)個のオ...
#comment
* noquote関数と文字列出力 [#q049f689]
cat関数とprint関数は文字列の出力形式(引用符の有無)が異な...
noquote関数は文字列ベクトルにクラス属性 "noquote" を加え...
> cat(letters[1:3],"\n") # cat関数は文字列を引用符なしで...
a b c
> print(letters[1:3]) # print関数では引用符付き
[1] "a" "b" "c"
> print(noquote(letters[1:3])) # noquote関数を使えばprin...
[1] a b c
# クラス'noquote'を持つ文字列ベクトルになる
> str(noquote(letters[1:3]))
Class 'noquote' chr [1:3] "a" "b" "c"
> noquote(letters)[1:3]
[1] a b c
#comment
* cat関数による出力に制御文字を使う [#a01fd840]
cat関数では制御文字 "\n" が改行の意味を持つことはよく知ら...
cat関数に使用でき，次の出力位置を制御できる制御文字(エス...
\a ベルを鳴らす
\b 一文字分戻る
\f 現在の位置から一行進む
\n 改行し行先頭に
\r 現在の行の先頭に戻る
\v 垂直タブ，一定行分進む
\t 水平タブ，次のタブ位置から行う
これらの制御文字は cat 関数によるコンソール出力で有効であ...
> cat("abc\n")
abc
# 途中で二回連続改行
> cat("abc\n\ndef\nijk\n")
abc
def
ijk
# 行頭に戻って次を出力し最後に改行
> cat("abc\rdef\n")
def
# abcの後，一文字戻ってdefを出力後改行
> cat("abc\bdef\n")
abdef
# print関数では無視される
> print("abc\bdef")
[1] "abc\bdef"
# 2文字戻って次を出力し改行
> cat("abc\b\b\n)
adef
# 3文字戻って次を出力
> cat("abc\b\b\b\n")
def
# abc出力後，行先頭に戻りdefを出力
> cat("abc\rdef\n")
def
# 水平タブを二回使用
> cat("abc\tdef\tijk\n")
abc def ijk
# 垂直タブを二回使用
> cat("abc\vdef\vijk\n")
abc
def
ijk
# 途中で二回行下げ
> cat("abc\fdef\fhij\n")
abc
def
hij
abdehij
> cat("abc\vdef\rijk\n")
abc
ijkdef
> cat("abc\vdef\v\rijk\n")
abc
def
ijk
# 垂直タブを二回連続使用
> cat("abc\vdef\v\v\rijk\n")
abc
def
ijk
画面への出力で行頭に戻って出力する制御文字 "\r" は，何回...
# 途中出力は上書きされ，最後の出力だけが残る
> for(i in 1:100)
{Sys.sleep(2)
cat("\r",i,"番目の出力")}; cat("\n")
100 番目の出力
# 更に空白を除く
> for(i in 1:100)
{Sys.sleep(1)
cat("\r",i,"\b番目の出力")}
cat("\n")
100番目の出力
#comment
* 行列・配列のコンパクトな表示 [#da387bd2]
Rオブジェクトのコンソールへの出力は多くの場合好ましい形式...
行列や配列のコンソールへの(print関数による)表示の際は暗黙...
> no.dimnames <- function(a) {
d <- list(); l <- 0
for(i in dim(a))
d[[l <- l + 1]] <- rep("", i)
dimnames(a) <- d
a }
# 既定の表示
> X <- matrix(1:16,4,4); X
[,1] [,2] [,3] [,4]
[1,] 1 5 9 13
[2,] 2 6 10 14
[3,] 3 7 11 15
[4,] 4 8 12 16
# 次元名なしに表示
> no.dimnames(X)
1 5 9 13
2 6 10 14
3 7 11 15
4 8 12 16
# 配列にも使える
> X <- array(1:8,c(2,2,2)); X
, , 1
[,1] [,2]
[1,] 1 3
[2,] 2 4
, , 2
[,1] [,2]
[1,] 5 7
[2,] 6 8
> no.dimnames(X)
, ,
1 3
2 4
, ,
5 7
6 8
- X <- matrix(1:16,4,4); X は (X <- matrix(1:16,4,4)) の...
- 念のため，このような機能は，一番最初 Version 0.49 Beta ...
#comment
* Lisp(Reduce)風の構文を持つ関数 Reduce, Filter, Find, Po...
Lisp(Reduce)にある機能を真似た関数が幾つか登場しています...
関数 Reduce は二項演算子をベクトル x 中の要素に逐次適用し...
以下で定義する関数 add, cadd はそれぞれ汎用加算関数と汎用...
> add <- function(x) Reduce("+",x)
> add(list(1,2,3))
[1] 6
> cadd <- function(x) Reduce("+",x,accumulate=TRUE)
> cadd(seq_len(7))
[1] 1 3 6 10 15 21 28
> cadd(list(1:3,2:4)) # ベクトルの累積和
[[1]]
[1] 1 2 3
[[2]]
[1] 3 5 7
> cadd(list(matrix(1,2,2), matrix(2,2,2))) # 行列の累...
[[1]]
[,1] [,2]
[1,] 1 1
[2,] 1 1
[[2]]
[,1] [,2]
[1,] 3 3
[2,] 3 3
連分数を計算する．
# 連分数を計算する関数
> cfrac <- function(x) Reduce(function(u,v) u+1/v,x,righ...
> cfrac(c(3,7,15,1,292)) # 円周率を近似する連分数
[1] 3.141593
> cfrac(c(2,1,2,1,1,4,1,1,6,1,1,8)) # exp(1)を近似する連...
[1] 2.718282
関数の重複適用．
> Funcall <- function(f,...) f(...)
# log(exp(acos(cos(0))を計算
> Reduce(Funcall,list(log,exp,acos,cos),0,right=TRUE)
[1] 0
# 黄金比の連分数近似
> cfrac(rep.int(1,31))
[1] 1.618034
# 関数 t |-> (t+x/t)/2 の不動点としてsqrt(x)を近似計算
> asqrt <- function(x,n) Iterate(function(t) (t+x/t)/2,n)
# 初期値を正の数とする
> asqrt(2,30)(10)
[1] 1.414214
# 初期値を負の数とする
> asqrt(2,30)(-1)
[1] -1.414214
#comment
* マニュアル，ヘルプ文章も進化している [#n7eb1874]
Rの魅力の一つが充実したマニュアルとオブジェクト毎のヘルプ...
#comment
* Rの商業的利用に関するR Foundationの見解 [#p46778ac]
Rの商業的利用もかっての，上司に内緒で会社のパソコンにイン...
2.11 Can I use R for commercial purposes?
R is released under the GNU General Public License (GPL)...
If you have any questions regarding the legality of usin...
in any particular situation you should bring it up with ...
We are in no position to offer legal advice.
It is the opinion of the R Core Team that one can use R ...
purposes (e.g., in business or in consulting). The GPL,...
licenses, permits all and any use of the package. It onl...
distribution of R or of other programs containing code f...
clear in clause 6 (“No Discrimination Against Fields of ...
the Open Source Definition:
The license must not restrict anyone from making use...
in a specific field of endeavor. For example, it may not...
from being used in a business, or from being used for g...
It is also explicitly stated in clause 0 of the GPL, whi...
Activities other than copying, distribution and modi...
covered by this License; they are outside its scope. The...
the Program is not restricted, and the output from the P...
only if its contents constitute a work based on the Prog...
Most add-on packages, including all recommended ones, al...
commercial use in this way. A few packages are restricte...
use”; you should contact the author to clarify whether t...
seek the advice of your legal counsel.
None of the discussion in this section constitutes legal...
Team does not provide legal advice under any circumstanc...
#comment
* パッケージ data.table [#wc16f06d]
データフレームはRのデータ形式の中心であり，多くの統計処理...
パッケージ data.table (Rの推奨パッケージではない)はデータ...
以下は同じ内容の大きなデータフレームとデータテーブルの一...
> x <- sample(1:10, 1e6,rep=TRUE)
> y <- sample(letters[1:10], 1e6,rep=TRUE)
> DF <- data.frame(x=x,y=y)
> DT <- as.data.table(DF)
# データフレーム版
> z <- subset(DF,y=="a", select=x)
# データテーブル版
> z <- DT[y=="a"]
# キー指定したデータテーブル版
> setkey(DT,y)
> z <- DT["a"]
次の例はデータテーブルの操作がデータフレームのそれとは異...
> x <- sample(0:9,1e5,replace=TRUE)
> y <- sample(letters[0:9],1e5,replace=TRUE))
> z <- runif(1e5)
> DF <- data.frame(x=x,y=y)
> tracemem(DT)
[1] "<0x9d195e0>"
# データフレームの連結操作．２回の内部コピー
> DF <- cbind(DF,z=z)
tracemem[0x9d195e0 -> 0xcc2e068]: data.frame cbind cbind
tracemem[0xcc2e068 -> 0xcc2c958]: data.frame cbind cbind
> DF <- data.frame(x=x,y=y)
> tracemem(DF)
[1] "<0xad5ad58>"
# 同じ事を別の操作で．4回の内部コピー
> DF <- transform(DF,z=z)
tracemem[0xb36b928 -> 0xb367928]: do.call transform.data...
tracemem[0xb367928 -> 0xb367a88]: do.call transform.data...
tracemem[0xb367a88 -> 0xb367d68]: data.frame do.call tra...
tracemem[0xb367d68 -> 0xb364658]: data.frame do.call tra...
> DT <- data.table(x=x,y=y)
> tracemem(DT)
[1] "<0x8354770>"
# データテーブルに新しい列を加える．コピーされていない
> DT[,z:=z]
#comment
* print関数のオプション zero.print [#m2fdd112]
Rのコンソールへの出力は，実はprint関数の無数のメソッド関...
> t1 <- round(abs(rt(200,df=1.8)))
> t2 <- round(abs(rt(200,df=1.4)))
# メソッドprint.table使用
> table(t1,t2)
t2
t1 0 1 2 3 4 5 6 7 8 10 17 21 30
0 21 22 14 4 1 0 1 1 1 0 1 0 0
1 25 21 7 3 4 2 1 1 1 1 0 0 0
--途中省略--
12 1 0 0 0 0 0 0 0 0 0 0 0 0
# 値0をドットで表現
> print(table(t1,t2),zero.print=".")
t2
t1 0 1 2 3 4 5 6 7 8 10 17 21 30
0 21 22 14 4 1 . 1 1 1 . 1 . .
1 25 21 7 3 4 2 1 1 1 1 . . .
--途中省略--
12 1 . . . . . . . . . . . .
- Version 0.60 Alpha (December 2, 1997) から使えるように...
#comment
* こんなのあったけ？ 行列・配列編 [#n71b1e69]
今回 R-fullrefman.pdf をあれこれ見ているうちに「こんなの...
関数 rowsum(x,group,reorder=TRUE,...) (rowSumsではありま...
> rowsum(1:5,c(1,2,2,2,2)) # 各列をグループ1,2に分けて...
[,1]
1 1 # グループ1の総和
2 14 # グループ2の総和
> rowsum(1:5,c(1,2,2,3,3)) # グループ数3
[,1]
1 1 # 1
2 5 # 2+3
3 9 # 4+5
> x
[,1] [,2] [,3] [,4]
[1,] 1 5 9 13
[2,] 2 6 10 14
[3,] 3 7 11 15
[4,] 4 8 12 16
> rowsum(x,c(1,2,2,3))
[,1] [,2] [,3] [,4]
1 1 5 9 13 # グループ番号1の行の列和 (第1列...
2 5 13 21 29 # グループ番号2の行の列和 (第2,3...
3 4 8 12 16 # グループ番号3の行の列和 (第4列...
> rowsum(x,c(1,2,2,2))
[,1] [,2] [,3] [,4]
1 1 5 9 13 # グループ番号1の行の列総和 (第1...
2 9 21 33 45 # グループ番号2の行の列総和(第2,3...
> rowsum(x, c("a","b","a","b")) # グループ変数は何でも...
[,1] [,2] [,3] [,4]
a 4 12 20 28
b 6 14 22 30
ベクトルは行列と異なり次元属性を持たないので nrow, ncol ...
> y <- 1:(3*4*5)
> c(NROW(y),NCOL(y)) # ベクトルに対するNROW,NCOL
[1] 60 1
配列 x とその次元番号 MARGIN に対する slice.index(x,MARGI...
> x <- array(1:24,c(2,3,4))
> x1 <- slice.index(x,1)
# x1[n,i,j]はn
> x1[1,,]
[,1] [,2] [,3] [,4]
[1,] 1 1 1 1
[2,] 1 1 1 1
[3,] 1 1 1 1
> x1[2,,]
[,1] [,2] [,3] [,4]
[1,] 2 2 2 2
[2,] 2 2 2 2
[3,] 2 2 2 2
> x <- array(1:8,c(2,2,2))
> x1 <- slice.index(x,1)
> x2 <- slice.index(x,2)
> x3 <- slice.index(x,3)
> x[x2 == x1 & x3 == x1] # 配列の一般化対角成分 x[1,1,...
[1] 1 8
- %%Version 1.8.1 (2003-11-21) には既にありました。%%&br;
rowsum は Version 0.63.1 (December 5, 1998) からです。&br;
slice.index は Version 1.7.0 (2003-04-16) には既にありま...
#comment
* こんなのあったけ？ apply関数族編 [#h30d367b]
vapply 関数は出力書式指定の sapply 関数だそうです．
関数 vapply(X,FUN,FUN.VALUE,...,USE.NAMES=TRUE) は sapply...
> i39 <- sapply(3:9,seq); i39 # ベクトルのリスト
[[1]]
[1] 1 2 3
[[2]]
[1] 1 2 3 4
--途中省略--
[[7]]
[1] 1 2 3 4 5 6 7 8 9
> vapply(i39,fivenum,c(0,0,0,0,0)) # 一つの結果は長さ5...
[,1] [,2] [,3] [,4] [,5] [,6] [,7]
[1,] 1.0 1.0 1 1.0 1.0 1.0 1
[2,] 1.5 1.5 2 2.0 2.5 2.5 3
[3,] 2.0 2.5 3 3.5 4.0 4.5 5
[4,] 2.5 3.5 4 5.0 5.5 6.5 7
[5,] 3.0 4.0 5 6.0 7.0 8.0 9
#comment
* こんなのあったけ？ 作表関数編 [#ded7d377]
関数 prop.table(x,margin=NULL) は分割表の各項目を margin ...
> m <- matrix(1:4,2)
> margin.table(m,1)
[1] 4 6
> margin.table(m,2)
[1] 3 7
> prop.table(m,1)
[,1] [,2]
[1,] 0.2500000 0.7500000
[2,] 0.3333333 0.6666667
> prop.table(m,2)
[,1] [,2]
[1,] 0.3333333 0.4285714
[2,] 0.6666667 0.5714286
- %%Version 1.8.1 (2003-11-21) には既にありました。%%&b...
function (x, margin)
sweep(x, margin, margin.table(x, margin), "/")
でした（今も大差ない） -- &new{2013-10-24 (木) 21:10:32};
> example(prop.table)
prp.tb> m <- matrix(1:4, 2)
prp.tb> m
[,1] [,2]
[1,] 1 3
[2,] 2 4
prp.tb> prop.table(m, 1)
[,1] [,2]
[1,] 0.2500000 0.7500000
[2,] 0.3333333 0.6666667
#comment
* こんなのあったけ？ タイマー編 [#u0512bbf]
setTimeLimit はトップレベルの各計算(つまり，コマンド行入...
> setSessionTimeLimit(10,10)
> for(i in seq(1e4)) {
cat(i); Sys.sleep(1)}
12345678910 以下にエラー Sys.sleep(1) : セッション時間が...
> {setTimeLimit(10,10);
for(i in seq(1e4))
{cat(i); Sys.sleep(1)}}
12345678910 以下にエラー Sys.sleep(1) : 時間が経過して上...
gc.time はガベージコレクションを行い，必要だった時間を返...
# 巨大なベクトルを作り，消してからガベージコレクション時...
> sapply(1:10, function(i) {x <- rnorm(1e8); rm(x); gc.t...
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,...
[1,] 1.574 1.584 1.602 1.612 1.622 1.644 1.654 1.664 1.6...
[2,] 1.026 1.072 1.122 1.168 1.218 1.268 1.318 1.364 1.4...
[3,] 1.097 1.138 1.190 1.231 1.272 1.326 1.367 1.407 1.4...
[4,] 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.0...
[5,] 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.0...
- R version 2.9.1 (2009-06-26) に既にありました.&br;なか...
- ご指摘ありがとうございます．ここで「あったけ？」と言っ...
- 全てのバージョンのRをインストールしておけばよい。あるバ...
#comment
*遅延コピー(とでもいうんだろうか?) [#id742d03]
help(tracemem)を見ていていまさらだが気づいたこと(当たり前...
> a <- 1:10
> tracemem(a) # aが指すメモリ領域がメモリ上でコピー...
[1] "<0x990d438>"
> b <- a
> b[1:3] # この段階ではa,bは同じメモリ領域を共...
[1] 1 2 3
> sum(b) # 同じく
[1] 55
> a[1] <- 10 # a が変更されたので実際のコピーが必...
tracemem[0x990d438 -> 0x990d558]: # 二度コピーされてい...
tracemem[0x990d558 -> 0x8ea59e8]:
> untracemem(a)
> a <- 1:10
> tracemem(a)
[1] "<0xa3e9cb0>"
> b <- a
> b[1] <- 10 # bが変更されればやはり実際のコピーが...
tracemem[0xa3e9cb0 -> 0xa3e9c68]:
tracemem[0xa3e9c68 -> 0xa2c5e08]:
> untracemem(a)
> tracemem(a)
[1] "<0xa3e9cb0>"
> b <- a
> rm(a) # コピー動作は無い，考えてみれば当然
#comment
*S3, S4クラス，メソッド [#p86303ff]
S言語本家の Chambers 氏が開発メンバーに加わったせいがある...
#comment
*豆知識 #line directive [#s3ee04a2]
Rでは，#記号から行末まではコメントとされ実行時には無視さ...
#line nn "ファイル名"
(#line は行頭から5文字にあるべき)はC言語のそれと同じく，p...
R news の記事の一部を引用
The #line directive
In some cases, R source code is written by a program,
not by a human being. For example, Sweave() extracts
lines of code from Sweave documents before sending
the lines to R for parsing and evaluation. To support
such preprocessors, the R 2.10.0 parser recognizes
a new directive of the form
#line nn "filename"
where nn is an integer. As with the same-named
directive in the C language, this tells the parser to
assume that the next line of source is line nn
from the given filename for the purpose of constructing
source references. The Sweave() function doesn’t
currently make use of this, but in the future, it (and
other preprocessors) could output #line directives
so that source references and syntax errors refer to
the original source location rather than to an inter-
mediate file.
The #line directive was a late addition to R
2.10.0. Support for this in Sweave() appeared in R 2.12.0
#comment
* ローマ数字の話 [#i55bfa5c]
Rには正整数をローマ数字に変換する関数があります．いまさら...
関数 as.roman(n) は正整数 n をローマ数字表記に変換する．...
# ローマ数字表記
> x <- as.roman(c(1,5,10,100)); x
[1] I V X C
# ローマ数字には0や負の数は無い
> as.roman(c(-1,0))
[1] <NA> <NA>
# クラス'roman'を持つ整数である
> str(x)
Class 'roman' int [1:6] 1 5 10 100
# ローマ数字に対する四則演算例
> as.roman(10)+1
[1] XI
> as.roman(10)-1
[1] IX
# ローマ数字では3,899までが一意的に表現可能
# 4,000以上の数は表現できない
> as.roman(3899); as.roman(3900)
[1] MMMDCCCXCIX
[1] <NA>
# ただの整数に戻る
> c(as.roman(3),as.roman(10))
[1] 3 10
- http://www.unicode.org/charts/nameslist/n_2150.html に...
- as.roman のヘルプが，http://en.wikipedia.org/w/index.ph...
#comment
終了行:
東京工業大学の[[間瀬茂]]です．拙著「Rプログラミングマニュ...
拙著の基本方針は「Rのプログラミング機能を網羅的に紹介する...
最後の仕事になります)．出版時期はまだ未定です．
改訂作業中に特に「!!!」と思った点が幾つかありました．これ...
もし，私の誤解，理解不足がありましたら，コメントいただけ...
注意：頁タイトルが誤解を与えるかもしれませんが，ここに書...
(独り言：Rバージョン4が出る頃，私はまだ生きているかな......
----
#contents
----
*「apply関数族は簡潔な既述ができるが，実体は for ループを...
lapply, vapply, sapply, rapply 関数は内部関数を用いて計算...
以下は sapply 関数と for ループを使った同値な処理をそれぞ...
foo <- function(x) x # つまらない関数
X <- numeric(1e6) # 念のため，予め変数を用意
# sapply関数を用いた処理を100回繰り返す
Y <- sapply(1:100, function(i) system.time(X <- sapply(1...
# そのユーザ時間の要約
> c(fivenum(Y[1,]),mean(Y[1,]),sd(Y[1,]))
[1] 2.048 2.240 2.280 2.388 2.604 2.303 0.09771
# forループを用いた処理を100回繰り返す
Z <- sapply(1:100, function(x) system.time(for(i in 1:1e...
# そのユーザ時間要約
> c(fivenum(Z[1,]),mean(Z[1,]),sd(Z[1,]))
[1] 2.548 2.592 2.606 2.624 2.684 2.609 0.02449
- OS やメモリーなどの環境によるのかもしれませんが，以下の...
示されたプログラムでは system.time(X <- sapply(1:1e6,foo)...
foo <- function(x) x # つまらない関数
sim <- function(n, m) {
X <- numeric(m) # 念のため，予め変数を用意
Y <- replicate(n, {gc(); gc(); system.time(X <- sapply(...
Z <- replicate(n, {gc(); gc(); system.time(for(i in 1:m...
return(list(s=Y, f=Z))
}
得られた測定データの最小値，最大値，平均値，標準偏差を以...
n=1e4 n=1e6
sapply for sapply for
min 0.00700 0.01600 2.11900 1.73000
max 0.01900 0.02300 2.91500 2.76100
mean 0.00868 < 0.01838 2.45032 > 1.84305
sd 0.00140 0.00183 0.12012 0.22043
n=1e4 のとき（左図）は for は sapply の倍くらい時間がかか...
&ref(fig1.png); &ref(fig2.png);
- 早速のコメントありがとうございます．実行環境により様子...
n=1e5 sapply
0.132 0.140 0.144 0.148 0.180 0.145 0.00770
n=1e5 for ループ
0.248 0.256 0.260 0.264 0.292 0.260 0.00780
n=1e4 sapply
0.0120 0.0120 0.0120 0.0140 0.0160 0.0130 0.00174
n=1e4 forループ
0.0240 0.0240 0.0240 0.0280 0.0320 0.0253 0.00205
- あ, sapplyはlapplyのただのラッパーです. 性能のバラツキ...
#comment
*「Rシステム」はあるが「R言語」は無い． [#x5862e1a]
現在「R言語」という呼び方がしばしばみられますが，私は二重...
例えば，公式マニュアル An Introduction to R には
”a well developed, simple and effective programming lang...
which includes conditionals, loops, user defined recursi...
and output facilities. (Indeed most of the system suppli...
themselves written in the S language.)”
公式FAQ集の記述では
"We can regard S as a language with three current implem...
or“engines”, the “old S engine” (S version 3; S-PLUS 3.x...
the “new S engine”(S version 4; S-PLUS 5.x and above), a...
Given this understanding, asking for “the differences be...
really amounts to asking for the specifics of the R impl...
the S language, i.e., the difference between the R and S...
#comment
* Tierny氏の compiler パッケージが推奨パッケージになった ...
R開発メンバーのTierny氏の compiler パッケージのおかげで，...
関数 enableJIT は実行時(Just-In-Time)コンパイルを指示しま...
関数 compilePKGS はパッケージを読み込む際に，含まれる関数...
# lapply関数の初期の定義
> la1 <- function(X,FUN,...) {
FUN <- match.fun(FUN)
if (!is.list(X))
X <- as.list(X)
rval <- vector("list",length(X))
for(i in seq(along=X))
rval[i] <- list(FUN(X[[i]],...))
names(rval) <- names(X)
return(rval) }
# 時間のかかるまずいコード
> tmp <- function(n) {
x <- numeric(0)
for(i in 1:n)
x <- c(x,runif(1))
return(x) }
# コンパイルする
> la1c <- cmpfun(la1)
> tmpc <- cmpfun(tmp)
> la1c
function(X,FUN, ...) {
FUN <- match.fun(FUN)
if (!is.list(X)) X <- as.list(X)
rval <- vector("list",length(X))
for(i in seq(along=X)) rval[i] <- list(FUN(X[[i]],...))
names(rval) <- names(X)
return(rval) }
<bytecode: 0xa9a3388> # バイトコンパイルされている証拠
> tmpc
function(n) {x <- numeric(0); for(i in 1:n) x <- c(x,run...
<bytecode: 0xa64e254>
# 以下は la1(1:100,tmp) を1万回実行したユーザ時間の要約
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.016 0.020 0.020 0.021 0.020 0.024 # JITコンパ...
0.004 0.012 0.012 0.011 0.012 0.024 # enableJIT(1)
0.004 0.012 0.012 0.011 0.012 0.016 # enableJIT(2)
0.008 0.012 0.012 0.011 0.012 0.016 # enableJIT(3)
#comment
* バイト型データと関連関数の導入 [#o50edeef]
Rの原子オブジェクトに新しく(何時から？)バイト(raw)型デー...
# 255は16進法ではff
> x <- as.raw(255); x
[1] ff
> str(x)
raw ff
> x <- "A test string"
# xのバイト列表現
> y <- charToRaw(x); y
[1] 41 20 74 65 73 74 20 73 74 72 69 6e 67
# 文字列に戻す
> rawToChar(y)
[1] "A test string"
# 整数8のビット列表現
> intToBits(8L)
[1] 00 00 00 01 00 00 00 00 00 00 00 00 00 00 00 00 00 ...
[26] 00 00 00 00 00 00 00
#comment
* 隠し変数 .Last.value [#oa840bee]
Rのトップレベル(つまりRコンソールのプロンプト > に対して...
続けて二回実行することが勧められています．
> x <- log(10); x
[1] 2.302585
> .Last.value
[1] 2.302585
# library関数の返り値は不可視
> library(splines)
# .Last.valueには記録
> .Last.value
[1] "splines" "stats"
--以下省略--
# 巨大なオブジェクトx
> x <- rnorm(1e6)
# .Last.valueにも記録
> identical(x,.Last.value)
[1] TRUE
# この時点で.Last.valueはすでに更新されている
> str(.Last.value)
TRUE
ユーザレベルでも .Last.value 変数の便利な使い方があります...
> sin(1)
[1] 0.841471
> res <- .Last.value; res
[1] 0.841471
#comment
* マルチコア並列処理パッケージ parallel が推奨パッケージ...
Rの推奨パッケージに新しく付け加わった並列処理パッケージ p...
# クラスタを作り並列処理する例
> require(parallel)
# 4コアのクラスタを生成
> cl <- makeCluster(4)
# 変数xxと関数fooを定義
> xx <- 1; foo <- function(y) xx + y
# 小プロセスにそれを移出，オブジェクト名は文字列で与える
> clusterExport(cl,c("xx","foo"))
# sapply関数の並列版parSapplyを用いfoo(1:10)をクラスタに...
> parSapply(cl,1:10,foo)
[1] 2 3 4 5 6 7 8 9 10 11
# 2変数関数を定義し，小プロセスにそれを移出
> bar <- function(x,y) x + y
> clusterExport(cl,"bar")
# 多変数apply関数であるmapply関数の並列版clusterMap．以...
> clusterMap(cl,bar,c(a=1,b=2,c=3),
c(10,0,-10))
$a # つまりbar(a=1,10)
[1] 11
$b
[1] 2
$c
[1] -7
> stopCluster(cl) # 最後にクラスタを閉じる
関数 pvec はクラスタの生成，apply 関数風の並列計算と結果...
次の結果は日付文字列一万個を as.POSIXct 関数で日付オブジ...
# 日付文字列一万個のベクトルdatesをas.POSIXctで日付オブ...
# すべて経過時間を示す
> dates <- sprintf('%04d-%02d-%02d',as.integer(2000+rnor...
as.integer(runif(1e4,1,12)),as.intege...
> N <- 1:1000
# 非並列化版
> x <- sapply(N, function(i) system.time(as.POSIXct(date...
> summary(x)
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.1750 0.1770 0.1780 0.1795 0.1810 0.2530
# 以下子プロセス数n=1,2,...,10,100で並列計算の経過時間を...
> x <- sapply(N, function(i) system.time(pvec(dates,as.P...
# summary(x)の結果
最小値 第1四分偏差 中央値 平均 第3四分偏差 ...
0.175 0.177 0.178 0.180 0.181 0....
n=1 0.176 0.177 0.178 0.180 0.182 0....
n=2 0.094 0.099 0.100 0.101 0.102 0....
n=3 0.066 0.071 0.072 0.072 0.073 0....
n=4 0.052 0.057 0.059 0.064 0.068 0....
n=5 0.059 0.070 0.072 0.071 0.073 0....
n=6 0.051 0.060 0.061 0.062 0.063 0....
n=7 0.051 0.053 0.054 0.056 0.056 0....
n=8 0.048 0.050 0.052 0.055 0.059 0....
n=9 0.053 0.058 0.061 0.061 0.063 0....
n=10 0.052 0.057 0.059 0.060 0.061 0....
n=100 0.148 0.181 0.192 0.190 0.201 0....
- 推奨パッケージ内のbootはsnow及びmulticoreを使えるように...
#comment
* long vectors 機構の登場 [#w5a197ef]
大規模データ処理の必要性から，バージョン3.0のRよりlong ve...
しかし，これは十分なメモリ(long vectors オブジェクトは結...
一次元ベクトルの長さは依然として 2^{32}-1 以下である．
そうしたベクトルに対する数値演算には，対応関数が long vec...
必要がある．詳細は help(”long vectors”) を参照せよ．
long vectors に対する添字は(整数を表す)倍精度実数であり，...
実際のところ，私はそんなにハイスペックの計算機を持ってい...
#comment
* 国際化の進展 [#qd936bd4]
中間さんのRの日本語化に端を発し，Rの国際化が始まったこと...
現在のRでは特に UTF エンコーディングのサポートが充実して...
\unnnn, \u{nnnn}, \Unnnnnnnn, \U{nnnnnnnn} は
4桁もしくは8桁の16進数を用いたUnicodeによる文字表現を与え...
# キリル(ロシア)文字の例
> c("\u0414","\u{0411}")
[1] "Д" "Б"
Rの推奨パッケージ tools 中には，非アスキー文字や多バイト...
# intToUtf8(946)でもよい
> intToUtf8(0x03B2L)
[1] "β"
# 元に戻る
> utf8ToInt(intToUtf8(0x03B2L))
[1] 946
> intToUtf8(945:950)
[1] "αβγδεζ"
> intToUtf8(945:950,TRUE)
[1] "α" "β" "γ" "δ" "ε" "ζ"
# エキゾチックな半角文字がいっぱい
> x <- intToUtf8(as.integer(
c(160:383,0x0192,0x02C6,0x02C7,0x02CA,
0x02D8,0x02D9,0x02DD,0x200C,0x2018,
0x2019,0x201C,0x201D,0x2020,0x2022,
0x2026,0x20AC)),multiple=TRUE)
# 見やすくする
> matrix(x,ncol=16,byrow=TRUE)
---出力省略---
#comment
* 群盲象を撫でる [#ac32d51b]
R3.0の基本パッケージ中には総計1,264(内部関数は488)個のオ...
#comment
* noquote関数と文字列出力 [#q049f689]
cat関数とprint関数は文字列の出力形式(引用符の有無)が異な...
noquote関数は文字列ベクトルにクラス属性 "noquote" を加え...
> cat(letters[1:3],"\n") # cat関数は文字列を引用符なしで...
a b c
> print(letters[1:3]) # print関数では引用符付き
[1] "a" "b" "c"
> print(noquote(letters[1:3])) # noquote関数を使えばprin...
[1] a b c
# クラス'noquote'を持つ文字列ベクトルになる
> str(noquote(letters[1:3]))
Class 'noquote' chr [1:3] "a" "b" "c"
> noquote(letters)[1:3]
[1] a b c
#comment
* cat関数による出力に制御文字を使う [#a01fd840]
cat関数では制御文字 "\n" が改行の意味を持つことはよく知ら...
cat関数に使用でき，次の出力位置を制御できる制御文字(エス...
\a ベルを鳴らす
\b 一文字分戻る
\f 現在の位置から一行進む
\n 改行し行先頭に
\r 現在の行の先頭に戻る
\v 垂直タブ，一定行分進む
\t 水平タブ，次のタブ位置から行う
これらの制御文字は cat 関数によるコンソール出力で有効であ...
> cat("abc\n")
abc
# 途中で二回連続改行
> cat("abc\n\ndef\nijk\n")
abc
def
ijk
# 行頭に戻って次を出力し最後に改行
> cat("abc\rdef\n")
def
# abcの後，一文字戻ってdefを出力後改行
> cat("abc\bdef\n")
abdef
# print関数では無視される
> print("abc\bdef")
[1] "abc\bdef"
# 2文字戻って次を出力し改行
> cat("abc\b\b\n)
adef
# 3文字戻って次を出力
> cat("abc\b\b\b\n")
def
# abc出力後，行先頭に戻りdefを出力
> cat("abc\rdef\n")
def
# 水平タブを二回使用
> cat("abc\tdef\tijk\n")
abc def ijk
# 垂直タブを二回使用
> cat("abc\vdef\vijk\n")
abc
def
ijk
# 途中で二回行下げ
> cat("abc\fdef\fhij\n")
abc
def
hij
abdehij
> cat("abc\vdef\rijk\n")
abc
ijkdef
> cat("abc\vdef\v\rijk\n")
abc
def
ijk
# 垂直タブを二回連続使用
> cat("abc\vdef\v\v\rijk\n")
abc
def
ijk
画面への出力で行頭に戻って出力する制御文字 "\r" は，何回...
# 途中出力は上書きされ，最後の出力だけが残る
> for(i in 1:100)
{Sys.sleep(2)
cat("\r",i,"番目の出力")}; cat("\n")
100 番目の出力
# 更に空白を除く
> for(i in 1:100)
{Sys.sleep(1)
cat("\r",i,"\b番目の出力")}
cat("\n")
100番目の出力
#comment
* 行列・配列のコンパクトな表示 [#da387bd2]
Rオブジェクトのコンソールへの出力は多くの場合好ましい形式...
行列や配列のコンソールへの(print関数による)表示の際は暗黙...
> no.dimnames <- function(a) {
d <- list(); l <- 0
for(i in dim(a))
d[[l <- l + 1]] <- rep("", i)
dimnames(a) <- d
a }
# 既定の表示
> X <- matrix(1:16,4,4); X
[,1] [,2] [,3] [,4]
[1,] 1 5 9 13
[2,] 2 6 10 14
[3,] 3 7 11 15
[4,] 4 8 12 16
# 次元名なしに表示
> no.dimnames(X)
1 5 9 13
2 6 10 14
3 7 11 15
4 8 12 16
# 配列にも使える
> X <- array(1:8,c(2,2,2)); X
, , 1
[,1] [,2]
[1,] 1 3
[2,] 2 4
, , 2
[,1] [,2]
[1,] 5 7
[2,] 6 8
> no.dimnames(X)
, ,
1 3
2 4
, ,
5 7
6 8
- X <- matrix(1:16,4,4); X は (X <- matrix(1:16,4,4)) の...
- 念のため，このような機能は，一番最初 Version 0.49 Beta ...
#comment
* Lisp(Reduce)風の構文を持つ関数 Reduce, Filter, Find, Po...
Lisp(Reduce)にある機能を真似た関数が幾つか登場しています...
関数 Reduce は二項演算子をベクトル x 中の要素に逐次適用し...
以下で定義する関数 add, cadd はそれぞれ汎用加算関数と汎用...
> add <- function(x) Reduce("+",x)
> add(list(1,2,3))
[1] 6
> cadd <- function(x) Reduce("+",x,accumulate=TRUE)
> cadd(seq_len(7))
[1] 1 3 6 10 15 21 28
> cadd(list(1:3,2:4)) # ベクトルの累積和
[[1]]
[1] 1 2 3
[[2]]
[1] 3 5 7
> cadd(list(matrix(1,2,2), matrix(2,2,2))) # 行列の累...
[[1]]
[,1] [,2]
[1,] 1 1
[2,] 1 1
[[2]]
[,1] [,2]
[1,] 3 3
[2,] 3 3
連分数を計算する．
# 連分数を計算する関数
> cfrac <- function(x) Reduce(function(u,v) u+1/v,x,righ...
> cfrac(c(3,7,15,1,292)) # 円周率を近似する連分数
[1] 3.141593
> cfrac(c(2,1,2,1,1,4,1,1,6,1,1,8)) # exp(1)を近似する連...
[1] 2.718282
関数の重複適用．
> Funcall <- function(f,...) f(...)
# log(exp(acos(cos(0))を計算
> Reduce(Funcall,list(log,exp,acos,cos),0,right=TRUE)
[1] 0
# 黄金比の連分数近似
> cfrac(rep.int(1,31))
[1] 1.618034
# 関数 t |-> (t+x/t)/2 の不動点としてsqrt(x)を近似計算
> asqrt <- function(x,n) Iterate(function(t) (t+x/t)/2,n)
# 初期値を正の数とする
> asqrt(2,30)(10)
[1] 1.414214
# 初期値を負の数とする
> asqrt(2,30)(-1)
[1] -1.414214
#comment
* マニュアル，ヘルプ文章も進化している [#n7eb1874]
Rの魅力の一つが充実したマニュアルとオブジェクト毎のヘルプ...
#comment
* Rの商業的利用に関するR Foundationの見解 [#p46778ac]
Rの商業的利用もかっての，上司に内緒で会社のパソコンにイン...
2.11 Can I use R for commercial purposes?
R is released under the GNU General Public License (GPL)...
If you have any questions regarding the legality of usin...
in any particular situation you should bring it up with ...
We are in no position to offer legal advice.
It is the opinion of the R Core Team that one can use R ...
purposes (e.g., in business or in consulting). The GPL,...
licenses, permits all and any use of the package. It onl...
distribution of R or of other programs containing code f...
clear in clause 6 (“No Discrimination Against Fields of ...
the Open Source Definition:
The license must not restrict anyone from making use...
in a specific field of endeavor. For example, it may not...
from being used in a business, or from being used for g...
It is also explicitly stated in clause 0 of the GPL, whi...
Activities other than copying, distribution and modi...
covered by this License; they are outside its scope. The...
the Program is not restricted, and the output from the P...
only if its contents constitute a work based on the Prog...
Most add-on packages, including all recommended ones, al...
commercial use in this way. A few packages are restricte...
use”; you should contact the author to clarify whether t...
seek the advice of your legal counsel.
None of the discussion in this section constitutes legal...
Team does not provide legal advice under any circumstanc...
#comment
* パッケージ data.table [#wc16f06d]
データフレームはRのデータ形式の中心であり，多くの統計処理...
パッケージ data.table (Rの推奨パッケージではない)はデータ...
以下は同じ内容の大きなデータフレームとデータテーブルの一...
> x <- sample(1:10, 1e6,rep=TRUE)
> y <- sample(letters[1:10], 1e6,rep=TRUE)
> DF <- data.frame(x=x,y=y)
> DT <- as.data.table(DF)
# データフレーム版
> z <- subset(DF,y=="a", select=x)
# データテーブル版
> z <- DT[y=="a"]
# キー指定したデータテーブル版
> setkey(DT,y)
> z <- DT["a"]
次の例はデータテーブルの操作がデータフレームのそれとは異...
> x <- sample(0:9,1e5,replace=TRUE)
> y <- sample(letters[0:9],1e5,replace=TRUE))
> z <- runif(1e5)
> DF <- data.frame(x=x,y=y)
> tracemem(DT)
[1] "<0x9d195e0>"
# データフレームの連結操作．２回の内部コピー
> DF <- cbind(DF,z=z)
tracemem[0x9d195e0 -> 0xcc2e068]: data.frame cbind cbind
tracemem[0xcc2e068 -> 0xcc2c958]: data.frame cbind cbind
> DF <- data.frame(x=x,y=y)
> tracemem(DF)
[1] "<0xad5ad58>"
# 同じ事を別の操作で．4回の内部コピー
> DF <- transform(DF,z=z)
tracemem[0xb36b928 -> 0xb367928]: do.call transform.data...
tracemem[0xb367928 -> 0xb367a88]: do.call transform.data...
tracemem[0xb367a88 -> 0xb367d68]: data.frame do.call tra...
tracemem[0xb367d68 -> 0xb364658]: data.frame do.call tra...
> DT <- data.table(x=x,y=y)
> tracemem(DT)
[1] "<0x8354770>"
# データテーブルに新しい列を加える．コピーされていない
> DT[,z:=z]
#comment
* print関数のオプション zero.print [#m2fdd112]
Rのコンソールへの出力は，実はprint関数の無数のメソッド関...
> t1 <- round(abs(rt(200,df=1.8)))
> t2 <- round(abs(rt(200,df=1.4)))
# メソッドprint.table使用
> table(t1,t2)
t2
t1 0 1 2 3 4 5 6 7 8 10 17 21 30
0 21 22 14 4 1 0 1 1 1 0 1 0 0
1 25 21 7 3 4 2 1 1 1 1 0 0 0
--途中省略--
12 1 0 0 0 0 0 0 0 0 0 0 0 0
# 値0をドットで表現
> print(table(t1,t2),zero.print=".")
t2
t1 0 1 2 3 4 5 6 7 8 10 17 21 30
0 21 22 14 4 1 . 1 1 1 . 1 . .
1 25 21 7 3 4 2 1 1 1 1 . . .
--途中省略--
12 1 . . . . . . . . . . . .
- Version 0.60 Alpha (December 2, 1997) から使えるように...
#comment
* こんなのあったけ？ 行列・配列編 [#n71b1e69]
今回 R-fullrefman.pdf をあれこれ見ているうちに「こんなの...
関数 rowsum(x,group,reorder=TRUE,...) (rowSumsではありま...
> rowsum(1:5,c(1,2,2,2,2)) # 各列をグループ1,2に分けて...
[,1]
1 1 # グループ1の総和
2 14 # グループ2の総和
> rowsum(1:5,c(1,2,2,3,3)) # グループ数3
[,1]
1 1 # 1
2 5 # 2+3
3 9 # 4+5
> x
[,1] [,2] [,3] [,4]
[1,] 1 5 9 13
[2,] 2 6 10 14
[3,] 3 7 11 15
[4,] 4 8 12 16
> rowsum(x,c(1,2,2,3))
[,1] [,2] [,3] [,4]
1 1 5 9 13 # グループ番号1の行の列和 (第1列...
2 5 13 21 29 # グループ番号2の行の列和 (第2,3...
3 4 8 12 16 # グループ番号3の行の列和 (第4列...
> rowsum(x,c(1,2,2,2))
[,1] [,2] [,3] [,4]
1 1 5 9 13 # グループ番号1の行の列総和 (第1...
2 9 21 33 45 # グループ番号2の行の列総和(第2,3...
> rowsum(x, c("a","b","a","b")) # グループ変数は何でも...
[,1] [,2] [,3] [,4]
a 4 12 20 28
b 6 14 22 30
ベクトルは行列と異なり次元属性を持たないので nrow, ncol ...
> y <- 1:(3*4*5)
> c(NROW(y),NCOL(y)) # ベクトルに対するNROW,NCOL
[1] 60 1
配列 x とその次元番号 MARGIN に対する slice.index(x,MARGI...
> x <- array(1:24,c(2,3,4))
> x1 <- slice.index(x,1)
# x1[n,i,j]はn
> x1[1,,]
[,1] [,2] [,3] [,4]
[1,] 1 1 1 1
[2,] 1 1 1 1
[3,] 1 1 1 1
> x1[2,,]
[,1] [,2] [,3] [,4]
[1,] 2 2 2 2
[2,] 2 2 2 2
[3,] 2 2 2 2
> x <- array(1:8,c(2,2,2))
> x1 <- slice.index(x,1)
> x2 <- slice.index(x,2)
> x3 <- slice.index(x,3)
> x[x2 == x1 & x3 == x1] # 配列の一般化対角成分 x[1,1,...
[1] 1 8
- %%Version 1.8.1 (2003-11-21) には既にありました。%%&br;
rowsum は Version 0.63.1 (December 5, 1998) からです。&br;
slice.index は Version 1.7.0 (2003-04-16) には既にありま...
#comment
* こんなのあったけ？ apply関数族編 [#h30d367b]
vapply 関数は出力書式指定の sapply 関数だそうです．
関数 vapply(X,FUN,FUN.VALUE,...,USE.NAMES=TRUE) は sapply...
> i39 <- sapply(3:9,seq); i39 # ベクトルのリスト
[[1]]
[1] 1 2 3
[[2]]
[1] 1 2 3 4
--途中省略--
[[7]]
[1] 1 2 3 4 5 6 7 8 9
> vapply(i39,fivenum,c(0,0,0,0,0)) # 一つの結果は長さ5...
[,1] [,2] [,3] [,4] [,5] [,6] [,7]
[1,] 1.0 1.0 1 1.0 1.0 1.0 1
[2,] 1.5 1.5 2 2.0 2.5 2.5 3
[3,] 2.0 2.5 3 3.5 4.0 4.5 5
[4,] 2.5 3.5 4 5.0 5.5 6.5 7
[5,] 3.0 4.0 5 6.0 7.0 8.0 9
#comment
* こんなのあったけ？ 作表関数編 [#ded7d377]
関数 prop.table(x,margin=NULL) は分割表の各項目を margin ...
> m <- matrix(1:4,2)
> margin.table(m,1)
[1] 4 6
> margin.table(m,2)
[1] 3 7
> prop.table(m,1)
[,1] [,2]
[1,] 0.2500000 0.7500000
[2,] 0.3333333 0.6666667
> prop.table(m,2)
[,1] [,2]
[1,] 0.3333333 0.4285714
[2,] 0.6666667 0.5714286
- %%Version 1.8.1 (2003-11-21) には既にありました。%%&b...
function (x, margin)
sweep(x, margin, margin.table(x, margin), "/")
でした（今も大差ない） -- &new{2013-10-24 (木) 21:10:32};
> example(prop.table)
prp.tb> m <- matrix(1:4, 2)
prp.tb> m
[,1] [,2]
[1,] 1 3
[2,] 2 4
prp.tb> prop.table(m, 1)
[,1] [,2]
[1,] 0.2500000 0.7500000
[2,] 0.3333333 0.6666667
#comment
* こんなのあったけ？ タイマー編 [#u0512bbf]
setTimeLimit はトップレベルの各計算(つまり，コマンド行入...
> setSessionTimeLimit(10,10)
> for(i in seq(1e4)) {
cat(i); Sys.sleep(1)}
12345678910 以下にエラー Sys.sleep(1) : セッション時間が...
> {setTimeLimit(10,10);
for(i in seq(1e4))
{cat(i); Sys.sleep(1)}}
12345678910 以下にエラー Sys.sleep(1) : 時間が経過して上...
gc.time はガベージコレクションを行い，必要だった時間を返...
# 巨大なベクトルを作り，消してからガベージコレクション時...
> sapply(1:10, function(i) {x <- rnorm(1e8); rm(x); gc.t...
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,...
[1,] 1.574 1.584 1.602 1.612 1.622 1.644 1.654 1.664 1.6...
[2,] 1.026 1.072 1.122 1.168 1.218 1.268 1.318 1.364 1.4...
[3,] 1.097 1.138 1.190 1.231 1.272 1.326 1.367 1.407 1.4...
[4,] 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.0...
[5,] 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.0...
- R version 2.9.1 (2009-06-26) に既にありました.&br;なか...
- ご指摘ありがとうございます．ここで「あったけ？」と言っ...
- 全てのバージョンのRをインストールしておけばよい。あるバ...
#comment
*遅延コピー(とでもいうんだろうか?) [#id742d03]
help(tracemem)を見ていていまさらだが気づいたこと(当たり前...
> a <- 1:10
> tracemem(a) # aが指すメモリ領域がメモリ上でコピー...
[1] "<0x990d438>"
> b <- a
> b[1:3] # この段階ではa,bは同じメモリ領域を共...
[1] 1 2 3
> sum(b) # 同じく
[1] 55
> a[1] <- 10 # a が変更されたので実際のコピーが必...
tracemem[0x990d438 -> 0x990d558]: # 二度コピーされてい...
tracemem[0x990d558 -> 0x8ea59e8]:
> untracemem(a)
> a <- 1:10
> tracemem(a)
[1] "<0xa3e9cb0>"
> b <- a
> b[1] <- 10 # bが変更されればやはり実際のコピーが...
tracemem[0xa3e9cb0 -> 0xa3e9c68]:
tracemem[0xa3e9c68 -> 0xa2c5e08]:
> untracemem(a)
> tracemem(a)
[1] "<0xa3e9cb0>"
> b <- a
> rm(a) # コピー動作は無い，考えてみれば当然
#comment
*S3, S4クラス，メソッド [#p86303ff]
S言語本家の Chambers 氏が開発メンバーに加わったせいがある...
#comment
*豆知識 #line directive [#s3ee04a2]
Rでは，#記号から行末まではコメントとされ実行時には無視さ...
#line nn "ファイル名"
(#line は行頭から5文字にあるべき)はC言語のそれと同じく，p...
R news の記事の一部を引用
The #line directive
In some cases, R source code is written by a program,
not by a human being. For example, Sweave() extracts
lines of code from Sweave documents before sending
the lines to R for parsing and evaluation. To support
such preprocessors, the R 2.10.0 parser recognizes
a new directive of the form
#line nn "filename"
where nn is an integer. As with the same-named
directive in the C language, this tells the parser to
assume that the next line of source is line nn
from the given filename for the purpose of constructing
source references. The Sweave() function doesn’t
currently make use of this, but in the future, it (and
other preprocessors) could output #line directives
so that source references and syntax errors refer to
the original source location rather than to an inter-
mediate file.
The #line directive was a late addition to R
2.10.0. Support for this in Sweave() appeared in R 2.12.0
#comment
* ローマ数字の話 [#i55bfa5c]
Rには正整数をローマ数字に変換する関数があります．いまさら...
関数 as.roman(n) は正整数 n をローマ数字表記に変換する．...
# ローマ数字表記
> x <- as.roman(c(1,5,10,100)); x
[1] I V X C
# ローマ数字には0や負の数は無い
> as.roman(c(-1,0))
[1] <NA> <NA>
# クラス'roman'を持つ整数である
> str(x)
Class 'roman' int [1:6] 1 5 10 100
# ローマ数字に対する四則演算例
> as.roman(10)+1
[1] XI
> as.roman(10)-1
[1] IX
# ローマ数字では3,899までが一意的に表現可能
# 4,000以上の数は表現できない
> as.roman(3899); as.roman(3900)
[1] MMMDCCCXCIX
[1] <NA>
# ただの整数に戻る
> c(as.roman(3),as.roman(10))
[1] 3 10
- http://www.unicode.org/charts/nameslist/n_2150.html に...
- as.roman のヘルプが，http://en.wikipedia.org/w/index.ph...
#comment
ページ名:
![[RjpWiki]](image/../Rlogo.png "[RjpWiki]")
