裏 RjpWiki

文字通り,RjpWiki の裏を行きます
R プログラム コンピュータ・サイエンス 統計学

ダメ出し:オンラインヘルプをよく読もう

2012年01月31日 | ブログラミング

Rな予測 で以下のようにいっているが,predict には,level という引数がある。これで level=0.90 とすれば,90%推定区間が計算される。

予測の話」にあったSEやSD(的なもの)を回帰分析においても計算してやって、信頼区間とか予測区間とかを出しているわけです。

SE、SDの使い分けは

A<-predict(model,new,se.fit=T,interval="confidence")    #推定平均の95%推定区間付き
B<-predict(model,new,se.fit=T,interval="prediction")    #推定データの95%推定区間付き

↑に注目

 predict関数とはその名の通り予測値を返してくれる関数です。se.fit=Tにすると誤差の範囲も出してくれます。ここで、interval="confidence"とinterval="prediction" とに分かれているのですが、これで「データの95%誤差範囲」か「平均の95%誤差範囲」を出してくれるわけです。confidenceのほうは平均の誤差(信頼区間)で、 predictionの方はデータの範囲(予測区間)を表しています。

 もっと詳しい説明です。

 上記のやりかたはとても便利なのですが、95%以外の範囲を指定してやることはできません。そこで、もっと柔軟に計算するためには interval=~~を使わずに、直接SEとかSD(的なもの)を計算することになります。

コメント

基本を勉強しよう

2012年01月31日 | ブログラミング

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で評価されるようで表示されません。

そもそも他の言語のように

> a <- 1 ; if (0 < a < 2) print(a)
エラー: 予想外の '<' です ( " if (0 < a <" の)

のように書いてエラーとなり、()をつけて&として解釈されないかなーと
上の式をなんとなく書いてみたんですが、
Rでどう評価されてるのかよくわかりませんでした。(素直に0<a & a<2って書けよですね)

よくわからないときは,部分部分を評価してみる。

a <- 1 のとき,(0 < a) は TRUE,(0 > a) は FALSE
R では TRUE/FALSE は数値としては 1/0 として認識される
よって,(0 < a) < 2 は TRUE < 2 そして 1 < 2 なのだから TRUE になる
同様に,(0 > a) < 2 は FALSE < 2 そして 0 < 2 なのだから,この場合も TRUE になる
(0 < a) > 2 は TRUE > 2 そして 1 > 2 なのだから FALSE
(0 > a) > 2 は FALSE > 2 そして 0 > 2 なのだから FALSE

そもそも「他の言語のように (0 < a < 2) のように書いてエラーとなり、()をつけて&として解釈されないかなーと」いうけど,そんなことを許す言語ってあるか?ないと思うぞ。

コメント (1)

ダメ出し:計算の無駄,プログラムの無駄をなくする

2012年01月31日 | 裏 RjpWiki

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]), 0.5 + 0.5 * sin(seq(rad[i], rad[i+1], length.out=div)), 0.5 + 0.5 * sin(rad[i+1]))
  y <- c(0.5 + dy, 0.5 + 0.5 * cos(rad[i]), 0.5 + 0.5 * cos(seq(rad[i], rad[i+1], length.out=div)), 0.5 + 0.5 * cos(rad[i+1]))
  grid.polygon(x=x, y=y, gp=gpar(fill=colors[i]))
  grid.text(引数省略)
}

grid.polygon を描くところは,以下のようにすればよい。
0.5 * sin(rad[i]),0.5 * sin(rad[i+1]) などを付け加える必要はない。
theta も一度計算しておけばよい。

for (i in 1:length(num)) {
  theta <- seq(rad[i], rad[i+1], length.out = div)
  x <- c(dx, 0.5 * sin(theta))
  y <- c(dy, 0.5 * cos(theta))
  grid.polygon(x = 0.5 + x, y = 0.5 + y, gp = gpar(fill = colors[i]))
  grid.text(引数省略)
}

