R Under development (unstable) (2024-09-18 r87177 ucrt) -- "Unsuffered Consequences" Copyright (C) 2024 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > > library("mlt") Loading required package: basefun Loading required package: variables > library("mlt.docreg") > library("truncreg") Loading required package: maxLik Loading required package: miscTools Please cite the 'maxLik' package as: Henningsen, Arne and Toomet, Ott (2011). maxLik: A package for maximum likelihood estimation in R. Computational Statistics 26(3), 443-458. DOI 10.1007/s00180-010-0217-1. If you have questions, suggestions, or comments regarding the 'maxLik' package, please use a forum or 'tracker' at maxLik's R-Forge site: https://r-forge.r-project.org/projects/maxlik/ > ### make M1mac happy > options(digits = 4) > > > ## Left-truncated > set.seed(29) > n <- 1000 > x <- runif(n, max = 2 * pi) > y <- rnorm(n, 2*x + 1, .25) > d <- data.frame(y = y, x = x) > ## truncated response > tr <- 2.5 > d$yt <- ifelse(d$y > tr, d$y, NA) > d$trunc_left <- tr > > tmod <- truncreg(yt ~ x, data = d, point = tr, direction = "left") > coef(tmod) (Intercept) x sigma 1.0193 1.9918 0.2575 > logLik(tmod) 'log Lik.' -33.99 (df=3) > vcov(tmod) (Intercept) x sigma (Intercept) 4.972e-04 -1.147e-04 -7.445e-06 x -1.147e-04 3.146e-05 1.637e-06 sigma -7.445e-06 1.637e-06 3.941e-05 > > ## MLT > fm <- as.formula("~ x") > yb <- polynomial_basis(numeric_var("yt", support = range(d$yt, na.rm = TRUE)), + coef = c(TRUE, TRUE), ci = c(-Inf, 0)) > m <- ctm(yb, shifting = fm, todistr = "Normal", data = d) > dfull <- d[!is.na(d$yt),,drop = FALSE] > dfull$yt <- R(dfull$yt, tleft = dfull$trunc_left) > mltmod <- mlt(m, data = dfull) > (cf <- coef(mltmod)) (Intercept) yt x -3.958 3.883 -7.734 > c(-cf[1] / cf[2], -cf[3] / cf[2], 1 / cf[2]) (Intercept) x yt 1.0193 1.9918 0.2575 > logLik(mltmod) 'log Lik.' -33.99 (df=3) > vcov(mltmod) (Intercept) yt x (Intercept) 0.017689 -0.009565 0.01723 yt -0.009565 0.008957 -0.01774 x 0.017226 -0.017744 0.03563 > > checkGH(mltmod) [1] TRUE > > ## right-truncated > > set.seed(29) > n <- 1000 > x <- runif(n, max = 2 * pi) > y <- rnorm(n, 2*x + 1, .25) > d <- data.frame(y = y, x = x) > ## truncated response > tr <- 11 > d$yt <- ifelse(d$y < tr, d$y, NA) > d$trunc_right <- tr > > tmod <- truncreg(yt ~ x, data = d, point = tr, direction = "right") > coef(tmod) (Intercept) x sigma 1.000 1.995 0.264 > logLik(tmod) 'log Lik.' -50.52 (df=3) > vcov(tmod) (Intercept) x sigma (Intercept) 3.312e-04 -1.022e-04 -2.841e-06 x -1.022e-04 4.353e-05 1.855e-06 sigma -2.841e-06 1.855e-06 4.526e-05 > > ## MLT > fm <- as.formula("~ x") > yb <- polynomial_basis(numeric_var("yt", support = range(d$yt, na.rm = TRUE)), + coef = c(TRUE, TRUE), ci = c(-Inf, 0)) > m <- ctm(yb, shifting = fm, todistr = "Normal", data = d) > dfull <- d[!is.na(d$yt),,drop = FALSE] > dfull$yt <- R(dfull$yt, tright = dfull$trunc_right) > mltmod <- mlt(m, data = dfull) > (cf <- coef(mltmod)) (Intercept) yt x -3.789 3.788 -7.557 > c(-cf[1] / cf[2], -cf[3] / cf[2], 1 / cf[2]) (Intercept) x yt 1.000 1.995 0.264 > logLik(mltmod) 'log Lik.' -50.52 (df=3) > vcov(mltmod) (Intercept) yt x (Intercept) 0.014383 -0.009473 0.01733 yt -0.009473 0.009316 -0.01849 x 0.017333 -0.018485 0.03730 > > checkGH(mltmod) [1] TRUE > > > proc.time() user system elapsed 2.17 0.37 2.53