R

線形回帰の結果から回帰係数とR2の値を取り出す(プログラムなどで使用)

s<-summary(lm(y~x1+x2+1))
b1<-s$coefficients[1, "Estimate"]
b2<-s$coefficients[2, "Estimate"]
r2<-s$r.squared

(行頭に半角スペースを入れておくと、こんな感じに表示されますよ。「□→」のボタンを押すと行頭にスペースがはいります。by津田)

エクセルのデータファイル(*.xls)の 読み込み

> library(RODBC)
>fm<-read.table("file.xls",header=T)


図をIllsutrator用epsファイルで保存

> hist(k,breaks=seq(0,120, by=5))
> dev.copy2eps(file="testfig.eps", family="Helvetica")


度数分布

breaks <- seq(-2, 4, by=2)               
result <- table(cut(x, breaks)) 


ディレクトリ内のファイル一覧を取得

files = list.files()

エラーバーや信頼区間をグラフに出力(RjpWiki?より)

x<-c(1,2,3,4,5) #横軸の値
y<-c(1.1, 2.3, 3.0, 3.9, 5.1) #平均値や中央値など
ucl<-c(1.3, 2.4, 3.5, 4.1, 5.3) #平均値+標準偏差の値
lcl<-c(.9, 1.8, 2.7, 3.8, 5.0) #平均値-標準偏差の値
plot(x,y, ylim=range(c(lcl,ucl))) #平均値の点をプロット
arrows(x,ucl,x,lcl,length=.05,angle=90,code=3) #(1)横棒付きのエラーバーを書く
#(1)の代わりに(2)を用いても良い、
segments(x,ucl,x,lcl) #(2) 縦棒のみのエラーバー

error_bar.png

比較演算に注意

計算機内では16進数なので、そのへんの近似の都合上こんなことが起こる

> (a = 0.3 - 0.1)
[1] 0.2
> (b = 0.2)
[1] 0.2
> a == b
[1] FALSE

[Mac] mi.appで選択中の文字列をR.appで実行するAppleScript?

下記のスクリプトを/Applications/AppleScript?/Script Editor.appにコピペ。
アプリケーションバンドルとして保存(例えばmi2R.app)。
/Applications以下の好きな場所に置く。
mi.appのモード設定でツールやキーバインドの挿入文字列に <<<LAUNCH(mi2R)>>> を入れる。
command + R とかに割り当てて呼ぶと幸せ。

try
	tell application "mi"
		tell document 1
			set CommandLines to selection
		end tell
	end tell
	
	using terms from application "R"
		tell application "R"
			activate
			with timeout of 3600 seconds
				cmd CommandLines
			end timeout
		end tell
	end using terms from
end try

3次元グラフ

test2.png

	#テストデータ生成
	num<-5 #num*numのデータ
	herbivore<-NULL
	for(tmp1 in 1:num){
		for(tmp in 1:num) herbivore<-c(herbivore,tmp)
	}
	plant<-NULL
	for(tmp in 1:num){
		for(tmp1 in 1:num) plant<-c(plant,tmp)
	}
	species<-c(1:(num*num)) 
	
	# rglパッケージを使う
	# http://cse.naro.affrc.go.jp/takezawa/r-tips/r/57.html
	library(rgl)
	open3d()
	rgl.bg(color=c("white", "black"))
	rgl.light(theta = 5, phi = 5)
	
	x<-c(1:num) #x軸
	y<-c(1:num) #y軸
	f<-function(x,y){ r<-species }
	z<-outer(x,y,f) #z軸
	
	#プロット
	plot3d(herbivore,plant,species)
	
	#面
	#テストデータを色分けする
	col2<- heat.colors(10)
	colo<-col2[species]
	terrain3d(x,y,z, color=colo)
	
	#各プロットを線でつなぐ
	for(tm in 0:(num-1)){
		lines3d(herbivore[c((num*tm+1):(num*(tm+1)))],plant[c((num*tm+1):(num*(tm+1)))],species[c((num*tm+1):(num*(tm+1)))], lwd=8)
		lines3d(herbivore[num*c(0:(num-1))+tm+1],plant[num*c(0:(num-1))+tm+1],species[num*c(0:(num-1))+tm+1], lwd=8)
	}
	
	rgl.viewpoint(0, -90) #図形の角度
	#for(i in 1:360) rgl.viewpoint(i,i/4)   # 図形の回転(私はうまく使えない)
	############ (1)pngに保存する場合 ############
	# まず図にしたい角度に調整します
	pngname<-paste("test.png",sep="")
	rgl.snapshot(pngname)
	rgl.quit()
	############ (2)画像を連続させてpngに保存する場合############ 
	M <- par3d("userMatrix") 
	#M <- par3dinterp( userMatrix=list(M,rotate3d(M, pi/2, 0, 0, 1) ), duration=2 ,movie="test1",dir=".") #(使い方を忘れました)
	movie3d( par3dinterp( userMatrix=list(M,rotate3d(M, pi/2, 0, 0, 1) )), duration=2 ,movie="test1",dir=".") #この場合、21個画像ファイルが作られるので注意。ただし何かエラーが出るので、どなたか直してもらえますか?

高速化Tips

計算にかかる時間を測定

> system.time( 測定したいコマンド )

演算のベクトル化:forを避ける

Rではvectorやmatrixの構造をそのままに、各要素に対して一気に演算できるようになっている。むき出しのforを自分で書かずにこれを利用したほうが高速なので、常に心がける。ifelseはちょっと分かりにくいけど感動的。

> (x = c(1:10))
[1]  1  2  3  4  5  6  7  8  9 10

> x * 2    #各要素かける2
[1]  2  4  6  8 10 12 14 16 18 20

> ifelse(x>5, x, 0)    #各要素>5をチェックし、Tならそのまま、Fなら0
[1]  0  0  0  0  0  6  7  8  9 10

data.frameは遅い

計算量が多い場合はmatrixに変換してからやったほうが断然速い。ただしmatrixは1種類のデータ型しか入れられないので、カテゴリ変数(フィールドの名前とか)などもすべてnumericで表せるように変換する必要がある。

JAGS

Just Another Gibbs Samplerの頭文字。MCMCのアルゴリズムのひとつであるギブズサンプリングをやるソフト。データをRで整えてJAGSに渡し、出て来た結果をRで受け取って解析する。

JAGS本体のインストール

  • JAGSのサイトからOSに合ったインストールファイルをダウンロードする。
    • Windowsなら *.exe
    • Mac(OS10.4以上)なら *.dmg
  • ダウンロードしたファイルを実行し、「はい」とか「次へ」で普通にインストールする。

Rから使えるようにする(JAGS本体をインストールしてから)

  • Rを立ち上げて、rjagsというパッケージを読み込む。これにはRからJAGSを動かすための橋渡しとなる関数が入っている。
    > library(rjags)
  • もし、rjagsなんて見つからないと怒られたら以下のコマンドでインストール。
    > install.packages("rjags")
  • それでもダメな場合、codaやlatticeというパッケージがインストールされてない可能性がある。
    あと、Rは最新版(2.7以上)じゃないとダメらしい。
    > install.packages("coda")

添付ファイル: filetest2.png 840件 [詳細] filetest.png 377件 [詳細] filegr.jpg 645件 [詳細] filegr.pdf 690件 [詳細] fileerror_bar.png 1578件 [詳細]

トップ   編集 凍結 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS
Last-modified: 2011-03-08 (火) 18:58:53 (2450d)