コメント

ダメ出し:Project Euler in R: Problem 25

2012年01月27日 | ブログラミング

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 (nchar(as.character(c)) >= 1000) break
        a <- b
        b <- c
    }
    print(i)

ところで,元のプログラムは,

fibon1 = list()
fibon1[[1]] = fibon1[[2]] = as.bigz(1)
digits = 1000
i=3
n=0

while(n < digits){
fibon1[[i]] = fibon1[[i-1]]+fibon1[[i-2]]
n = nchar(as.character(fibon1[[i]]))
i=i+1
}
i-1

のように,答えの表示で i-1 としている。
これは,カウンターの位置が不適切だからそうしないと答えが合わないからだが,以下のようにするのが自然。

fibon1 = list()
fibon1[[1]] = fibon1[[2]] = as.bigz(1)
digits = 1000
i=2
n=0

while(n < digits){
i=i+1
fibon1[[i]] = fibon1[[i-1]]+fibon1[[i-2]]
n = nchar(as.character(fibon1[[i]]))
}
i

プログラム書法の一つ:カウントのしそこないに注意

while ループの前に n=0 などと無駄な設定が必要なのは,while ループを使っているから。新しい項を計算して,その項の桁数を見ればよいのだから,repeat を使えば n=0 という設定は不要。r の repeat は繰り返し終了判定を含まないので repeat ループの中で if を使う。

fibon1 = list()
fibon1[[1]] = fibon1[[2]] = as.bigz(1)
digits = 1000
i=2

repeat {
i=i+1
fibon1[[i]] = fibon1[[i-1]]+fibon1[[i-2]]
if (nchar(as.character(fibon1[[i]])) >= digits) break
}
i

プログラム書法の一つ:適切なところで条件判定をしよう

コメント

プログラムの使い方が間違っているのではないか?

2012年01月25日 | 裏 RjpWiki

npパッケージのnpcmstestでbwsで具体的にbandwidthを指定すると「仮引数 ”bws” が複数の実引数にマッチしました」というエラーがでてしまいます

ななみ? (2012-01-21 (土) 21:55:03)

はじめまして。npパッケージでモデルの検定をしています。

model <- lm(Y ~ X, y = TRUE, x = TRUE)
X <- data.frame(X)
npcmstest(model = model, xdat = X, ydat = Y, bws = 1)

のようにバンド幅を指定すると

「仮引数 ”bws” が複数の実引数にマッチしました」

というエラーがでてしまいます。
どなたか、これに対する解決方法をご存じでしたらお教え頂きますようお願い致します。

  • 再現できる最小のサンプルデータがないと答えるのが難しいと思う。bwsが一意でないエラーに見えるけど、サンプルがないと何とも言えない。 -- 2012-01-23 (月) 15:58:06

少なくとも「サンプルがないと何とも言えない」というのは,この場合たぶんまとはずれ。事前に,X <- rnorm(10); Y <- rnorm(10) とでも定義しておくだけで問題は再現する。

問題の核心は,確かに npcmstest では,... 引数で指定できるものとして,bws と bwscaling があり,このような状況では bws を指定すると件のエラーメッセージが出るのは必至。

にもかかわらずこれが問題にならない・なっていないのは,そもそもこのような使い方はないのか,単なる放置されたバグなのか,判断にこまるところ。

コメント

ダメ出し:自分で

2012年01月25日 | ブログラミング

この一つ前の記事で。

二重の 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[j] ^ 2 - 2 * x[i]) ^ 3 - 4 * x[i] ^ 2 * y[j] ^ 2
    }
}
contour(x, y, z, drawlabels = FALSE, levels = c(0, 1, 2, 5, 10, 20, 30, 50))
})

#   ユーザ   システム       経過 
#    11.962      0.119     12.111

10 倍速だ。

