CRAN Package Check Results for Package SemiParSampleSel

Last updated on 2019-12-05 11:50:21 CET.

Check Details

Version: 1.5
Check: examples
Result: ERROR
    Running examples in 'SemiParSampleSel-Ex.R' failed
    The error most likely occurred in:
    
    > base::assign(".ptime", proc.time(), pos = "CheckExEnv")
    > ### Name: SemiParSampleSel
    > ### Title: Semiparametric Sample Selection Modelling with Continuous or
    > ### Discrete Response
    > ### Aliases: SemiParSampleSel
    > ### Keywords: copula sample selection semiparametric sample selection
    > ### modelling sample selection model smooth regression spline shrinkage
    > ### smoother variable selection
    >
    > ### ** Examples
    >
    >
    > library(SemiParSampleSel)
    >
    > ######################################################################
    > ## Generate data
    > ## Correlation between the two equations and covariate correlation 0.5
    > ## Sample size 2000
    > ######################################################################
    >
    > set.seed(0)
    >
    > n <- 2000
    >
    > rhC <- rhU <- 0.5
    >
    > SigmaU <- matrix(c(1, rhU, rhU, 1), 2, 2)
    > U <- rmvnorm(n, rep(0,2), SigmaU)
    >
    > SigmaC <- matrix(rhC, 3, 3); diag(SigmaC) <- 1
    >
    > cov <- rmvnorm(n, rep(0,3), SigmaC, method = "svd")
    > cov <- pnorm(cov)
    >
    > bi <- round(cov[,1]); x1 <- cov[,2]; x2 <- cov[,3]
    >
    > f11 <- function(x) -0.7*(4*x + 2.5*x^2 + 0.7*sin(5*x) + cos(7.5*x))
    > f12 <- function(x) -0.4*( -0.3 - 1.6*x + sin(5*x))
    > f21 <- function(x) 0.6*(exp(x) + sin(2.9*x))
    >
    > ys <- 0.58 + 2.5*bi + f11(x1) + f12(x2) + U[, 1] > 0
    > y <- -0.68 - 1.5*bi + f21(x1) + + U[, 2]
    > yo <- y*(ys > 0)
    >
    >
    > dataSim <- data.frame(ys, yo, bi, x1, x2)
    >
    > ## CLASSIC SAMPLE SELECTION MODEL
    > ## the first equation MUST be the selection equation
    >
    > out <- SemiParSampleSel(list(ys ~ bi + x1 + x2,
    + yo ~ bi + x1),
    + data = dataSim)
     ----------- FAILURE REPORT --------------
     --- failure: the condition has length > 1 ---
     --- srcref ---
    :
     --- package (from environment) ---
    SemiParSampleSel
     --- call from context ---
    SemiParSampleSel(list(ys ~ bi + x1 + x2, yo ~ bi + x1), data = dataSim)
     --- call from argument ---
    if (class(X2s) == "try-error") stop("Check that the range of the covariate values of the selected sample is the same as that of the complete dataset.")
     --- R stacktrace ---
    where 1: SemiParSampleSel(list(ys ~ bi + x1 + x2, yo ~ bi + x1), data = dataSim)
    
     --- value of length: 2 type: logical ---
    [1] FALSE FALSE
     --- function from context ---
    function (formula, data = list(), weights = NULL, subset = NULL,
     start.v = NULL, start.theta = NULL, BivD = "N", margins = c("probit",
     "N"), fp = FALSE, infl.fac = 1, rinit = 1, rmax = 100,
     iterlimsp = 50, pr.tolsp = 1e-06, bd = NULL, parscale)
    {
     if (margins[2] == "GA")
     margins[2] <- "G"
     if (margins[1] == "probit")
     margins[1] <- "N"
     marg2 <- c("N", "G", "P", "NB", "D", "PIG", "S", "BB", "BI",
     "GEOM", "LG", "NBII", "WARING", "YULE", "ZIBB", "ZABB",
     "ZABI", "ZIBI", "ZALG", "ZANBI", "ZINBI", "ZAP", "ZIP",
     "ZIP2", "ZIPIG")
     marg2D <- c("P", "NB", "D", "PIG", "S", "BB", "BI", "GEOM",
     "LG", "NBII", "WARING", "YULE", "ZIBB", "ZABB", "ZABI",
     "ZIBI", "ZALG", "ZANBI", "ZINBI", "ZAP", "ZIP", "ZIP2",
     "ZIPIG")
     copul <- c("N", "FGM", "F", "AMH", "C0", "C90", "C180", "C270",
     "J0", "J90", "J180", "J270", "G0", "G90", "G180", "G270")
     if (!(margins[1] %in% c("N", "probit")))
     stop("Error in selection equation's margin name. It can currently only be probit or equivalently N.")
     if (!(margins[2] %in% marg2))
     stop("Error in outcome equation's margin name. It should be one of: N, GA, P, NB, D, PIG, S, BB, BI, GEOM, LG, NBII, WARING, YULE, ZIBB, ZABB, ZABI, ZIBI, ZALG, ZANBI, ZINBI, ZAP, ZIP, ZIP2, ZIPIG.")
     if (!(BivD %in% copul))
     stop("Error in parameter BivD value. It should be one of: N, FGM, F, AMH, C0, C90, C180, C270, J0, J90, J180, J270, G0, G90, G180, G270.")
     if (length(formula) > 2 && margins[2] %in% c("N", "G", "NB",
     "PIG", "BB", "NBII", "WARING", "ZABI", "ZIBI", "ZALG",
     "ZAP", "ZIP", "ZIP2")) {
     if (length(formula) != 4)
     stop("You need to specify four equations.")
     }
     if (length(formula) > 2 && margins[2] %in% c("P", "BI", "GEOM",
     "LG", "YULE")) {
     if (length(formula) != 3)
     stop("You need to specify three equations.")
     }
     if (length(formula) > 2 && margins[2] %in% c("D", "S", "ZIBB",
     "ZABB", "ZANBI", "ZINBI", "ZIPIG")) {
     if (length(formula) != 5)
     stop("You need to specify five equations.")
     }
     qu.mag <- sp <- nu <- gam2.1 <- KendTau <- gamlss.fit <- edf <- edf1 <- NULL
     sp3 <- gp3 <- gam3 <- X3 <- sp.theta <- NULL
     sp4 <- gp4 <- gam4 <- X4 <- NULL
     sp5 <- gp5 <- gam5 <- X5 <- NULL
     X3.d2 <- X4.d2 <- X5.d2 <- l.sp3 <- l.sp4 <- l.sp5 <- 0
     if (length(formula) > 2) {
     f3t <- try(formula[[3]][[3]], silent = TRUE)
     if (class(f3t) != "try-error")
     stop("The third equation does not require a response.")
     if (length(formula) > 3) {
     f4t <- try(formula[[4]][[3]], silent = TRUE)
     if (class(f4t) != "try-error")
     stop("The fourth equation does not require a response.")
     if (length(formula) > 4) {
     f5t <- try(formula[[5]][[3]], silent = TRUE)
     if (class(f5t) != "try-error")
     stop("The fifth equation does not require a response.")
     }
     }
     }
     ig <- interpret.gam(formula)
     mf <- match.call(expand.dots = FALSE)
     if (length(formula) == 2)
     pred.n <- union(ig[[1]]$pred.names, c(ig[[2]]$pred.names,
     ig[[2]]$response))
     if (length(formula) == 3)
     pred.n <- union(ig[[1]]$pred.names, c(ig[[2]]$pred.names,
     ig[[3]]$pred.names, ig[[2]]$response))
     if (length(formula) == 4)
     pred.n <- union(ig[[1]]$pred.names, c(ig[[2]]$pred.names,
     ig[[3]]$pred.names, ig[[4]]$pred.names, ig[[2]]$response))
     if (length(formula) == 5)
     pred.n <- union(ig[[1]]$pred.names, c(ig[[2]]$pred.names,
     ig[[3]]$pred.names, ig[[4]]$pred.names, ig[[5]]$pred.names,
     ig[[2]]$response))
     if (length(formula) == 6)
     pred.n <- union(ig[[1]]$pred.names, c(ig[[2]]$pred.names,
     ig[[3]]$pred.names, ig[[4]]$pred.names, ig[[5]]$pred.names,
     ig[[6]]$pred.names, ig[[2]]$response))
     fake.formula <- paste(ig[[1]]$response, "~", paste(pred.n,
     collapse = " + "))
     environment(fake.formula) <- environment(ig$fake.formula)
     mf$formula <- fake.formula
     mf$start.v <- mf$method <- mf$start.theta <- mf$BivD <- mf$margins <- mf$infl.fac <- mf$fp <- mf$rinit <- mf$rmax <- mf$iterlimsp <- mf$pr.tolsp <- mf$parscale <- mf$bd <- NULL
     mf$drop.unused.levels <- TRUE
     mf$na.action <- na.pass
     mf[[1]] <- as.name("model.frame")
     data <- eval(mf, parent.frame())
     indS <- as.logical(data[, ig[[1]]$response]) == FALSE
     indS <- ifelse(is.na(indS), FALSE, indS)
     data[indS, ig[[2]]$response] <- ifelse(is.na(data[indS, ig[[2]]$response]),
     0, data[indS, ig[[2]]$response])
     data <- na.omit(data)
     if (is.null(weights)) {
     weights <- rep(1, dim(data)[1])
     data$weights <- weights
     names(data)[length(names(data))] <- "(weights)"
     }
     else weights <- data[, "(weights)"]
     formula.eq1 <- formula[[1]]
     formula.eq2 <- formula[[2]]
     gam1 <- eval(substitute(gam(formula.eq1, binomial(link = "probit"),
     gamma = infl.fac, weights = weights, data = data), list(weights = weights)))
     y1 <- gam1$y
     inde <- y1 > 0
     if (margins[2] == "N")
     gam2 <- eval(substitute(gam(formula.eq2, gamma = infl.fac,
     data = data, weights = weights, subset = inde), list(weights = weights,
     inde = inde)))
     if (margins[2] == "G")
     gam2 <- eval(substitute(gam(formula.eq2, gamma = infl.fac,
     data = data, weights = weights, subset = inde, family = Gamma(link = "log")),
     list(weights = weights, inde = inde)))
     if (margins[2] %in% c("P", "NB", "D", "PIG", "S", "GEOM",
     "NBII", "WARING", "YULE", "ZANBI", "ZINBI", "ZAP", "ZIP",
     "ZIP2", "ZIPIG", "LG", "ZALG", "BB", "BI", "ZIBB", "ZABB",
     "ZABI", "ZIBI")) {
     gam2 <- eval(substitute(gam(formula.eq2, gamma = infl.fac,
     data = data, weights = weights, subset = inde, family = poisson(link = "log")),
     list(weights = weights, inde = inde)))
     }
     n <- length(y1)
     X1 <- model.matrix(gam1)
     X1.d2 <- dim(X1)[2]
     X2.d2 <- length(coef(gam2))
     X2 <- matrix(0, n, X2.d2, dimnames = list(c(1:n), c(names(coef(gam2)))))
     X2[inde, ] <- model.matrix(gam2)
     if (margins[2] == "LG") {
     y2 <- rep(1, n)
     }
     else {
     y2 <- rep(0, n)
     }
     y2[inde] <- gam2$y
     l.sp1 <- length(gam1$smooth)
     if (l.sp1 != 0)
     sp1 <- gam1$sp
     else sp1 <- NULL
     l.sp2 <- length(gam2$smooth)
     if (l.sp2 != 0)
     sp2 <- gam2$sp
     else sp2 <- NULL
     gp1 <- gam1$nsdf
     gp2 <- gam2$nsdf
     if (length(formula) == 3)
     theta <- rnorm(n)
     nad1 <- "theta"
     if (length(formula) == 4) {
     sigma <- rnorm(n)
     nad1 <- "sigma"
     theta <- rnorm(n)
     nad2 <- "theta"
     }
     if (length(formula) == 5) {
     sigma <- rnorm(n)
     nad1 <- "sigma"
     nu <- rnorm(n)
     nad2 <- "nu"
     theta <- rnorm(n)
     nad3 <- "theta"
     }
     if (length(formula) > 2) {
     formula.eq3 <- formula[[3]]
     formula.eq3 <- as.formula(paste(nad1, "~", formula.eq3[2],
     sep = ""))
     gam3 <- eval(substitute(gam(formula.eq3, weights = weights,
     data = data)))
     X3 <- model.matrix(gam3)
     X3.d2 <- llM <- dim(X3)[2]
     l.sp3 <- length(gam3$sp)
     if (l.sp3 != 0)
     sp3 <- sp.theta <- gam3$sp
     else sp3 <- NULL
     environment(gam3$formula) <- environment(gam2$formula)
     gp3 <- gam3$nsdf
     a.theta <- coef(gam3)
     }
     if (length(formula) > 3) {
     formula.eq4 <- formula[[4]]
     formula.eq4 <- as.formula(paste(nad2, "~", formula.eq4[2],
     sep = ""))
     gam4 <- eval(substitute(gam(formula.eq4, weights = weights,
     data = data)))
     X4 <- model.matrix(gam4)
     X4.d2 <- llM <- dim(X4)[2]
     l.sp4 <- length(gam4$sp)
     if (l.sp4 != 0)
     sp4 <- sp.theta <- gam4$sp
     else sp4 <- NULL
     environment(gam4$formula) <- environment(gam2$formula)
     gp4 <- gam4$nsdf
     a.theta <- coef(gam4)
     }
     if (length(formula) > 4) {
     formula.eq5 <- formula[[5]]
     formula.eq5 <- as.formula(paste(nad3, "~", formula.eq5[2],
     sep = ""))
     gam5 <- eval(substitute(gam(formula.eq5, weights = weights,
     data = data)))
     X5 <- model.matrix(gam5)
     X5.d2 <- llM <- dim(X5)[2]
     l.sp5 <- length(gam5$sp)
     if (l.sp5 != 0)
     sp5 <- sp.theta <- gam5$sp
     else sp5 <- NULL
     environment(gam5$formula) <- environment(gam2$formula)
     gp5 <- gam5$nsdf
     a.theta <- coef(gam5)
     }
     if ((l.sp1 != 0 || l.sp2 != 0 || l.sp3 != 0 || l.sp4 != 0 ||
     l.sp5 != 0) && fp == FALSE) {
     qu.mag <- S.m(gam1, gam2, gam3, gam4, gam5, l.sp1, l.sp2,
     l.sp3, l.sp4, l.sp5, margins, eq1 = "yes")
     sp <- c(sp1, sp2, sp3, sp4, sp5)
     qu.mag1 <- S.m(gam1, gam2, gam3, gam4, gam5, l.sp1, l.sp2,
     l.sp3, l.sp4, l.sp5, margins, eq1 = "no")
     sp1 <- c(sp2, sp3, sp4, sp5)
     }
     if (length(formula) < 3) {
     p.g1 <- predict(gam1)
     imr <- data$imr <- dnorm(p.g1)/pnorm(p.g1)
     data$bd <- bd
     formula.eq2.1 <- update.formula(formula.eq2, ~. + imr)
     if (margins[2] %in% c("P", "NB", "D", "PIG", "S", "NBII",
     "ZANBI", "ZINBI", "ZAP", "ZIP", "ZIP2", "ZIPIG")) {
     formula.eq2.1 <- update.formula(formula.eq2, round((. +
     mean(.[. >= 0]))/2) ~ . + imr)
     }
     else if (margins[2] %in% c("BB", "BI")) {
     formula.eq2.1 <- update.formula(formula.eq2, cbind(round(bd *
     (. + 0.5)/(bd + 1)), bd - round(bd * (. + 0.5)/(bd +
     1))) ~ . + imr)
     }
     else if (margins[2] %in% c("GEOM", "YULE", "WARING",
     "LG", "ZALG")) {
     formula.eq2.1 <- update.formula(formula.eq2, ~. +
     imr)
     }
     else if (margins[2] %in% c("ZABI", "ZIBI", "ZIBB", "ZABB")) {
     formula.eq2.1 <- update.formula(formula.eq2, cbind(.,
     bd - .) ~ . + imr)
     }
     if (margins[2] == "N")
     gam2.1 <- eval(substitute(gam(formula.eq2.1, gamma = infl.fac,
     data = data, weights = weights, subset = inde),
     list(weights = weights, inde = inde)))
     if (margins[2] == "G")
     gam2.1 <- eval(substitute(gam(formula.eq2.1, gamma = infl.fac,
     data = data, weights = weights, subset = inde,
     family = Gamma(link = "log")), list(weights = weights,
     inde = inde)))
     if (margins[2] %in% c("P", "NB", "D", "PIG", "S", "GEOM",
     "NBII", "WARING", "YULE", "ZANBI", "ZINBI", "ZAP",
     "ZIP", "ZIP2", "ZIPIG", "LG", "ZALG")) {
     gam2.1 <- eval(substitute(gam(formula.eq2.1, gamma = infl.fac,
     data = data, weights = weights, subset = inde,
     family = poisson(link = "log")), list(weights = weights,
     inde = inde)))
     }
     if (margins[2] %in% c("BB", "BI", "ZIBB", "ZABB", "ZABI",
     "ZIBI")) {
     gam2.1 <- eval(substitute(gam(formula.eq2.1, gamma = infl.fac,
     data = data, weights = weights, subset = inde,
     family = binomial(link = "logit")), list(weights = weights,
     inde = inde)))
     }
     sigma <- sqrt(mean(residuals(gam2.1, type = "deviance")^2) +
     mean(imr[inde] * (imr[inde] + p.g1[inde])) * gam2.1$coef["imr"]^2)[[1]]
     names(sigma) <- "sigma"
     co <- (gam2.1$coef["imr"]/sigma)[[1]]
     if (margins[2] == "G") {
     sigma <- k <- (summary(gam2)$dispersion)^(-1)
     names(k) <- "sigma"
     }
     a.theta <- st.theta.star(start.theta, co, BivD)
     names(a.theta) <- "theta.star"
     if (margins[2] %in% c("N", "NB", "PIG", "WARING", "BB",
     "NBII")) {
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(sigma), a.theta)
     names(start.v)[length(start.v) - 1] <- "log.sigma"
     }
     if (margins[2] %in% c("ZIP", "ZAP", "ZIP2", "ZALG", "ZIBI",
     "ZABI")) {
     sigma.star.start1 <- sum(gam2.1$y == 0)/length(gam2.1$y)
     sigma.star.start1 <- ifelse(sigma.star.start1 > 0.999,
     0.999, sigma.star.start1)
     sigma.star.start1 <- ifelse(sigma.star.start1 < 0.001,
     0.001, sigma.star.start1)
     sigma.star.start1 <- qlogis(sigma.star.start1)
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], sigma.star.start1, a.theta)
     names(start.v)[length(start.v) - 1] <- "logi.sigma"
     }
     if (margins[2] == "G") {
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(k), a.theta)
     names(start.v)[length(start.v) - 1] <- "log.shape"
     }
     if (margins[2] %in% c("P", "YULE", "GEOM", "BI", "LG")) {
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], a.theta)
     }
     if (margins[2] %in% c("D")) {
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(sigma), qlogis(0.51), a.theta)
     names(start.v)[length(start.v) - 2] <- "log.sigma"
     names(start.v)[length(start.v) - 1] <- "logi.nu"
     }
     if (margins[2] %in% c("ZINBI")) {
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(sigma), qlogis(((sum(gam2$y ==
     0)/length(gam2$y)) + 0.01)/2), a.theta)
     names(start.v)[length(start.v) - 2] <- "log.sigma"
     names(start.v)[length(start.v) - 1] <- "logi.nu"
     }
     if (margins[2] %in% c("ZANBI")) {
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(sigma), qlogis(max(sum(gam2$y ==
     0)/length(gam2$y), 0.01)), a.theta)
     names(start.v)[length(start.v) - 2] <- "log.sigma"
     names(start.v)[length(start.v) - 1] <- "logi.nu"
     }
     if (margins[2] %in% c("ZIPIG")) {
     nu.star.start1 <- sum(gam2.1$y == 0)/length(gam2.1$y)
     nu.star.start1 <- ifelse(nu.star.start1 > 0.999,
     0.999, nu.star.start1)
     nu.star.start1 <- ifelse(nu.star.start1 < 0.001,
     0.001, nu.star.start1)
     nu.star.start1 <- qlogis(nu.star.start1)
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(sigma), nu.star.start1, a.theta)
     names(start.v)[length(start.v) - 2] <- "log.sigma"
     names(start.v)[length(start.v) - 1] <- "logi.nu"
     }
     if (margins[2] %in% c("ZIBB", "ZABB")) {
     nu.star.start1 <- sum(gam2.1$y == 0)/length(gam2.1$y)
     nu.star.start1 <- ifelse(nu.star.start1 > 0.999,
     0.999, nu.star.start1)
     nu.star.start1 <- ifelse(nu.star.start1 < 0.001,
     0.001, nu.star.start1)
     nu.star.start1 <- qlogis(nu.star.start1)
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(sigma), nu.star.start1, a.theta)
     names(start.v)[length(start.v) - 2] <- "log.sigma"
     names(start.v)[length(start.v) - 1] <- "logi.nu"
     }
     if (margins[2] == "S") {
     if (is.null(start.v))
     start.v <- c(coef(gam1), coef(gam2.1)[names(coef(gam2.1)) !=
     "imr"], log(sigma), -0.5, a.theta)
     names(start.v)[length(start.v) - 2] <- "log.sigma"
     names(start.v)[length(start.v) - 1] <- "nu"
     }
     start.v1 <- start.v[-c(1:X1.d2, length(start.v))]
     }
     else {
     start.v <- start.v <- c(coef(gam1), coef(gam2), coef(gam3),
     coef(gam4), coef(gam5))
     start.v1 <- start.v[-c(1:X1.d2, (length(start.v) - (llM -
     1)):length(start.v))]
     }
     if (missing(parscale))
     parscale <- 1
     if (margins[2] %in% c("N", "G"))
     funcop <- ghss
     else funcop <- ghssD
     respvec <- cbind(y1, y2)
     VC <- list(X1 = X1, X2 = X2, X3 = X3, X4 = X4, X5 = X5, X1.d2 = X1.d2,
     X2.d2 = X2.d2, X3.d2 = X3.d2, X4.d2 = X4.d2, X5.d2 = X5.d2,
     gp1 = gp1, gp2 = gp2, gp3 = gp3, gp4 = gp4, gp5 = gp5,
     n = n, bd = bd, l.sp1 = l.sp1, l.sp2 = l.sp2, l.sp3 = l.sp3,
     l.sp4 = l.sp4, l.sp5 = l.sp5, infl.fac = infl.fac, weights = weights,
     fp = fp, BivD = BivD, margins = margins, parscale = parscale)
     if (margins[2] %in% marg2D) {
     gamlss.fit <- fit.SemiParSampleSel(funcop = ghssDuniv,
     start.v = start.v1, dat = respvec, rinit = rinit,
     rmax = rmax, VC = VC, qu.mag = qu.mag1, sp = sp1,
     iterlimsp = iterlimsp, pr.tolsp = pr.tolsp, univariate = TRUE)
     start.v <- c(coef(gam1), gamlss.fit$fit$argument, a.theta)
     sp <- c(sp1, gamlss.fit$sp, sp.theta)
     gam2$coefficients <- gamlss.fit$fit$argument[1:X2.d2]
     gam2$Vp <- gamlss.fit$magpp$Vb[1:X2.d2, 1:X2.d2]
     gam2$sig2 <- 1
     gam2$edf <- gamlss.fit$magpp$edf[1:X2.d2]
     if (length(formula) > 3) {
     if (margins[2] %in% c("NB", "PIG", "D", "S", "BB",
     "NBII", "WARING", "ZABI", "ZIBI", "ZALG", "ZAP",
     "ZIP", "ZIP2", "ZIBB", "ZABB", "ZANBI", "ZINBI",
     "ZIPIG")) {
     gam3$coefficients <- gamlss.fit$fit$argument[X2.d2 +
     (1:X3.d2)]
     gam3$Vp <- gamlss.fit$magpp$Vb[X2.d2 + (1:X3.d2),
     X2.d2 + (1:X3.d2)]
     gam3$sig2 <- 1
     gam3$edf <- gamlss.fit$magpp$edf[X2.d2 + (1:X3.d2)]
     }
     if (margins[2] %in% c("D", "S", "ZIBB", "ZABB", "ZANBI",
     "ZINBI", "ZIPIG")) {
     gam4$coefficients <- gamlss.fit$fit$argument[X2.d2 +
     X3.d2 + (1:X4.d2)]
     gam4$Vp <- gamlss.fit$magpp$Vb[X2.d2 + X3.d2 +
     (1:X4.d2), X2.d2 + X3.d2 + (1:X4.d2)]
     gam4$sig2 <- 1
     gam4$edf <- gamlss.fit$magpp$edf[X2.d2 + X3.d2 +
     (1:X4.d2)]
     }
     }
     }
     fit <- fit.SemiParSampleSel(funcop = funcop, start.v = start.v,
     dat = respvec, rinit = rinit, rmax = rmax, VC = VC, qu.mag = qu.mag,
     sp = sp, iterlimsp = iterlimsp, pr.tolsp = pr.tolsp,
     univariate = FALSE)
     epsilon <- sqrt(.Machine$double.eps)
     He <- fit$fit$hessian
     He.eig <- eigen(He, symmetric = TRUE)
     if (min(He.eig$values) < epsilon)
     He.eig$values[which(He.eig$values < epsilon)] <- 1e-07
     Vb <- He.eig$vectors %*% tcrossprod(diag(1/He.eig$values),
     He.eig$vectors)
     if ((l.sp1 != 0 || l.sp2 != 0 || l.sp3 != 0 || l.sp4 != 0 ||
     l.sp5 != 0) && fp == FALSE) {
     HeSh <- He - fit$fit$S.h
     F <- Vb %*% HeSh
     F1 <- fit$magpp$edf1
     R <- fit$bs.mgfit$R
     Ve <- F %*% Vb
     }
     else {
     HeSh <- He
     Ve <- Vb
     F <- F1 <- diag(rep(1, dim(Vb)[1]))
     R <- fit$bs.mgfit$R
     }
     t.edf <- sum(diag(F))
     param <- fit$fit$argument[-c(1:X1.d2)]
     param <- param[1:X2.d2]
     X2s <- try(predict.gam(gam2, newdata = data, type = "lpmatrix"),
     silent = TRUE)
     if (class(X2s) == "try-error")
     stop("Check that the range of the covariate values of the selected sample is the same as that of the complete dataset.")
     eta2 <- X2s %*% param
     dimnames(fit$fit$hessian)[[1]] <- dimnames(fit$fit$hessian)[[2]] <- dimnames(Vb)[[1]] <- dimnames(Vb)[[2]] <- dimnames(HeSh)[[1]] <- dimnames(HeSh)[[2]] <- dimnames(F)[[1]] <- dimnames(F)[[2]] <- dimnames(He)[[1]] <- dimnames(He)[[2]] <- names(fit$fit$argument)
     if ((l.sp1 != 0 || l.sp2 != 0 || l.sp3 != 0 || l.sp4 != 0 ||
     l.sp5 != 0)) {
     edf <- edf1 <- list(0, 0, 0, 0, 0)
     for (i in 1:5) {
     if (i == 1) {
     mm <- l.sp1
     if (mm == 0)
     next
     }
     if (i == 2) {
     mm <- l.sp2
     if (mm == 0)
     next
     }
     if (i == 3) {
     mm <- l.sp3
     if (mm == 0)
     next
     }
     if (i == 4) {
     mm <- l.sp4
     if (mm == 0)
     next
     }
     if (i == 5) {
     mm <- l.sp5
     if (mm == 0)
     break
     }
     for (k in 1:mm) {
     if (i == 1) {
     gam <- gam1
     ind <- (gam$smooth[[k]]$first.para):(gam$smooth[[k]]$last.para)
     }
     if (i == 2) {
     gam <- gam2
     ind <- (gam$smooth[[k]]$first.para:gam$smooth[[k]]$last.para) +
     X1.d2
     }
     if (i == 3) {
     gam <- gam3
     ind <- (gam$smooth[[k]]$first.para:gam$smooth[[k]]$last.para) +
     X1.d2 + X2.d2
     }
     if (i == 4) {
     gam <- gam4
     ind <- (gam$smooth[[k]]$first.para:gam$smooth[[k]]$last.para) +
     X1.d2 + X2.d2 + X3.d2
     }
     if (i == 5) {
     gam <- gam5
     ind <- (gam$smooth[[k]]$first.para:gam$smooth[[k]]$last.para) +
     X1.d2 + X2.d2 + X3.d2 + X4.d2
     }
     edf[[i]][k] <- sum(diag(F)[ind])
     edf1[[i]][k] <- sum(F1[ind])
     }
     }
     if (l.sp1 != 0)
     names(edf[[1]]) <- names(edf1[[1]]) <- names(gam1$sp)
     if (l.sp2 != 0)
     names(edf[[2]]) <- names(edf1[[2]]) <- names(gam2$sp)
     if (l.sp3 != 0)
     names(edf[[3]]) <- names(edf1[[3]]) <- names(gam3$sp)
     if (l.sp4 != 0)
     names(edf[[4]]) <- names(edf1[[4]]) <- names(gam4$sp)
     if (l.sp5 != 0)
     names(edf[[5]]) <- names(edf1[[5]]) <- names(gam5$sp)
     }
     KendTau <- theta <- sigma <- k <- phi <- nu <- NULL
     if (length(formula) < 3) {
     theta <- theta.tau(BivD = BivD, theta.star = fit$fit$argument["theta.star"])
     KendTau <- theta[, 2]
     theta <- theta[, 1]
     names(theta) <- names(KendTau) <- NULL
     if (margins[2] %in% c("N", "NB", "PIG", "D", "S", "BB",
     "NBII", "WARING", "ZIBB", "ZABB", "ZANBI", "ZINBI",
     "ZIPIG")) {
     sigma <- exp(fit$fit$argument["log.sigma"])
     names(sigma) <- k <- NULL
     phi <- sigma^2
     }
     if (margins[2] %in% c("ZABI", "ZIBI", "ZALG", "ZAP",
     "ZIP", "ZIP2")) {
     sigma <- plogis(fit$fit$argument["logi.sigma"])
     names(sigma) <- k <- NULL
     phi <- sigma^2
     }
     if (margins[2] == "G") {
     sigma <- k <- exp(fit$fit$argument["log.shape"])
     names(k) <- names(sigma) <- NULL
     phi <- k^{
     -1
     }
     }
     if (margins[2] == "S") {
     nu <- fit$fit$argument["nu"]
     names(nu) <- NULL
     }
     if (margins[2] %in% c("D", "ZIBB", "ZABB", "ZANBI", "ZINBI",
     "ZIPIG")) {
     nu <- plogis(fit$fit$argument["logi.nu"])
     names(nu) <- NULL
     }
     theta.a <- theta
     KendTau.a <- KendTau
     sigma.a <- sigma
     phi.a <- phi
     nu.a <- nu
     }
     if (length(formula) > 2) {
     sigma <- k <- phi <- nu <- NULL
     etatheta <- fit$fit$etatheta
     etasqv <- fit$fit$etasqv
     etak <- fit$fit$etak
     etak <- fit$fit$etak
     theta <- theta.tau(BivD = BivD, theta.star = fit$fit$etatheta)
     KendTau <- theta[, 2]
     theta <- as.matrix(theta[, 1])
     if (margins[2] %in% c("N", "NB", "PIG", "D", "S", "BB",
     "NBII", "WARING", "ZIBB", "ZABB", "ZANBI", "ZINBI",
     "ZIPIG")) {
     sigma <- exp(fit$fit$etasqv)
     phi <- sigma^2
     }
     if (margins[2] %in% c("ZABI", "ZIBI", "ZALG", "ZAP",
     "ZIP", "ZIP2")) {
     sigma <- plogis(fit$fit$etasqv)
     phi <- sigma^2
     }
     if (margins[2] == "G") {
     sigma <- k <- exp(fit$fit$etak)
     phi <- k^{
     -1
     }
     }
     if (margins[2] == "S")
     nu <- fit$fit$etanu
     if (margins[2] %in% c("D", "ZIBB", "ZABB", "ZANBI", "ZINBI",
     "ZIPIG"))
     nu <- plogis(fit$fit$etanu)
     theta.a <- mean(as.numeric(theta))
     KendTau.a <- mean(as.numeric(KendTau))
     sigma.a <- mean(as.numeric(sigma))
     phi.a <- mean(as.numeric(phi))
     nu.a <- mean(as.numeric(nu))
     }
     L <- list(fit = fit$fit, gam1 = gam1, gam2 = gam2, gam3 = gam3,
     gam4 = gam4, gam5 = gam5, gam2.1 = gam2.1, coefficients = fit$fit$argument,
     weights = weights, sp = fit$sp, iter.if = fit$iter.if,
     iter.sp = fit$iter.sp, iter.inner = fit$iter.inner, start.v = start.v,
     formula = formula, phi = phi, sigma = sigma, nu = nu,
     theta = theta, tau = KendTau, n = n, n.sel = length(gam2$y),
     bd = bd, X1 = X1, X2 = X2, X1.d2 = X1.d2, X2.d2 = X2.d2,
     X3 = X3, X3.d2 = X3.d2, X4 = X4, X4.d2 = X4.d2, X5 = X5,
     X5.d2 = X5.d2, l.sp1 = l.sp1, l.sp2 = l.sp2, l.sp3 = l.sp3,
     l.sp4 = l.sp4, l.sp5 = l.sp5, He = He, HeSh = HeSh, Vb = Vb,
     F = F, BivD = BivD, margins = margins, t.edf = t.edf,
     bs.mgfit = fit$bs.mgfit, conv.sp = fit$conv.sp, wor.c = fit$wor.c,
     eta1 = fit$fit$eta1, eta2 = eta2, y1 = y1, y2 = y2, logLik = -fit$fit$l,
     fp = fp, gp1 = gp1, gp2 = gp2, gp3 = gp3, gp4 = gp4,
     gp5 = gp5, X2s = X2s, magpp = fit$magpp, sigma.a = sigma.a,
     phi.a = phi.a, theta.a = theta.a, nu.a = nu.a, tau.a = KendTau.a,
     gamlss.fit = gamlss.fit, edf = edf, edf11 = edf1, R = R,
     Ve = Ve)
     class(L) <- "SemiParSampleSel"
     L
    }
    <bytecode: 0xb92a600>
    <environment: namespace:SemiParSampleSel>
     --- function search by body ---
    Function SemiParSampleSel in namespace SemiParSampleSel has this body.
     ----------- END OF FAILURE REPORT --------------
    Fatal error: the condition has length > 1
Flavor: r-devel-linux-x86_64-debian-clang