0

I was trying to do simulation on mle for my model which is a parallel exponential model to estimate the covariate which is b0 and b1. t and x I got generate using random numbers.

My pdf is f(t) = 2 * (lambda * exp(-lambda * t)) * (1- exp(-lambda * t))^2 The lambda is the covariate. lambda = b0 + b1*x

 n=20
 u<- runif(n,0,1)
 x<- rnorm(n,0,1)
 
 b0=2 ; b1=4
   t1= -log(1-sqrt(u))/(b0+b1*x) #inverse method
 
 c1<- rexp(n,0.01)
 c<- 1*(t1<c1)
 t= pmin(t1,c1)
 
 data1<-data.frame(x,t,t1,c1,c)
 
 
 #mle
 LLF<- function(para){
   set.seed(1)
   
   b0=para[1]
   b1=para[2]
   
   z1= log(2*(b0+b1*x)) - (b0+b1*x)*t + log(1-exp(-(b0+b1*x)*t))
   
   j=sum(z1)
   return(j)
 }
 mle<-maxLik(LLF,start = c(2,4))
Warning message:
In log(2 * (b0 + b1 * x)) : NaNs produced

How can I estimate MLE when I am using random numbers that include negative and log cannot compute negative numbers.

zira
  • 15
  • 3
  • 1
    Thoughts: (1) don't name objects "c" in R because `c()` is a common function and this often leads to confusion, (2) you should set a seed before you randomly generate numbers (ie, before rnorm), you should not set a seed in your LLF function, (3) as currently specified, f(t,lambda) is not a density since lambda can be negative (a negative lambda yields a negative value for f(t) and thus f(t) is not a density). If you want to force lambda to be positive, you should parameterize lambda=exp(b0+b1x). (4) when posting on SO, you need to specify the library for non-base R functions like `maxLik()` – DanY Nov 25 '20 at 18:31
  • actually, this is only a segment of my coding. I already put maxLik() and set.seed before the coding. But the lambda=exp(b0+b1x) works! Thank you so much no more NaNs produced. – zira Nov 26 '20 at 06:12

0 Answers0