system.time({
x <- seq(-0.5, 3, length = 1000)
y <- seq(-2.5, 2.5, length = 1000)
z <- outer(x, y, function(xi, yj) (xi ^ 2 + yj ^ 2 - 2 * xi) ^ 3 - 4 * xi ^ 2 * yj ^ 2)
contour(x, y, z, drawlabels = FALSE, levels = c(0, 1, 2, 5, 10, 20, 30, 50))
})
#   ユーザ   システム       経過 
#     0.923      0.037      0.956


下記のような cxxfunction を使うのと遜色ない実行速度

library(inline)
src <- '
    Rcpp::NumericVector x(X);
    Rcpp::NumericVector y(Y);
    int nx = x.size();
    int ny = y.size();
    Rcpp::NumericMatrix z(nx, ny);
    for (int i = 0; i < nx; i++) {
        for (int j = 0; j < ny; j++) {
            z(i, j) = pow(x(i) * x(i) + y(j) * y(j) - 2 * x(i), 3) - 4 * x(i) * x(i) * y(j) * y(j);
        }
    }
    return wrap(z);
'
calcZ <- cxxfunction(signature(X="numeric", Y="numeric"), src, plugin="Rcpp")
x <- seq(-0.5, 3, length = 1000)
y <- seq(-2.5, 2.5, length = 1000)
z <- calcZ(x, y)
contour(x, y, z, drawlabels = FALSE, levels = c(0, 1, 2, 5, 10, 20, 30, 50))

もっとも,これを以下のようにするとちょっとだけ速くなる。ループ内で不変な計算はループの外でやることと,何回も同じ計算をするときには一回計算して変数に代入しておくこと。

src <- '
    Rcpp::NumericVector x(X);
    Rcpp::NumericVector y(Y);
    int nx = x.size();
    int ny = y.size();
    Rcpp::NumericMatrix z(nx, ny);
    for (int i = 0; i < nx; i++) {
        double xi2 = x(i) * x(i);
        for (int j = 0; j < ny; j++) {
            double yj2 = y(j) * y(j);
            z(i, j) = pow(xi2 + yj2 - 2 * x(i), 3) - 4 * xi2 * yj2;
        }
    }
    return wrap(z);
'

コメント

陰関数で表された曲線の描画

2012年01月24日 | 裏 RjpWiki

陰関数のプロットについて のコメントはちょっと外れているような気がする。

ちょっと,裏技的であるが,3次元関数として計算し,任意の等高線を描くというようにすることができることもある。
f(x, y) = (x^2+y^2-2x)^3- 4(xy)^2 = c のグラフの例。

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[j] ^ 2 - 2 * x[i]) ^ 3 - 4 * x[i] ^ 2 * y[j] ^ 2
    }
}
contour(x, y, z, drawlabels = FALSE, levels = c(0, 1, 2, 5, 10, 20, 30, 50))

f(x, y) = x^2-xy+y^2 = c のグラフの例。

x <- seq(-3, 3, length = 1000)
y <- seq(-3, 3, length = 1000)
z <- matrix(0, 1000, 1000)
for (i in 1:1000) {
    for (j in 1:1000) {
        z[i, j] <- x[i] ^ 2 - x[i]*y[j] + y[j] ^ 2
    }
}
contour(x, y, z, drawlabels = FALSE, levels = c(0, 0.2, 0.5, 1:25))

コメント

GPSデータの解析・処理

2012年01月20日 | ブログラミング

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)
    N <- ne$N
    E <- ne$E
    lat <- range(N)
    lon <- range(E)
    center <- c(mean(lat), mean(lon))
    if (is.null(size)) {
        size <- c(diff(lon), diff(lat))
        size <- round(640 * size / max(size))
    }
    if (is.null(zoom)) {
        zoom <- min(MaxZoom(lat, lon), size = size)
    }
    destfile <- paste(fn, "-orig.png", sep = "")
    map <- GetMap(center = center, zoom = zoom, destfile = destfile,

                  maptype = maptype, size = size)
    size <- dim(map$myTile)
    if (format == "png") {
        png(paste(fn, "-res.png", sep=""), width=size[2],
            height=size[1])
    }
    else {
        pdf(paste(fn, "-res.pdf", sep = ""), width = size[2] / 72 * 3,
                  height = size[1] / 72 *  3)
    }
    PlotOnStaticMap(map, ne$N, ne$E, destfile = destfile,
                    size = c(size[2], size[1]), FUN = lines,
                    col = col, lwd = lwd)
    dev.off()
}

