R
与えられたベクトルxから，ランダムに要素をn個抜き出す関数． sample(x, n, replace=TRUE, prob=p) 重複したものを出す，出さないはreplaceオプションで決める．(デフォルトはFalse) 抽出する確率はprobオプションで決める．(デフォルトは等確率=1/length(x…
Quantile点の取り出し > q <- quantile(x) > q 0% 25% 50% 75% 100% 5.0 7.0 10.0 27.5 63.0 > names(q) [1] "0%" "25%" "50%" "75%" "100%" > q[1] 5.0 > q[2] 7.0
nls実行時に停止されると，その後ろの処理は一切してくれない． 途中までの結果もそのまま捨てられるので注意． 以下にエラー nls(y ~ B1 * exp(-b * x), dat, start = list(b = 0.01, B1 = y[1]/exp(0.01))) : 繰り返し数が最大値 50 を超えました Calls: n…
ある配列yの最小値を取る添字(i0)を知りたい． その添字に対応する，別の配列xの値を取得する． N <- 100 x <- rexp(N,rate=1) y <- runif(N) i0 <- which(y==min(y)) x[i0]
頻度が指数分布する乱数の発生． rm(list=ls()) set.seed(1) N <- 1000 #発生させる個数 La <- 1/140 #指数分布のパラメータλ x <- round(rexp(N, rate=La),0) #Rの組み込み関数を使う x2 <- round(-(log(1-runif(N)))/La,0) #一様乱数から変換する mean(x) …
Rのトップに表示している図をPDFでそのまま保存する． dev.copy2eps(file="hoge.eps")だそうです． PDFやらEPSに出力，保存しようと思ったら，なぜか ずれたり，天地逆になったり，対数にしたら保存されなかったりと， 全然できませんでしたので，この方法で…
Rで微分法定式を解くには、odesolveというパッケージが必要らしい。 ということで、Rの画面で install.packages("odesolve",depend=T) library(odesolve)とする。depend=Tは依存性のある他のやつらもインストールしますよ、ということらしい。 sample関数を…
- バッヂモードで外部引数を数値として受け取る 通常の引数は文字列として扱われるので、as.numericをつけて受け取る。 x<-read.table(commandArgs()[5]) alpha <- as.numeric(commandArgs()[6]) y<-x[,1]*(-alpha) cat(y);cat("\n")R --vanilla --slave --a…
べき分布のべき指数推定、Hill estimatorによる場合の計算 Inputdataの中の数値が対象。 x<-read.table("Inputdata") x1<-x[,1] xmin<-min(x1) alpha=1+length(x1)*(sum(log(x1/xmin)))^-1引用元Title:A Simple General Approach to Inference About the Tai…
分布の値を返す、pなんとか。 例えば、最小値が0,最大値が1の一様分布において、0.5より小さい値が出る確率 punif(0.5,0,1)0.5になる。 カイ二乗分布では。 例えば、自由度が2、カイ二乗値が9を与えた場合、カイ二乗分布における、 0から9の値のところまでの…
平均が5のポアソン乱数を700個発生させる。 po<-rpois(700,5)
Rでヒストグラムを出力して、外部ファイルにも出力。 x<-read.table("Inputfile.csv") y<-hist(x[,1],breaks=seq(min(x),max(x)+0.01,0.01)) z<-paste(y$breaks,y$count) write.table(z,row.names=F,col.names=F,file="Outputfile.csv")このままだとダブルコ…
変数の名前を確認する「names」 x<-read.table("TMP") y<-acf(x) names(y) [1] "acf" "type" "n.used" "lag" "series" "snames" #結果「こんな変数の名前使ってるよ」 y$n.used #「n.used」は使ったデータポイント数 [1] 137184 nrow(x) [1] 137184 #元のデ…
library(igraph) g<-read.graph(commandArgs()[5], format="edgelist",directed=FALSE) cs<-clusters(g) csh<-hist(cs$csize,breaks=seq(min(cs$csize),max(cs$csize),1)) S<-0 for(i in 1:max(csh$breaks)){ cat(i,csh$counts[i],"\n") if(!is.na(csh$count…
Rによるネットワーク描画、解析ソフトライブラリ。 igraph Rへのインストールと、ライブラリ読み込み。参考 install.packages("igraph") library(igraph) ランダムグラフ g <- erdos.renyi.game(100, 3/100) plot(g, layout=layout.fruchterman.reingold,ver…
hist(x[,1], breaks=seq(-0.5, 1, 0.01),col="blue")
異時刻相関 cross correlation function。 xとyというデータがあったとして、その異時刻相関。 ccf(x,y) option : lag.max=10など。
Output result of acf function. x s write.table(data.frame(s$lag,s$acf),row.names=F,col.names=F,file="output.txt")
計算の実行時間を調べる。 100万件の計算をするときには、ラフにどのくらい計算時間がかかるのか 知っておくことも大事。1日は86400秒なので、1件1秒の計算だと、10日以上かかるし、 1件10秒の計算だと、3か月以上かかってしまうことに・・・。 timer<-Sys.t…
use R from unix terminal Simple example; echo "x Test.tmp echo "y > Test.tmp echo "cor(x,y,method=\"spearman\")" >> Test.tmp R --vanilla --slave (or R --vanilla --slave More tricky way; R --vanilla --slave --args 1 2 Corr1.R as follows; x …
x列y行の時系列データからデンドログラムを作成する。 x方向がデータの種類、y方向が時刻とするファイル。 x y z plot(z)
