2016-08-23 9 views
3

data.frameのデータ操作に関する質問があります。各行の要素に関数を適用して要約する

本質的に、私は大規模なデータセットを持っている - 以下の略記バージョン:私は

は同様に、私がしたいパラメータとしてベータ分布とαおよびβを用いて、確率を生成するrbetaを使用する

structure(list(nm_mean = c(194213914.326, 194213914.326, 194213914.326, 
194213914.326, 194213914.326, 217947112.739), nm_se = c(9984735.05918367, 
9984735.05918367, 9984735.05918367, 9984735.05918367, 9984735.05918367, 
11010386.0760204), alpha = c(193.197697846336, 214.592588477741, 
240.246557258741, 258.116959355425, 282.560024775668, 306.610038660465 
), beta = c(61526.2664158025, 57950.9563448233, 56085.1512614369, 
52919.4794239927, 51483.4591654126, 50405.8186695088)), .Names = c("nm_mean", 
"nm_se", "alpha", "beta"), row.names = c(NA, 6L), class = "data.frame") 

nmnとnm_seを平均とsdとする正規分布を使用して乱数を生成するには、rnormを使用します。

Iは、行本質的に1

x <- rbeta(1000,193.1977,61526.27) 
y <- rnorm(1000,194213914,9984735) 
z <- x*y 

dat$ce <- quantile(z,0.5) 
dat$ll <- quantile(z,0.25) 
dat$ul <- quantile(z,0.975) 

ための例として

だからバックデータフレームにrnorm値によって生成rbeta値を乗算し、第50、第25及び第75分位数を抽出します私は、lbetaとrnormの積のためにce、llとulをデータベースに追加します。

+1

SOで '[r] data frame apply rows'を検索すると、多くの良い答えが得られます。 TL DR: 'apply'の使用は、data.frameを数値を文字に変換する'配列 '(' matrix')に変換するので問題になります。ベクトル化された計算、ループ、 'plyr'、' dplyr'はあなたの友達です。 – r2evans

+1

'quantile'は確率のベクトルをとることができます。 – shayaa

答えて

1

これは@thelatemailと私の会話に基づいてベクトル化ソリューションです:

n <- 1000 
grp <- nrow(dat) 
z <- with(dat, rnorm(grp*n, nm_mean, nm_se) * rbeta(grp*n, alpha, beta)) 
m <- 1 

for(i in 1:nrow(dat)){ 
    dat$ce[i] <- quantile(z[m:(i*1000)],0.5) 
    dat$ll[i] <- quantile(z[m:(i*1000)],0.25) 
    dat$ul[i] <- quantile(z[m:(i*1000)],0.975) 
    m <- m + 1000 
} 

少ないベクトル化ソリューションは、次のとおりです。

HackRのコード@が動機
for(i in 1:nrow(dat)){ 
    x <- rbeta(1000, shape1 = dat$alpha[i], shape2 = dat$beta[i]) 
    y <- rnorm(n=1000,dat$nm_mean[i],dat$nm_se[i]) 
    z <- x*y 

    dat$ce[i] <- quantile(z,0.5) 
    dat$ll[i] <- quantile(z,0.25) 
    dat$ul[i] <- quantile(z,0.975) 
} 

dat 
nm_mean nm_se alpha  beta  ce  ll  ul 
1 194213914 9984735 193.1977 61526.27 607563.9 573229.9 713057.2 
2 194213914 9984735 214.5926 57950.96 712268.5 674826.3 836950.8 
3 194213914 9984735 240.2466 56085.15 823322.9 777482.8 981156.7 
4 194213914 9984735 258.1170 52919.48 937331.2 884945.0 1095876.3 
5 194213914 9984735 282.5600 51483.46 1059980.4 1003596.4 1225615.6 
6 217947113 11010386 306.6100 50405.82 1316733.1 1250190.1 1515185.0 
+0

これは確かにうまくいくでしょうが、おそらく、 'rbeta'と' rnorm'をベクトル化し、アルファ/ベータ/平均/ sdで繰り返し実行するだけでよいのです。 – thelatemail

+0

@thelatemail良いアイデア。私はそれについて考えます。しかし、もし私がインデックスを取り除き、x、y、xの部分をループの外側に動かすと、1000行の結果のどれが 'dat 'のどの行に対応するのかをどうやって知ることができますか? –

+1

3列... 1,2,3,1,2,3と同じように繰り返されますので、1,2,3で分けたりタップしたりすることができます。 'n < - 10;と同じです。 grp < - nrow(dat); tmp < - with(dat、 rnorm(grp * n、nm_mean、nm_se)* rbeta(grp * n、alpha、beta) )「私はそうだと思います。 – thelatemail

6

、何私は機能的なベクトル化バージョンだと思う:

set.seed(42) 
n <- 1000 
nrows <- nrow(dat) 
rn <- matrix(rnorm(nrows * n, dat$nm_mean, dat$nm_se), ncol = nrows, byrow = TRUE) 
rb <- matrix(rbeta(nrows * n, shape1 = dat$alpha, shape2 = dat$beta), 
      ncol = nrows, byrow = TRUE) 
cbind(dat, 
     structure(t(apply(rn * rb, 2, function(z) quantile(z, c(0.5, 0.25, 0.975)))), 
       .Dimnames = list(NULL, c("ce", "ll", "ul")))) 
#  nm_mean nm_se alpha  beta  ce  ll  ul 
# 1 194213914 9984735 193.1977 61526.27 608455.3 570100.5 710373.6 
# 2 194213914 9984735 214.5926 57950.96 715305.0 677754.3 856570.7 
# 3 194213914 9984735 240.2466 56085.15 825143.7 778351.2 979361.1 
# 4 194213914 9984735 258.1170 52919.48 943261.4 895832.6 1091899.3 
# 5 194213914 9984735 282.5600 51483.46 1054514.3 995640.8 1226176.4 
# 6 217947113 11010386 306.6100 50405.82 1312325.0 1247030.8 1515630.5 
+0

ああ、コメント形式よりはるかに明確になっています。私はそれが好きです。 +1 –

+1

もちろん、私の中の[CDO](http://www.urbandictionary.com/define.php?term=CDO)は、それらを '' ll ''、 '' ce "'、 '' ' ul "' :-) – r2evans

関連する問題