decodeGPS <- function(fn) {
    base.dir <- "/foo"
    file.name <- paste(base.dir, "/", fn, ".log", sep = "")
    a <- readLines(file.name)
    a <- a[-1]
    n <- length(a)/3
    quality <- hh <- mm <- ss <- integer(n)
    km <- m <- N <- E <- numeric(n)
    cnt <- 0
    for (i in 1:n * 3 - 2) {
        d <- unlist(strsplit(a[i], ","))
        stopifnot(d[1] == "$GPGGA")
        cnt <- cnt + 1
        hh[cnt] <- as.integer(substr(d[2], 1, 2))
        mm[cnt] <- as.integer(substr(d[2], 3, 4))
        ss[cnt] <- as.integer(substr(d[2], 5, 10))
        N[cnt] <- as.double(substr(d[3], 1, 2)) +
                  as.double(substr(d[3],  3, 9)) / 60
        E[cnt] <- as.double(substr(d[5], 1, 3)) +
                  as.double(substr(d[5],  4, 10)) / 60
        quality[cnt] <- as.integer(d[7])
        m[cnt] <- as.double(d[10])
        d <- unlist(strsplit(a[i + 1], ","))
        stopifnot(d[1] == "$GPRMC")
        d <- unlist(strsplit(a[i + 2], ","))
        stopifnot(d[1] == "$GPVTG")
        km[cnt] <- as.double(d[8])
    }
    m[quality == 0] <- NA
    km[quality == 0] <- NA
    d <- data.frame(hh, mm, ss, E, N, m, km, quality)
    layout(matrix(c(1, 1, 2, 3), 2))
    par(mar = c(3, 3, 0.5, 0.5), mgp = c(1.8, 0.6, 0))
    plot(E, N, type = "l", asp = 1)
    points(E[km == 0], N[km == 0], col = 2, pch = 19, cex = 0.2)
    plot(m, type = "l", ylab = "標高(m)")
    plot(km, type = "l", ylab = "対地速度(km)")
    layout(1)
    return(list(N=N, E=E))
}

コメント

octave にある vander 関数

2012年01月18日 | ブログラミング

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]
[1,]    1
[2,]    1
[3,]    1
[4,]    1
[5,]    1

> vander(1:5, 2)
     [,1] [,2]
[1,]    1    1
[2,]    2    1
[3,]    3    1
[4,]    4    1
[5,]    5    1

> vander(1:5)
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    1    1    1    1
[2,]   16    8    4    2    1
[3,]   81   27    9    3    1
[4,]  256   64   16    4    1
[5,]  625  125   25    5    1

> vander(runif(5))
             [,1]         [,2]        [,3]      [,4] [,5]
[1,] 4.578474e-02 0.0989784306 0.213973690 0.4625729    1
[2,] 8.531321e-01 0.8876917796 0.923651483 0.9610679    1
[3,] 4.391669e-02 0.0959339217 0.209563093 0.4577806    1
[4,] 3.281008e-05 0.0004335163 0.005728008 0.0756836    1
[5,] 6.778974e-01 0.7470900815 0.823345243 0.9073837    1


コメント

ダメ出し:シミュレーションの条件設定

2012年01月18日 | ブログラミング

http://d.hatena.ne.jp/foo22222/20120115

