x<-c(0.53495, 0.33676, 0.59671, 1.28840, 0.28603, 0.18338, 0.15082, 0.42654, 0.81693, 0.69833, 0.93841, 0.06836, 0.92658, 0.39314, 0.08227, 0.83003, 0.05288, 0.60613, 0.60410, 1.21632, 0.62962, 0.03493, 1.48802, 0.51257, 1.03245, 0.15208, 0.12643, 0.17467, 0.41677, 0.66373, 1.05877, 0.24880, 0.20516, 0.70610, 0.01899, 0.45714, 0.51985, 0.73178, 0.09021, 0.26240, 0.24517, 0.06909, 0.06497, 0.57614, 0.52521, 0.21695, 0.01568, 0.76230, 0.84435, 0.33877) loglike <-function(lam,x,c) { length(x)*log(lam)-lam*sum(x)-length(x)*lam*c } dlike<-function(lam,x,c) { (length(x)/lam) - sum(x)-length(x)*c } dlike2<-function(lam,x,c) { -length(x)/lam^2 } modified.newton<- function(fun, derf, derf2, x0,n,x,eps){ k<- 0 repeat { k<- k + 1 x1<- x0 - derf(x0,n,x)/derf2(x0,n,x) if(abs(x0 - x1) < eps) break x0<- x1 cat("Iteration No: ", k, "Current iterate = ",x1 ,"Current Likelihood = ",fun(x0,n,x),"dl/dx = ",derf(x0,n,x),"dl2/dx2=",derf2(x0,n,x),fill=T) } x1 } # produce a data vector x n<-length(data) c<-1.5 # plot the log-likelihood function ltemp<- seq(0.01,2,0.01) lhood<- loglike(ltemp,n,x) plot(ltemp,lhood,type="l",xlab="parameter values", ylab="log-likelihood", main="Log-likelihood function of the truncated exponential distribution") # initial value est<- mean(x) final.est<-modified.newton(loglike,dlike,dlike2,est,n,x,0.000001) final.est #truncexp.mle <-function(x,c) { #xbar <- mean(x) #lam0<-(1/length(x))*sum(x)+c*(1/(exp(c/xbar)-1)) #res<-nlm(nloglike,lam0,x=x) #list(lam = res$estimate[1],loglik = -res$minimum) #} #random number generation n.desired<-length(x) rtexp<-rexp(n.desired,1.0) %% 1.5