0

R の経験はあまりありません。次のような for ループがある Gibbs サンプラーを作成しようとしています。

for (iNum in 1:totNum) {
    rateNum <- Y3[iNum]
    if(Y3[iNum] > 0) {
        yStar3[iNum] <- rtnorm(1, mean = Mean3[iNum], sd = sqrt(Var3), 
                 lower = gz[rateNum], upper = gz[rateNum + 1])
    } else if(Y3[iNum] == 0) {
    yStar3[iNum] <- rtnorm(1, mean = Mean3[iNum], sd = sqrt(Var3), 
                  lower = -Inf, upper = Inf);
    }
}

これには時間がかかりすぎています。私は使用しようとしましlapplyたが、それも十分に高速ではありません。このループをベクトル化する方法はありますか?

ありがとうございます。

4

3 に答える 3

2

最も簡単な方法は、条件付きの2つの半分を生成し、どちらを使用するかを選択することです。meanパラメータはベクトル平均を取るので、次のようになります。

yStar3 <- ifelse(
  Y3 > 0,
    rtnorm(totNum, mean=Mean3, sd=sqrt(Var3), lower=gz[ratenum], upper=gz[rateNum+1]),
    rtnorm(totNum, mean=Mean3, sd=sqrt(Var3), lower=-Inf, upper=Inf))

ifelseY3がゼロ未満の場合は、追加の条件を使用して、を調整する必要がありますが、これが一般的な考え方です。

更新：@hadleyは、ifelseをrtnorm内に移動することを提案しています：

yStar3 <- rtnorm(totNum, mean=Mean3, sd=sqrt(Var3),
  lower=ifelse(Y3>0,gz[rateNum], -Inf),
  upper=ifelse(Y3>0,gz[rateNum+1], Inf))

現在、不要な計算は基本的にゼロです。

更新：もちろん、コメント投稿者が指摘したように、1は間違っています。代わりにtotNumにする必要があります。

于 2013-02-21T20:02:28.223 に答える
2

そのため、反復間に依存関係があるようには見えないため、ベクトル化が非常に簡単になります

  lhs = rtnorm(length(Y3), mean = Mean3, sd = sqrt(Var3), lower = gz[Y3],
              upper = gz[Y3 + 1])
  rhs = rtnorm(length(Y3), mean=Mean3, sd = sqrt(Var3), lower=-Inf, upper=Inf)

  ifelse(Y3 > 0, lhs, rhs)

ここでの問題は、rtnorm をその入力パラメーター、平均、下限、および上限でベクトル化する必要があることです。そうでない場合もあります。その場合は、さらに作業を行う必要があります。

于 2013-02-21T20:00:24.827 に答える
1

これは、変数にいくつかの値がないと少し問題になりますが、やりたいことはかなり簡単です。このためにベクトル化されたすべてのステートメントに固執し、メモリを大量に消費しないようにします。これが基本的な戦略です。

ステップ 1 : すべての数値を計算する方法を理解します。

# The number of values you need from 'rtnorm'
sum(Y3 > 0)
sum(Y3 == 0)

# The means you need from the 'Mean3' array
Mean3[Y3 > 0]
Mean3[Y3 == 0]

# Lower and upper limits for Y3 > 0
gz[Y3[Y3 > 0]]
gz[Y3[Y3 > 0] + 1]

ステップ 2 : これらの値を yStar3 のベクトル フィルターで使用します。サンプル データと変数値がなければ、すべての構文が完全であるとは断言できませんが、次のようになります。

yStar3[Y3 > 0] <- rtnorm(
  sum(Y3 > 0), 
  mean = Mean3[Y3 > 0], 
  sd = sqrt(Var3), 
  lower = gz[Y3[Y3 > 0]], 
  upper = gz[Y3[Y3 > 0] + 1])

yStar3[Y3 == 0] <- rtnorm(
  sum(Y3 == 0), 
  mean = Mean3[Y3 == 0], 
  sd = sqrt(Var3), 
  lower = -Inf, 
  upper = Inf)
于 2013-02-21T20:25:46.543 に答える