Nt <- 1000
df1 <- 1
set.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, replace = FALSE)
p1 <- pchisq(t1, df1, lower.tail = FALSE)
p2 <- pchisq(t2, df1, lower.tail = FALSE)
plot(p1, p2)

としているが,
ts1 と ts2 は全く等しい(等しいようにしているから等しいのだけど)
ならば,ts1, ts2 の区別はする必要はないので余計なことは不要で

Nt <- 1000
df1 <- 1
set.seed(101)
ts1 <- rchisq(Nt, df1)
t1 <- sample(ts1, 500, replace = FALSE)
t2 <- sample(ts1, 500, replace = FALSE)

でよいし,そもそも rchisq( ) は無限母集団からの乱数抽出なのだから,そこから Nt 個取り出してそれを母集団として sample( ) により無作為抽出をするというような二段階の乱数生成をする必要はない。つまり,

Nt <- 1000
df1 <- 1
set.seed(101)
t1 <- rchisq(500, df1)
t2 <- rchisq(500, df1)
p1 <- pchisq(t1, df1, lower.tail = FALSE)
p2 <- pchisq(t2, df1, lower.tail = FALSE)
plot(p1, p2)

で十分。

コメント (1)

良いプログラムとはどういうものか

2012年01月11日 | ブログラミング

このブログで実行速度を比較していることがあるが,実行速度至上主義でないことはご理解頂きたい。

100 秒かかるものが1秒で終わっても,なんにもうれしくない。100 秒かかるものが 1.1 秒で終わろうが, 0.9 秒で終わろうが,問題ではない。

100 時間かかるものが 1 時間で終わればうれしいだろう。

そういうこと。

ただ,100 時間かかるものが 1 時間で終わるプログラムというのは,いろいろな意味で示唆に富むものであろうということ。

たとえば,もともと 100 時間も掛かるような問題ではなく,プログラムがまずいために 100 時間も掛かるような場合。

なぜそのようになっているかを解析しておくことは,将来に備えることになる。そもそも,いつもそのようにプログラムしていれば,そもそも「いざと言うことが起こらない」ということ。

コメント

ダメ出し:既存の関数(かなり有名)があれば,それを使おう

2012年01月10日 | ブログラミング

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

コメント (1)

cxxfunction array の実装 その2

2012年01月10日 | ブログラミング

return のときに次元を決めるのを前の記事の最後に書いたけど,やはり宣言する場所で決めるのが普通なので,そちらを。
しかし,いずれにしても,添え字は一個。

> library(inline)
> src <- '
+     int n1 = as<int>(N1);
+     int n2 = as<int>(N2);
+     int n3 = as<int>(N3);
+     Rcpp::NumericVector x(Dimension(n1, n2, n3));
+     x[1+n1*2+n1*n2*3] = 123;
+     x[0+n1*1+n1*n2*1] = 11;
+     x[0] = 0;
+     x[1] = 1;
+     x[5] = 5;
+     x[23] = 23;
+     return wrap(x);
+ '
> arry <- cxxfunction(signature(N1="integer", N2="integer", N3="integer"), src, plugin="Rcpp")
> arry(2, 3, 4)
, , 1

     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    1    0    5

, , 2

     [,1] [,2] [,3]
[1,]    0   11    0
[2,]    0    0    0

, , 3

     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0

, , 4

     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0   23

コメント

cxxfunction array の実装

2012年01月09日 | ブログラミング

array とはいっても,C 言語では単なる一次元ベクトルを三個以上の添え字でアクセスするだけ。

C 言語と R では,行と列の優先順序が正反対という点のみ。

「Rcpp::Arrayがないので」 では,本当の意味での array を操作していることにはなっていない。

つまり,添え字で array をアクセスするためには,以下のようにすべしということ。

