Tobit: Ошибка в lm.fit(X.vlm, y = z.vlm, ...): NA/NaN/Inf в "x" (смоделированные данные)
Я провожу исследование Монте-Карло. У меня есть линейная модель с гетероскедастичностью и левой цензурой зависимой переменной на 0. Среднее значение цензуры составляет 25,9.
Я получаю ошибку
Error in lm.fit(X.vlm, y = z.vlm, ...) : NA/NaN/Inf in 'x'
после попытки оценить модель тобита.
vglm(y[i,]~X[1,i,]+X[2,i,]+X[3,i,]+X[4,i,],family=tobit(Lower=0))
Мои данные смоделированы из стандартного распределения, поэтому проблема не должна быть связана с нечетными переменными.
Я нашел два других вопроса, которые имели ту же проблему с реальными данными: ошибка lm() NA/NaN/Inf, ошибка lm() NA/NaN/Inf Но, похоже, не было удовлетворительных ответов. Кроме того, мои данные легко воспроизводимы, поэтому это должно помочь выявить проблему
Вот коды:
library(VGAM)
set.seed(12345)
nobs=100
nsim=100
b=c(2,-2,-3,3)
g=c(1,0.2)
y=matrix(rep(0,nobs*nsim),ncol=nobs,nrow=nsim)
X=array(0,dim=c(4,nsim,nobs))
res=matrix(rep(0,nobs*nsim),ncol=nobs,nrow=nsim)
tobit=vector(mode="list",length=nsim)
for(i in 1:nsim){
# generate covariates :
X[1,i,]=rlnorm(n=nobs)
X[2,i,]=runif(n=nobs)<=.75
X[3,i,]=rnorm(mean = 3,n=nobs)
X[4,i,]=runif(n=nobs,min=0,max=10)
res[i,]=(g[1]+g[2]*X[4,i,])*rnorm(n=nobs)
# generate censored dependent variable
y[i,]=b[1]*X[1,i,]+b[2]*X[2,i,]+b[3]*X[3,i,]+b[4]*X[4,i,]+res[i,]
y[i,]=sapply(y[i,],FUN=function(x){max(0,x)}) #apply censoring
tobit[[i]]<-vglm(y[i,]~X[1,i,]+X[2,i,]+X[3,i,]+X[4,i,],
family = tobit(Lower=0))
}
Вот обратная связь
traceback()
5: lm.fit(X.vlm, y = z.vlm, ...)
4: vlm.wfit(xmat = X.vlm.save, z, Hlist = NULL, U = U, matrix.out =FALSE,
is.vlmX = TRUE, qr = qr.arg, xij = NULL)
3: vglm.fitter(x = x, y = y, w = w, offset = offset, Xm2 = Xm2,
Ym2 = Ym2, etastart = etastart, mustart = mustart, coefstart =coefstart,
family = family, control = control, constraints = constraints,
criterion = control$criterion, extra = extra, qr.arg = qr.arg,
Terms = mt, function.name = function.name, ...)
2: vglm(y[1, ] ~ X[1, 1, ] + X[2, i, ] + X[3, i, ] + X[4, i, ],
family = tobit(Lower = 0))
1: traceback(vglm(y[1, ] ~ X[1, 1, ] + X[2, i, ] + X[3, i, ] + X[4,
i, ], family = tobit(Lower = 0)))
*** Редактировать:
Удалив одну ковариату (я пробовал с X[3,i,] и X[4,i,]) и установив нижнюю цензуру на -0,001, как предполагает BondedDust, она работает нормально, и я даже увеличиваю количество повторений до 1000 без Основные проблемы.
Просто установив нижнюю цензуру на -0.001 и сохранив все ковариаты, я получаю две ошибки из 100 итераций. Стоит отметить, что ошибка сейчас
Error in lm.fit(X.vlm, y = z.vlm, ...) : NA/NaN/Inf in 'y'
Кроме того, я получаю эти предупреждения
In vglm.fitter(x = x, y = y, w = w, offset = offset, Xm2 = Xm2, ... :
iterations terminated because half-step sizes are very small
1 ответ
Я заметил, что это воспроизводимо не работает при i=1, поэтому подумал, что может быть проблема с vglm позвони сам. Глядя на примеры в ?tobit Я добавил некоторые параметры, связанные с цензурированными дистрибутивами, и начал получать несколько дополнительных итераций. Затем я попытался сузить диапазон цензуры и добился большего успеха с провалом только в 10% случаев. Поэтому я наконец добавил оболочку try(), чтобы цикл мог повторяться без остановки вычислений, и получил большинство успешных запусков с:
for(i in 1:nsim){
X[1,i,]=rlnorm(n=nobs)
X[2,i,]=runif(n=nobs)<=.75
X[3,i,]=rnorm(mean = 3,n=nobs)
X[4,i,]=runif(n=nobs,min=0,max=10)
res[i,]=(g[1]+g[2]*X[4,i,])*rnorm(n=nobs)
y[i,]=b[1]*X[1,i,]+b[2]*X[2,i,]+b[3]*X[3,i,]+b[4]*X[4,i,]+res[i,]
y[i,]=pmax(0,y[i,])
tobit[[i]]<-try( vglm(y[i,]~X[1,i,]+X[2,i,]+X[3,i,]+X[4,i,], crit = "coeff",
family = tobit(Lower=-.001, Upper=30, type.f = "cens")) )
}
Обратите внимание, что я заменяю вашу неуклюжую и, возможно, неэффективную sapply( ... max) с эквивалентным pmax,
> table( sapply(tobit, class))
try-error vglm
12 88
Вы можете перебрать успешные результаты с:
sapply( tobit[ sapply(tobit, class) == "vglm"], coefficients)
Топ результатов:
[,1] [,2] [,3] [,4] [,5] [,6]
(Intercept):1 2.8460081 1.910137 1.672237 1.2888827 2.4970536 1.0006290
(Intercept):2 0.9183935 1.042424 1.094658 0.9767228 0.9263946 0.9250609
X[1, i, ] 1.7777788 1.880506 1.662835 1.6204394 1.4412304 1.6275208
X[2, i, ] -3.0847792 -0.453110 -1.152709 -0.9900163 -2.4705355 -0.9651577
X[3, i, ] -2.4272169 -2.094114 -2.314748 -2.4628501 -1.9001385 -2.1076416
X[4, i, ] 2.6225234 2.245107 2.460182 2.7027493 2.3653673 2.3841989
[,7] [,8] [,9] [,10] [,11] [,12]
(Intercept):1 0.9520376 1.6319010 1.572563 1.4709517 1.616158 2.4992492
(Intercept):2 0.8698777 0.9005506 1.147485 0.9285724 1.012186 0.9229233
X[1, i, ] 1.6483879 1.6789573 1.718641 1.6544123 1.599116 1.7204001
X[2, i, ] -0.3718720 -1.8690782 -2.408657 -1.7278915 -1.208939 -2.0037999
X[3, i, ] -2.2601637 -1.9118288 -2.359274 -1.7828438 -2.257556 -2.3778443
X[4, i, ] 2.5381367 2.3091630 2.583869 2.3582418 2.333988 2.4389336
После достижения этой скромной степени успеха я попытался установить нижний уровень обратно на 0 и получил все ошибки. Увеличение значения Upper, похоже, не влияло на показатели успеха в ограниченном тестировании. Я не могу объяснить эти выводы, но, возможно, можно посоветоваться с автором пакета.