lambda1 = exp(1)
lambda2 = exp(3)
(phi = logit(-1, inverse=TRUE))
y1 = rexp(n <- 1000, lambda1)
y2 = rexp(n, lambda2)
y = ifelse(runif(n) < phi, y1, y2)
fit = vglm(y ~ 1, mix2exp)
coef(fit, matrix=TRUE)
# Compare the results with the truth
round(rbind('Estimated'=Coef(fit),
'Truth'=c(phi, lambda1, lambda2)), dig=2)
# Plot the results
hist(y, prob=TRUE, main="red=estimate, blue=truth")
abline(v=1/Coef(fit)[c(2,3)], lty=2, col="red", lwd=2)
abline(v=1/c(lambda1, lambda2), lty=2, col="blue", lwd=2)
Run the code above in your browser using DataLab