# 必要なパッケージ
library(Rcpp)
library(inline)
# C++ ソース本体
src<-'
    int Nx=as<int>(N1);
    int Ny=as<int>(N2);
    int Nz=as<int>(N3);
    Rcpp::NumericVector a(Nx*Ny*Nz);
    for (int k = 0; k < Nz; k++) {
        for (int j = 0; j < Ny; j++) {
            for (int i = 0; i < Nx; i++) {
                a[i+Nx*j+Nx*Ny*k]= (i + 1)*(j + 1)*(k + 1)*1.0;
            }
        }
    }
    return a;
'
# R 内で、関数をコンパイル
fx.test1 <- cxxfunction(signature(N1 = "integer", N2 = "integer", N3 = "integer"),
                        include = c(""), src, plugin="Rcpp")
# 使う
library(rgl)
N1 <- 10
N2 <- 15
N3 <- 20
test.out <- fx.test1(N1, N2, N3)

test.out2 <- array(1, c(N1, N2, N3))
addresses <- which(test.out2 > 0, arr.ind = TRUE)

plot3d(addresses, col = gray((test.out - min(test.out)) /
                             (max(test.out) - min(test.out))))

A <- array(test.out, c(N1, N2, N3))

plot3d(slice.index(A, 1), slice.index(A, 2), slice.index(A, 3),
          col = rainbow(256)[A / (max(A)) * 256 + 1], alpha = ifelse(A != 0, 1, 0))

以下に,基本的例を

> library(inline)
> src <- '
+     int n1 = as<int>(N1);
+     int n2 = as<int>(N2);
+     int n3 = as<int>(N3);
+     Rcpp::NumericVector x(n1*n2*n3);
+     int i, j, k;
+     for (k = 0; k < n3; k++) {
+         for (j = 0; j < n2; j++) {
+             for (i = 0; i < n1; i++) {
+                 x[i+n1*j+n1*n2*k] = (i+1)*100+(j+1)*10+k+1;
+             }
+         }
+     }
+     return wrap(x);
+ '
>
> arry <- cxxfunction(signature(N1="integer", N2="integer", N3="integer"), src, plugin="Rcpp")
> a <- arry(2, 3, 4)
> array(a, dim=c(2, 3, 4))
, , 1

     [,1] [,2] [,3]
[1,]  111  121  131
[2,]  211  221  231

, , 2

     [,1] [,2] [,3]
[1,]  112  122  132
[2,]  212  222  232

, , 3

     [,1] [,2] [,3]
[1,]  113  123  133
[2,]  213  223  233

, , 4

     [,1] [,2] [,3]
[1,]  114  124  134
[2,]  214  224  234

配列の次元は,Rcpp の内部で .attr("dim") = Dimension(n1, n2, n3); によって定義することも出来る。

> library(inline)
> src <- '
+     int n1 = as<int>(N1);
+     int n2 = as<int>(N2);
+     int n3 = as<int>(N3);
+     Rcpp::NumericVector x(n1*n2*n3);
+     int i, j, k;
+     for (k = 0; k < n3; k++) {
+         for (j = 0; j < n2; j++) {
+             for (i = 0; i < n1; i++) {
+                 x[i+n1*j+n1*n2*k] = (i+1)*100+(j+1)*10+k+1;
+             }
+         }
+     }
+     x.attr("dim") = Dimension(n1, n2, n3);
+     return wrap(x);
+ '
>
> arry <- cxxfunction(signature(N1="integer", N2="integer", N3="integer"), src, plugin="Rcpp")
> ( a <- arry(2, 3, 4) )
, , 1

     [,1] [,2] [,3]
[1,]  111  121  131
[2,]  211  221  231

, , 2

     [,1] [,2] [,3]
[1,]  112  122  132
[2,]  212  222  232

, , 3

     [,1] [,2] [,3]
[1,]  113  123  133
[2,]  213  223  233

, , 4

     [,1] [,2] [,3]
[1,]  114  124  134
[2,]  214  224  234

コメント

R 用の indent その4

2012年01月07日 | ブログラミング

なぜか,これらの整形プログラムは,"/" などの前後には空白を「置けない」

> 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

コメント