Site admin: motohiro PukiWiki 1.5.1 © 2001-2016 PukiWiki Development Team. Powered by PHP 7.4.33
HTML convert time: 0.041 sec.
個人メモを整理中です．不適切な記述が多々あるかと思います．お気づきの際は，ishida-m.at.ias.tokushima-u.ac.jp. までご連絡ください． 個人メモを整理中です．不適切な記述が多々あるかと思います．お気づきの際は，ishida-m(plase put "atmark" here)ias.tokushima-u.ac.jp までご連絡ください． #contents * scatterplot.matrix 散布図行列の作成 2007 01 18 [#e6d68c3f] plot だと対角要素は変数名 対角要素をヒストグラムにしたければ car パッケージの scatterplot.matrix 例えば [[エヴェリット, RとS-PLUSによる多変量解析:http://www.amazon.co.jp/R%E3%81%A8S-PLUS%E3%81%AB%E3%82%88%E3%82%8B%E5%A4%9A%E5%A4%89%E9%87%8F%E8%A7%A3%E6%9E%90-B-%E3%82%A8%E3%83%B4%E3%82%A7%E3%83%AA%E3%83%83%E3%83%88/dp/4431713123/]] usair.dat<-source( "/everitt/Data/chap3usair.dat")$value # install.packages("car") library(car) scatterplot.matrix(usair.dat[,-1],diagonal= "histo" ) plot(usair.dat) * [[エヴェリット, RとS-PLUSによる多変量解析:http://www.amazon.co.jp/R%E3%81%A8S-PLUS%E3%81%AB%E3%82%88%E3%82%8B%E5%A4%9A%E5%A4%89%E9%87%8F%E8%A7%A3%E6%9E%90-B-%E3%82%A8%E3%83%B4%E3%82%A7%E3%83%AA%E3%83%83%E3%83%88/dp/4431713123/]] の累積分布密度グラフ 2007 04 05 [#r4379818] pdf("p9.pdf") par(mar=c(2,2,2,2)) plot(seq(-10, 10, 0.01), seq(0, 1, 0.0005), type = "n", axes =F, xlab = "", ylab = "", main = "Cumulative distribution function") lines(x[order(x)], xp, xlim=c(-10, 10), lwd = 3 ) lines(y[order(y)], yp, xlim=c(-10, 10), lwd = 3 ) lines(c(-9,9), c(1,1)) lines(c(-9,-9), c(1,0)) lines(c(-9,9), c(0,0)) lines(c(9,9), c(1,0)) arrows(-9, 0.8, x.p0.8, 0.8, length = 0.2, lty = 2, lwd = 2) arrows(-9, 0.8, y.p0.8, 0.8, length = 0.2, lty = 1, lwd = 2) arrows( x.p0.8, 0.8, x.p0.8, 0, length = 0.2, lty = 2, lwd = 2) arrows(y.p0.8, 0.8, y.p0.8, 0, length = 0.2, lty = 1, lwd = 2) text(3.8, 0.85, expression(paste("c.d.f ", F[x], "( )")), cex = 1.1) text(-2.0, 0.85, expression( paste("c.d.f ", F[y], "( )")), cex = 1.1) text(-10, 0.8, "p", cex = 1.1) text(x.p0.8+0.5, -.02, expression(q [x] (p)), cex = 1.1) text(y.p0.8, -.02, expression(q [y] (p)), cex = 1.1) arrows(x.p0.2, 0.0, x.p0.2, 0.2, length = 0.2, lty = 1, lwd = 2) arrows(x.p0.2, 0.0, x.p0.2, pnorm( x.p0.2, mean = -1.5, sd = 2, length = 0.2, lty = 2, lwd = 2) arrows(x.p0.2, 0.2, -9, 0.2, length = 0.2, lty = 1, lwd = 2) arrows(x.p0.2, pnorm( x.p0.2, mean = -1.5, sd = 2) , -9, pnorm( x.p0.2, mean = -1.5, sd = 2), length = 0.2, lty = 2, lwd = 2) text(x.p0.2, -.02, "q", cex = 1.1) text(-10, pnorm( x.p0.2, mean = -1.5, sd = 2), expression(p [x] (q)), cex = 1.1) text(-10, 0.2, expression(p [y] (q)), cex = 1.1) text(-10,0, "0", cex = 1.2) text(-10,1, "1", cex = 1.2) ## text(-10, 1, pch = 64, cex = 1.2)# for test dev.off() * R の日本語化 2007 04 15 [#n8831b6f] [[Rjp Wiki:http://www.okada.jp.org/RWiki/?%C6%FC%CB%DC%B8%EC%B2%BD%B7%C7%BC%A8%C8%C4]] を参照 * 正規分布の理論曲線の引きかた 2007 05 23 [#o911dcb4] x<- rnorm(1000, mean = 50, sd = 10) riron<- dnorm(sort(x), mean = 50, sd = 10) * sum(x) / sd(x) hist(x) curve(sort(x), riron) 他に[[リゲス,Rの基礎とプログラミング技法:http://www.amazon.co.jp/gp/product/4431712186/]]や [[Crawley Statistics Using R:http://www.amazon.co.jp/Statistics-Introduction-Michael-J-Crawley/dp/0470022981/]] p. 56 - 57. * 画像内の日本語フォント Illustrator 2006 06 06 [#t5a66a15] [[Rjp Wiki:http://www.okada.jp.org/RWiki/?%BD%E9%B5%E9%A3%D1%A1%F5%A3%C1%20%A5%A2%A1%BC%A5%AB%A5%A4%A5%D6(6)]] より * pdf の作成 2007 06 06 [#z6170668] options(X11fonts = c("-alias-mincho-%s-%s-*-*-%d-*-*-*-*-*-*-*", "-adobe-symbol-*-*-*-*-%d-*-*-*-*-*-*-*")) ps.options(family= "Japan1Ryumin") # これを指定しない場合， # 例えば Illustrator で文字化けする(開くことはできる) pdf("/home/ishida/test.pdf") plot(sin,-pi,pi,main="正弦") dev.off() * フォントの埋め込み embedFonts 2006 06 06 [#nc4d8299] options(X11fonts = c("-alias-mincho-%s-%s-*-*-%d-*-*-*-*-*-*-*", "-adobe-symbol-*-*-*-*-%d-*-*-*-*-*-*-*")) ps.options(family= "Japan1Ryumin") plot(sin,-pi,pi,main="正弦") dev.copy2eps(file = "/home/ishida/test.eps") このファイルは，例えばIllustrator では開けない． embedFonts(file = "/home/ishida/test.eps", outfile = "/home/ishida/test2.eps") このファイルは，Illustrator では開けるが，文字はビットマップ化されて美しくない． ちなみに options(X11fonts = c("-alias-mincho-%s-%s-*-*-%d-*-*-*-*-*-*-*", "-adobe-symbol-*-*-*-*-%d-*-*-*-*-*-*-*")) # ps.options(family= "Japan1Ryumin") plot(sin,-pi,pi,main="正弦") dev.copy2eps(file="/home/ishida/test3.eps") では日本語がそのまま入るが，例えばIllustratorでは開けない． ps.options(family= "Japan1Ryumin") を実行すると，日本語部分はコード化されて記号となるが，Illustrator では開けない． * Ghostscript 用の環境変数 R_GSCMD 2007 06 12 [#v034c83d] system("which gs") /usr/local/TeX/bin/gs R_GSCMD=${R_GSCMD='/usr/local/TeX/bin/gs'} Sys.putenv("R_GSCMD"="/usr/local/TeX/bin/gs") * CarbonEL for Mac Carbon Emacs 2007 06 12 [#e5274d86] How to install: install.packages("CarbonEL",,"http://rforge.net/") to load: library(CarbonEL) Make sure CarbonEL package gets loaded by ESS such that Quartz is always responsive. This also works for R run from a shell such that Terminal or xterm. Cheers, Simon PS: The above doesn't fix other problems with the old Quartz device (R.app uses a different, newer Quartz device! E.g. quartz.save(..) works only in R.app). It is just a hack for people that don't mind using the old Quartz device ;). http://login.stat.ucla.edu/pipermail/statcompute/2003-November/000629.html * 文字ベクトルの一部を取り除く 2007 06 14 [#cf95914c] writers[writers!="Kajii"] writers[c(writers!= "Kajii" & writers != "Koda")] writers[writers!= "Kajii" | writers != "Koda"] * プロットを重ねる par new = T 2007 06 15 [#e1d7192b] par(new = T) * テーブルの区間ラベルの操作．行列名を取り出すなど 2007 06 16 [#ld5529b9] as.numeric(dimnames(table(buntyo))$buntyo) x<- c(1,1,1,1,2, 5,5,5, 6,6) table(cut(x, breaks = 0:max(x))) table(cut(x, breaks = 0:max(x), lab = F)) # 頻度表を作る as.numeric(labels(table(cut(x, breaks = 0:max(x), lab = F)))[[1]]) table(cut(x, breaks = 0:max(x), lab = 0:max(x))) # 頻度表を作る * file open close ファイルを一行ずつ読み込む 2007 07 10 [#fa05e4f5] con<- file("test.txt") open(con) # 最初の 2 行は con<- file("file", "r") でも良い while(1){ x<- readLines(con, n=1) if(length(x)> 0){ print(x) }else{ break; } } close(con) * UTF-8 で R で Rstem を実行する方法 2007 07 11 [#sad16559] 始めに [ishida@amd64 tmp]$ export LANG=ja_JP.UTF-8 ちなみにこの時，コンソールのメニュー設定も変更すること コンソールから sessionInfo() R version 2.5.1 (2007-06-27) x86_64-unknown-linux-gnu locale: LC_CTYPE=ja_JP.UTF-8;LC_NUMERIC=C; LC_TIME=ja_JP.UTF-8;LC_COLLATE=ja_JP.UTF-8; LC_MONETARY=ja_JP.UTF-8;LC_MESSAGES=ja_JP.UTF-8; LC_PAPER=ja_JP.UTF-8;LC_NAME=C; LC_ADDRESS=C;LC_TELEPHONE=C; LC_MEASUREMENT=ja_JP.UTF-8;LC_IDENTIFICATION=C kafka <- readLines("tmp/kafkaUTF8.txt") library(Rstem) kafka.vec<- unlist(strsplit(kafka[5], "[[:space:]]")) kafka.vec [1] "" "Zu" "diesem" [16] "Gesetz." Warning message: Currently, only 'english' is tested. You will need support for UTF characters. in: wordStem(kafka.vec, language = "german") * 正規表現 使用例 2007 07 11 [#uccbaf1a] 日本語文字コードチェック perl での日本語処理は以下を参照いたしました． http://www.din.or.jp/~ohzaki/perl.htm# # or ~/progSource/Perl/perl.html regex.html # Perl で間違ってマッチしてしまう例 $str = 'これはテストです'; $pattern = '好'; if ($str =~ /$pattern/) { print "マッチした\n"; } 日本語文字コードチェック str.jp<- "これはテストです" length(str.jp) gsub("です","--", str.jp) # "これはテスト--" gsub("^(\\w)","J", str.jp)#"Jれはテストです" regexpr("です",str.jp) # [1] 7 # attr(,"match.length") # [1] 2 gregexpr("です",str.jp) # [[1]] # [1] 7 # attr(,"match.length") # [1] 2 gsub("好","--", str.jp)#"これはテストです" regexpr("好",str.jp) # [1] -1 # attr(,"match.length") #[1] -1 gregexpr("好",str.jp) # [[1]] # [1] -1 # attr(,"match.length") # [1] -1 茶筅の解析結果を取り込んで読み込むたびにスペースで区切る con open(con) while(1){ if(length(x)> 0){ for(j in 1:length(z)){ print(z[j]) } }else{ break; } } close(con) 読み込むたびにスペースで区切る con<- file("chasen.txt") x<- readLines(con, n=1) close(con) x [1] "ここ\tココ\tここ\t名詞-代名詞-一般\t\t" y <- gsub("^(.+)[[:space:]]+(.+) [[:space:]]+(.+)[[:space:]]+(.+)[[:space:]]+", "\\1_\\2_\\3_\\4", x) y [1] "ここ_ココ_ここ_名詞-代名詞-一般\t" # y<- gsub("^(.*)[[:space:]]+(.*) [[:space:]]+(.*)[[:space:]]+(.*)[[:space:]]+", "\\1_\\2_\\3_\\4", x) # y<- gsub("([[:alpha:]]*) [[:space:]]+([[:alpha:]]*) [[:space:]]+([[:alpha:]]*) [[:space:]]+([[:alpha:]]*)[[:space:]]+", "\\1_\\2_\\3_\\4", x) # y<- gsub("(.*)\\t(.*)\\t(.*)\\t(.*)", "\\1_\\2_\\3_\\4", x, perl = T) 文字列の長さ mixed.text<- "花a子" nchar(mixed.text ) # 5 nchar(mixed.text, type = "chars") # 3 sub('[[:alnum:]]', '', mixed.text) # "a子" gsub('[[:alnum:]]', '', mixed.text) # "" gsub('[[:lower:]]', '', mixed.text) # "花子" gsub('[[:upper:]]', '', mixed.text) # "花a子" * 箱鬚図のノッチについて 2007 07 30 [#m1b1f908] ''以下からの引用です！（石田）'' http://d.hatena.ne.jp/syou6162/20070702/1183349416 より ボックスプロットは非常に分かりやすいプロットであるため、 使い易いツールである。 しかし、medianが統計的に有意に違うかどうかなどを 顧客に説明するには足りない部分がある。 そこでボックスプロットを用いてmedianが 統計的に有意に異なるかを調べるためのツールが、 このノッチの付いたボックスプロットである。 ノッチのついたボックスプロットは以下のようにして 使うことができる。 boxplot(lograifall.control,lograinfall.seeded,notch=TRUE) 図にてへこんだ部分がノッチである。 これがかぶっていなければ統計的に有意にmedianが 違いということが言える。 また、ノッチのついたボックスプロットはふたつの違いが 有意かどうかしか検定できないことに注意した。 みっつ以上のノッチについては使えないということである (おそらくt検定だけではなく、F検定が必要な理由と同じと思う)。 * モデル診断 model.diagnostics (cook, residuals) プロットの見方 2007 07 31 [#s7b28df7] [[行動計量学会チュートリアル:http://150.59.18.68]]より 最後に紹介したプロットは，モデルの適合度を確認する補助となるプロットを作成した例である． 詳細は Faraway\citep[pp.14 -- 17]{faraway05}, Crawley\citep[pp.357 -- 359]{crawley07} Wood, p.33 参照されたい． ここでは左上のプロットは，残差と適合値を対応させたもので， データの変動になんらかの規則性がなければ， この散布図上で各ポイントはランダムに散っているはずである． この図の場合，体重の変化には何らかのシステマティックな要因が働いている， あるいは分散に変動があると思われる．左下の図は残渣を標準化し (平均 0，標準偏差 1 とする)， その絶対値の平方根の取った図である(絶対値のままだと skew が生じるため)． 後者の図では点が三角形の形になる場合，heteroscedascity が疑われる． この種の図が得られた場合は，例えば変数の二次の項を導入するなどの変換を検討することになる． 右上の図は QQ プロットと呼ばれ，残差が正規分布に従っているかをチェックするのに使われる． 最後の図はクックの距離で，モデルの当てはめに影響のあるはずれ値を検出するのに利用される． ここで点線がクックの距離を表し，0.5 を越えるデータは影響が少なくない， 1 を越える場合はかなり影響があると見なす． それぞれのグラフで，問題となるデータ点にはデータ番号が振られる．