Rな予測 で以下のようにいっているが，predict には，level という引数がある。これで level=0.90 とすれば，90%推定区間が計算される。 予測の話」にあったSEやSD（的なもの）を回帰分析においても計算してやって、信頼区間とか予測区間とかを出しているわ…

ifの条件式に関して > a <- 1 ; if ( (0 < a) < 2) print(a)[1] 1> a <- 1 ; if ( (0 > a) < 2) print(a)[1] 1> a <- 1 ; if ( (0 < a) > 2) print(a)> a <- 1 ; if ( (0 > a) > 2) print(a)上2つの式はa2で評価されるようで表示されません。そもそも他の言…

CX's 'reality distortion' pie-graph (フジTV のステマ円グラフ) 元のプログラムは，無駄なことをしているので，プログラムが長く煩雑で，何をやっているのかわかりにくくなっている。 for (i in 1:length(num)) { x <- c(0.5 + dx, 0.5 + 0.5 * sin(rad[i]…

Project Euler in R: Problem 25 の中の aatrujilloJan 25, 2012 03:02 PM なぜ list なんかを使うかなあ。途中の項も不要なんだから，以下のようにすれば，計算時間は半分になる。 a <- as.bigz(1) b <- as.bigz(1) i <- 2 repeat { i <- i+1 c <- a+b if (…

npパッケージのnpcmstestでbwsで具体的にbandwidthを指定すると「仮引数 ”bws” が複数の実引数にマッチしました」というエラーがでてしまいます † ななみ? (2012-01-21 (土) 21:55:03) はじめまして。npパッケージでモデルの検定をしています。 model <- lm(…

この一つ前の記事で。 二重の for は効率が悪い。outer を使うのがよい。 system.time({x <- seq(-0.5, 3, length = 1000)y <- seq(-2.5, 2.5, length = 1000)z <- matrix(0, 1000, 1000)for (i in 1:1000) { for (j in 1:1000) { z[i, j] <- (x[i] ^ 2 + y[…

陰関数のプロットについて のコメントはちょっと外れているような気がする。 追記 2022/07/24以下を，R で書いてみたOctave で陰関数のプロット ちょっと，裏技的であるが，3次元関数として計算し，任意の等高線を描くというようにすることができることもあ…

Google Map 上に描く fn <- "new"ne <- decodeGPS(fn) genMap(fn, ne, size = c(640, 640)) 使う関数は以下の2つ genMap <- function(fn, ne, col = "red", lwd = 3, maptype = "roadmap", size = NULL, zoom = NULL, format = "png") { library(RgoogleMaps…

vander <- function(x, n = length(x)) { stopifnot(as.integer(n) == n && n > 0) z <- matrix(rep(x, n), length(x)) for (i in 1:n) { z[, i] <- z[, i] ^ (i - 1) } return(z[, n:1, drop = FALSE])} のように定義する。実行例は， > vander(1:5, 1) [,1…

http://d.hatena.ne.jp/foo22222/20120115Nt <- 1000df1 <- 1set.seed(101)ts1 <- rchisq(Nt, df1)set.seed(101) # set.seed()で生成される乱数を同じものにする。ts2 <- rchisq(Nt, df1)t1 <- sample(ts1, 500, replace = FALSE)t2 <- sample(ts2, 500, rep…

このブログで実行速度を比較していることがあるが，実行速度至上主義でないことはご理解頂きたい。 100 秒かかるものが1秒で終わっても，なんにもうれしくない。100 秒かかるものが 1.1 秒で終わろうが， 0.9 秒で終わろうが，問題ではない。 100 時間かかる…

ぱらぱらめくる『Rによる計算機統計学』と『生物学のための計算機統計学』の「多変量正規乱数」なんでまあ，ぱらぱらめくるのかな？それはともかく，MASSパッケージの mvrnorm 関数を使うのがよろしかろうと。

return のときに次元を決めるのを前の記事の最後に書いたけど，やはり宣言する場所で決めるのが普通なので，そちらを。しかし，いずれにしても，添え字は一個。 > library(inline)> src <- '+ int n1 = as<int>(N1);+ int n2 = as<int>(N2);+ int n3 = as<int>(N3);+ Rcpp::</int></int></int>…

array とはいっても，C 言語では単なる一次元ベクトルを三個以上の添え字でアクセスするだけ。 C 言語と R では，行と列の優先順序が正反対という点のみ。 「Rcpp::Arrayがないので」 では，本当の意味での array を操作していることにはなっていない。 つま…

なぜか，これらの整形プログラムは，"/" などの前後には空白を「置けない」 > library(formatR)> tidy.source(text = "result <- a+b-c*d/e^f%%g%/%h%in%i%o%j", width.cutoff = 150)result <- a + b - c * d/e^f%%g%/%h %in% i %o% j

n<-3ds<-runif(n)us<-sample(1:10,n)Niter<-100000dt<-0.01m<-matrix(0,Niter,n)p<-rep(0,Niter)m[1,]<-runif(n)p[1]<-0for(i in 1:n){ p[1]<-p[1]+us[i]*m[1,i]}#pが1からスタートするように補正m[1,]<-m[1,]/p[1]p[1]<-0for(i in 1:n){ p[1]<-p[1]+us[i]*m…

コンピュータ上の実数は近似値である。ただし，2 のべき乗数の和で表される数は正確な値を保持している。整数は全て2のべき乗数の和であるから正確な値を保持している。 > options(digits=20)> x <- 0.001 # 0.001 は 10 進数では正確な値であるが，2 進数で…

たぶん，最終的に用いているのは皆同じプログラムなんだろうけど，Sweave を使う例も挙げておこう。 <<>>=test<-function(x,y,z) {n<-length(x) # 日本語のコメントresult<-numeric(n)for (i in 1:n) {result[i]<-(x[i]+y[i])*z[i]}return(result)}@を含む *…

formatR パッケージの tidy.souce 関数を使う。この方法では，コメントはちゃんと残る。 元のプログラム（encoding: utf-8） test<-function(x,y,z) {n<-length(x) # 日本語のコメントresult<-numeric(n)for (i in 1:n) {result[i]<-(x[i]+y[i])*z[i]}return…

source("プログラムファイル.R", keep.source=FALSE)dump("関数名", file="新プログラムファイル.R") 難点 元のプログラムについている注釈が全部なくなる 元のプログラム testsource.R test<-function(x,y,z){n<-length(x)result<-numeric(n)for(i in 1:n){…

定数の引用 その言語処理システムで定義されているなら，それを使う。 R では，円周率は pi で引用できる。 もし定義されていない場合には，pi <- atan(1)*4 で正確に定義できる。 pi <- 3.14159 などと，いい加減に定義しないこと。 ネイピア定数は定義され…

Rプログラミングのレクチャーなどするお方は，標準的なプログラミング書法に従わないと，受講生は皆，悪しきお手本に従って，悪しきプログラミングをする ようになるという実例があるわけだから，早々に改める方がよいと思うけどなあ。なぜ，悪しきプログラ…