階層型多項ロジットに ChoiceModelR を使用しています。外部財 (正規分布に従う) の効用を推定したいと考えています。外側の商品には内側の商品のような共変量はありません。たとえば、価格やブランドのダミーを持つことはできません。 ChoiceModelR) にのみ適用されますが、y (選択) データのみに適用されます。
反復は正常に開始し、ある時点で停止して次のように言います
"Error in betadraw[good, ] = newbeta[good, ] : NAs are not allowed in subscripted assignments".
これは、関数「choicemodelr」の行 388 で、「良い」添え字が NA であるために発生する可能性があります。
私は Choicemodelr ( this、this、this ) と添え字の NA ( this、this ) に関するいくつかの質問を見ましたが、私の推測では、おそらく反復のいくつかの入力がちょうど「良い」がNAになるように大きく/小さくなります。
以下は非常に簡単な例です。さまざまな属性を持つ 3 つの製品でデータを生成します。期間の半分では、製品 3 は提供されません。2000 人の消費者には、通常分布する 3 つの属性 (および外部財に対する嗜好) の好みがあります。モデルと一致するようにロジット エラーが追加されました。社外品は製品 4 として索引付けされます (選択セットに 3 つの製品と 2 つの製品があった場合の両方)。
NA エラーを回避するにはどうすればよいですか? 私は何か間違ったことをしていますか、それとも関数の一般的なバグですか?
また、オプション none=TRUE を設定してオンラインで例を検索しましたが、再現可能なものは見つかりませんでした。none=FALSE に設定した場合に真のパラメーターを回復するのに問題はなく、顧客に外部オプションを選択させないため、おそらくこのオプションは問題のあるものにすぎません。
したがって、NA バグが発生するコードは次のとおりです。
library("ChoiceModelR")
library("MASS")
set.seed(36)
# Set demand pars
beta_mu = c(-3,4,1)
beta_sigma = diag(c(1,1,1))
alfa_mu = 5 #outside good mean utility
alfa_sigma = 2 #outside good sd
# Three/two products, 3 vars (2 continuous,1 dummy)
threeprod <- list()
twoprod <- list()
purchase <- list()
for (t in 1:1000){
threeprod[[t]] = cbind(rep(t,3),c(1,1,1),c(1,2,3),runif(3),runif(3),ceiling(runif(3,-0.5,0.5)))
purchase[[t]] = which.max(rbind(threeprod[[t]][,c(4,5,6)]%*%mvrnorm(1,beta_mu,beta_sigma) +
matrix( -log(-log(runif(3))), 3, 1),rnorm(1,alfa_mu,alfa_sigma)) )
threeprod[[t]] = cbind(threeprod[[t]],c(purchase[[t]],0,0))
}
for (t in 1001:2000){
twoprod[[t]] = cbind(rep(t,2),c(1,1),c(1,2),runif(2),runif(2),ceiling(runif(2,-0.5,0.5)))
purchase[[t]] = which.max(rbind(twoprod[[t]][,c(4,5,6)]%*%mvrnorm(1,beta_mu,beta_sigma) +
matrix( -log(-log(runif(2))), 2, 1),rnorm(1,alfa_mu,alfa_sigma)) )
if (purchase[[t]] == 3) {purchase[[t]] <- 4}
twoprod[[t]] = cbind(twoprod[[t]],c(purchase[[t]],0))
}
X <- rbind(do.call(rbind,threeprod),do.call(rbind,twoprod))
xcoding <- c(1,1,1)
mcmc = list(R = 5000, use = 2000)
options = list(none=TRUE, save=TRUE, keep=5)
out = choicemodelr(X, xcoding, mcmc = mcmc,options = options)