Неограниченная ошибка конечных разностей все еще существует даже при использовании mle2
Для предъявления по месту требования,
Я выполнил анализ MLE для отдельного набора данных задачи двоичного выбора (каждый набор данных включает в себя несколько испытаний, в которых конкретный субъект должен выбрать вариант A или B) с помощью R. Сначала я попробовал optim (method=BFGS). Для некоторых людей я не могу получить оценки параметров из-за "конечной разностной ошибки". Я искал в Интернете и обнаружил, что набор инструментов bbmle может исправить эту ошибку, предоставив "опцию инвертирования матрицы Гессе" во время вычислений. Однако ошибка все еще появлялась, даже когда я принял функцию mle2 и установил "skip.hessian=T". Я знаю, что могут быть некоторые другие решения, такие как изменение первоначального предположения; однако, действительно раздражает, что иногда обновленное начальное предположение работает для этого предмета, но не работает для другого предмета. Может кто-нибудь сказать мне, как я могу решить эту проблему в mle2? Я с нетерпением жду вашего ответа. Спасибо и хорошего дня:-)
Бест, Ян
Извините, что не предоставил тестовые данные, так как у меня нет онлайн-диска. Вы можете получить тестовые данные по электронной почте (huyang@uni-bonn.de)
Мой код следующий:
##############################################################
### Step 1: Define the Log-likelihood Function of M1 ###
##############################################################
logl <- function(alpha,beta,delta1,delta2,delta3,epsilon){
data<-df_sbj
y <- data$Y
s1 <- data$s1
o1 <- data$o1
s2 <- data$s2
o2 <- data$o2
adv1<- ifelse(data$s1>data$o1,c(data$s1-data$o1),c(0))
adv2<- ifelse(data$s2>data$o2,c(data$s2-data$o2),c(0))
disadv1<- ifelse(data$o1>data$s1,c(data$o1-data$s1),c(0))
disadv2<- ifelse(data$o2>data$s2,c(data$o2-data$s2),c(0))
u1<-vector()
u2<-vector()
for (j in 1:nrow(data)) {
if (data$f1[j]=="human" & data$f2[j]=="greedy"){
u1[j]<-s1[j]-alpha*disadv1[j]-beta*adv1[j]+delta1*adv1[j]-delta1*disadv1[j]
u2[j]<-s2[j]-alpha*disadv2[j]-beta*adv2[j]+delta1*adv2[j]-delta1*disadv2[j]
} else if (data$f1[j]=="human" & data$f2[j]=="equal") {
u1[j]<-s1[j]-alpha*disadv1[j]-beta*adv1[j]-delta2*adv1[j]+delta2*disadv1[j]
u2[j]<-s2[j]-alpha*disadv2[j]-beta*adv2[j]-delta2*adv2[j]+delta2*disadv2[j]
} else if (data$f1[j]=="human" & data$f2[j]=="generous") {
u1[j]<-s1[j]-alpha*disadv1[j]-beta*adv1[j]-delta3*adv1[j]+delta3*disadv1[j]
u2[j]<-s2[j]-alpha*disadv2[j]-beta*adv2[j]-delta3*adv2[j]+delta3*disadv2[j]
} else {
u1[j]<-s1[j]-alpha*disadv1[j]-beta*adv1[j]
u2[j]<-s2[j]-alpha*disadv2[j]-beta*adv2[j]
}
}
#u1<-s1-alpha%*%disadv1-beta%*%adv1
#u2<-s2-alpha%*%disadv2-beta%*%adv2
udiff<-epsilon%*%(u2-u1)
# Use the log-likelihood of the Bernouilli distribution, where p is
# defined as the logistic transformation of a linear combination
# of predictors, according to logit(p)=(x%*%beta)
loglik <- sum(-y*log(1 + exp(-udiff)) - (1-y)*log(1 + exp(udiff)))
return(-loglik)
}
##############################################################
### Step 2: Prepare the data ###
##############################################################
df_sbj<-read.csv('./model/testdata_bbmle.csv',header=T)
fit <- mle2(logl, start = list(alpha = 0.2, beta=0.2, delta1=0.02, delta2=0.02, delta3=0.02, epsilon=0),method="BFGS",data=df_sbj,skip.hessian=T)