MLE ковариантного параметра

Привет, я сейчас занимаюсь кодированием для моделирования данных с использованием обратного метода. Я использую параллельную экспоненциальную модель, где я позволяю лямбда = b0 + b1x. Моя симуляция основана на анализе выживания.

#generate data
gen <- function(n,lambda,b0,b1){
   set.seed(1)
   
   u <- runif(n,0,1)
   c1 <- rexp(n,lambda)
   x <- rnorm(n,0,1)
   
   t1 = -log(1 - sqrt(u) ) / (b0 + b1*x) #inverse method
   
   c <- 1*(t1 < c1)
   t = pmin(t1, c1)
   
   data1 <- data.frame(x, t, t1, c1, c)
   return(data1)
 }

data2 <- gen(20,0.01,2,4)
data2
x = data2$x
t = data2$t
xsum = sum(x)
tsum = sum(t)

Проблема в том, что при запуске второго кода ниже он не покажет мои mle для b0 и b1

#Likelihood
library(maxLik)
LLF <- function(para){
  set.seed(1)
  
  b0 = para[1]
  b1 = para[2]
  
  n = 1
  
  z1 = (n*log(2)) + (n*log(b0+b1*xsum)) - ((b0+b1*xsum)*tsum) + (n*log(1-exp((-(b0 + b1*xsum)*tsum))))

  return(z1)
}

mle <- maxLik(LLF, start = c(2,4))

r estimation inverse mle
person zira    schedule 21.11.2020    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Проблема в том, что вы назначили n=1 в файле LLF. Поскольку мы обычно максимизируем параметры, учитывая все данные, n должно быть равно количеству наблюдений. Если вы обновите эту информацию, ваши mle сойдутся. Например,

n<-nrow(data2)

#Likelihood
library(maxLik)
LLF <- function(para){
set.seed(1)

b0 = para[1]
b1 = para[2]

#n = 1

z1 = (n*log(2)) + (n*log(b0+b1*xsum)) - ((b0+b1*xsum)*tsum) + (n*log(1-exp((-(b0 + b1*xsum)*tsum))))

return(z1)
}

mle <- maxLik(LLF, start = c(2,4))
summary(mle)

Maximum Likelihood estimation
Newton-Raphson maximisation, 3 iterations
Return code 1: gradient close to zero
Log-Likelihood: -22.7055 
2  free parameters
Estimates:
     Estimate Std. error t value Pr(> t)
[1,]    1.986         NA      NA      NA
[2,]    3.986         NA      NA      NA
person ssaha    schedule 21.11.2020
comment
Спасибо! Оно работает. Теперь единственная проблема — стандартное возвращаемое значение ошибки NA. z1 должно быть логарифмическим правдоподобием, верно? Я не уверен в этой части. - person zira; 22.11.2020
comment
@зира: Здорово! Пожалуйста, отметьте значок галочки рядом с ответом, чтобы пометить его как отвеченный, и проголосуйте за мой ответ. :) - person ssaha; 22.11.2020
comment
Мне сложно что-либо подсказать по поводу z1, так как информации о типе модели очень мало. - person ssaha; 22.11.2020