Nothing
#' rigorous Lasso for Linear Models: Inference
#'
#' Estimation and inference of (low-dimensional) target coefficients in a high-dimensional linear model.
#'
#' The functions estimates (low-dimensional) target coefficients in a high-dimensional linear model.
#' An application is e.g. estimation of a treatment effect \eqn{\alpha_0} in a
#' setting of high-dimensional controls. The user can choose between the so-called post-double-selection method and partialling-out.
#' The idea of the double selection method is to select variables by Lasso regression of
#' the outcome variable on the control variables and the treatment variable on
#' the control variables. The final estimation is done by a regression of the
#' outcome on the treatment effect and the union of the selected variables in
#' the first two steps. In partialling-out first the effect of the regressors on the outcome and the treatment variable is taken out by Lasso and then a regression of the residuals is conducted. The resulting estimator for \eqn{\alpha_0} is normal
#' distributed which allows inference on the treatment effect. It presents a wrap function for \code{rlassoEffect}
#' which does inference for a single variable.
#'
#' @param x matrix of regressor variables serving as controls and potential
#' treatments. For \code{rlassoEffect} it contains only controls, for \code{rlassoEffects} both controls and potential treatments. For \code{rlassoEffects} it must have at least two columns.
#' @param y outcome variable (vector or matrix)
#' @param index vector of integers, logicals or variables names indicating the position (column) of
#' variables (integer case), logical vector of length of the variables (TRUE or FALSE) or the variable names of \code{x} which should be used for inference / as treatment variables.
#' @param method method for inference, either 'partialling out' (default) or 'double selection'.
#' @param I3 For the 'double selection'-method the logical vector \code{I3} has same length as the number of variables in \code{x};
#' indicates if variables (TRUE) should be included in any case to the model and they are exempt from selection. These variables should not be included in the \code{index}; hence the intersection with \code{index} must be the empty set.
#' In the case of partialling out it is ignored.
#' @param post logical, if post Lasso is conducted with default \code{TRUE}.
#' @param \dots parameters passed to the \code{rlasso} function.
#' @return The function returns an object of class \code{rlassoEffects} with the following entries: \item{coefficients}{vector with estimated
#' values of the coefficients for each selected variable} \item{se}{standard error (vector)}
#' \item{t}{t-statistic} \item{pval}{p-value} \item{samplesize}{sample size of the data set} \item{index}{index of the variables for which inference is performed}
#' @references A. Belloni, V. Chernozhukov, C. Hansen (2014). Inference on
#' treatment effects after selection among high-dimensional controls. The
#' Review of Economic Studies 81(2), 608-650.
#' @export
#' @rdname rlassoEffects
#' @examples
#' library(hdm); library(ggplot2)
#' set.seed(1)
#' n = 100 #sample size
#' p = 100 # number of variables
#' s = 3 # number of non-zero variables
#' X = matrix(rnorm(n*p), ncol=p)
#' colnames(X) <- paste("X", 1:p, sep="")
#' beta = c(rep(3,s), rep(0,p-s))
#' y = 1 + X%*%beta + rnorm(n)
#' data = data.frame(cbind(y,X))
#' colnames(data)[1] <- "y"
#' fm = paste("y ~", paste(colnames(X), collapse="+"))
#' fm = as.formula(fm)
#' lasso.effect = rlassoEffects(X, y, index=c(1,2,3,50))
#' lasso.effect = rlassoEffects(fm, I = ~ X1 + X2 + X3 + X50, data=data)
#' print(lasso.effect)
#' summary(lasso.effect)
#' confint(lasso.effect)
#' plot(lasso.effect)
# library(hdm)
# ## DGP
# n <- 250
# p <- 100
# px <- 10
# X <- matrix(rnorm(n*p), ncol=p)
# beta <- c(rep(2,px), rep(0,p-px))
# intercept <- 1
# y <- intercept + X %*% beta + rnorm(n)
# ## fit rlassoEffects object with inference on three variables
# rlassoEffects.reg <- rlassoEffects(x=X, y=y, index=c(1,7,20))
# ## methods
# summary(rlassoEffects.reg)
# confint(rlassoEffects.reg, level=0.9)
rlassoEffects <- function(x, ...)
UseMethod("rlassoEffects") # definition generic function
#' @export
#' @rdname rlassoEffects
rlassoEffects.default <- function(x, y, index = c(1:ncol(x)), method = "partialling out",
I3 = NULL, post = TRUE, ...) {
checkmate::checkChoice(method, c("partialling out", "double selection"))
if (is.logical(index)) {
k <- p1 <- sum(index)
} else {
k <- p1 <- length(index)
}
n <- dim(x)[1]
x <- as.matrix(x)
# preprocessing index numerischer Vektor
if (is.numeric(index)) {
index <- as.integer(index)
stopifnot(all(index <= ncol(x)) && length(index) <= ncol(x))
} else {
# logical Vektor
if (is.logical(index)) {
stopifnot(length(index) == ncol(x))
index <- which(index == T)
} else {
# character Vektor
if (is.character(index)) {
stopifnot(all(is.element(index, colnames(x))))
index <- which(is.element(colnames(x), index))
} else {
stop("argument index has an invalid type")
}
}
}
if (method == "double selection") {
# check validity of I3
I3ind <- which(I3 == T)
if (length(intersect(index, I3ind) != 0))
stop("I3 and index must not overlap!")
}
if (is.null(colnames(x)))
colnames(x) <- paste("V", 1:dim(x)[2], sep = "")
coefficients <- as.vector(rep(NA_real_, k))
se <- rep(NA_real_, k)
t <- rep(NA_real_, k)
pval <- rep(NA_real_, k)
lasso.regs <- vector("list", k)
reside <- matrix(NA, nrow = n, ncol = p1)
residv <- matrix(NA, nrow = n, ncol = p1)
coef.mat <- list()
selection.matrix <- matrix(NA, ncol = k, nrow = dim(x)[2])
names(coefficients) <- names(se) <- names(t) <- names(pval) <- names(lasso.regs) <- colnames(reside) <- colnames(residv) <- colnames(selection.matrix) <- colnames(x)[index]
rownames(selection.matrix) <- colnames(x)
for (i in 1:k) {
d <- x[, index[i], drop = FALSE]
Xt <- x[, -index[i], drop = FALSE]
I3m <- I3[-index[i]]
lasso.regs[[i]] <- try(col <- rlassoEffect(Xt, y, d, method = method,
I3 = I3m, post = post, ...), silent = TRUE)
if (is(lasso.regs[[i]], "try-error")) {
next
} else {
coefficients[i] <- col$alpha
se[i] <- col$se
t[i] <- col$t
pval[i] <- col$pval
reside[, i] <- col$residuals$epsilon
residv[, i] <- col$residuals$v
coef.mat[[i]] <- col$coefficients.reg
selection.matrix[-index[i],i] <- col$selection.index
}
}
#colnames(coef.mat) <- colnames(x)[index]
residuals <- list(e = reside, v = residv)
res <- list(coefficients = coefficients, se = se, t = t, pval = pval,
lasso.regs = lasso.regs, index = index, call = match.call(), samplesize = n,
residuals = residuals, coef.mat = coef.mat, selection.matrix = selection.matrix)
class(res) <- "rlassoEffects"
return(res)
}
#' @rdname rlassoEffects
#' @param formula An element of class \code{formula} specifying the linear model.
#' @param I An one-sided formula specifying the variables for which inference is conducted.
#' @param included One-sided formula of variables which should be included in any case (only for method="double selection").
#' @param data an optional data frame, list or environment (or object coercible by as.data.frame to a data frame) containing the variables in the model.
#' If not found in data, the variables are taken from environment(formula), typically the environment from which the function is called.
#' @export
rlassoEffects.formula <- function(formula, data, I, method = "partialling out",
included = NULL, post = TRUE, ...) {
cl <- match.call()
if (missing(data)) data <- environment(formula)
mf <- match.call(expand.dots = FALSE)
m <- match(c("formula", "data"), names(mf), 0L)
mf <- mf[c(1L, m)]
mf$drop.unused.levels <- TRUE
mf[[1L]] <- quote(stats::model.frame)
mf <- eval(mf, parent.frame())
mt <- attr(mf, "terms")
attr(mt, "intercept") <- 1
y <- model.response(mf, "numeric")
n <- length(y)
x <- model.matrix(mt, mf)[,-1, drop=FALSE]
cn <- attr(mt, "term.labels")
try(if (is.matrix(eval(parse(text=cn)))) cn <- colnames(eval(parse(text=cn))), silent=TRUE)
I.c <- check_variables(I, cn)
#I.c <- grep(cn[I.c],colnames(X))
I.c <- which(colnames(x) %in% cn[I.c])
I3 <- check_variables(included, cn)
#I3 <- grep(cn[I.c],colnames(X))
I3 <- which(colnames(x) %in% cn[I.c])
#if (length(intersect(I.c, I3) != 0))
# stop("I and included should not contain the same variables!")
est <- rlassoEffects(x, y, index = I.c, method = method,
I3 = I3, post = post, ...)
est$call <- cl
return(est)
}
#' @rdname rlassoEffects
#' @param d variable for which inference is conducted (treatment variable)
#' @export
rlassoEffect <- function(x, y, d, method = "double selection", I3 = NULL,
post = TRUE, ...) {
d <- as.matrix(d, ncol = 1)
y <- as.matrix(y, ncol = 1)
kx <- dim(x)[2]
n <- dim(x)[1]
if (is.null(colnames(d)))
colnames(d) <- "d1"
if (is.null(colnames(x)) & !is.null(x))
colnames(x) <- paste("x", 1:kx, sep = "")
if (method == "double selection") {
I1 <- rlasso(d ~ x, post = post, ...)$index
I2 <- rlasso(y ~ x, post = post, ...)$index
if (is.logical(I3)) {
I <- I1 + I2 + I3
I <- as.logical(I)
} else {
I <- I1 + I2
I <- as.logical(I)
names(I) <- union(names(I1),names(I2))
}
if (sum(I) == 0) {
I <- NULL
}
x <- cbind(d, x[, I, drop = FALSE])
reg1 <- lm(y ~ x)
alpha <- coef(reg1)[2]
names(alpha) <- colnames(d)
xi <- reg1$residuals * sqrt(n/(n - sum(I) - 1))
if (is.null(I)) {
reg2 <- lm(d ~ 1)
}
if (!is.null(I)) {
reg2 <- lm(d ~ x[, -1, drop = FALSE])
}
v <- reg2$residuals
var <- 1/n * 1/mean(v^2) * mean(v^2 * xi^2) * 1/mean(v^2)
se <- sqrt(var)
tval <- alpha/sqrt(var)
pval <- 2 * pnorm(-abs(tval))
if (is.null(I)) {
no.selected <- 1
} else {
no.selected <- 0
}
res <- list(epsilon = xi, v = v)
# results <- list(alpha=unname(alpha), se=drop(se), t=unname(tval),
# pval=unname(pval), no.selected=no.selected,
# coefficients=unname(alpha), coefficient=unname(alpha),
# coefficients.reg=coef(reg1), residuals=res, call=match.call(),
# samplesize=n)
se <- drop(se)
names(se) <- colnames(d)
results <- list(alpha = alpha, se = se, t = tval, pval = pval,
no.selected = no.selected, coefficients = alpha, coefficient = alpha,
coefficients.reg = coef(reg1), selection.index = I, residuals = res, call = match.call(),
samplesize = n)
}
if (method == "partialling out") {
reg1 <- rlasso(y ~ x, post = post, ...)
yr <- reg1$residuals
reg2 <- rlasso(d ~ x, post = post, ...)
dr <- reg2$residuals
reg3 <- lm(yr ~ dr)
alpha <- coef(reg3)[2]
var <- vcov(reg3)[2, 2]
se <- sqrt(var)
tval <- alpha/sqrt(var)
pval <- 2 * pnorm(-abs(tval))
res <- list(epsilon = reg3$residuals, v = dr)
I1 <- reg1$index
I2 <- reg2$index
I <- as.logical(I1 + I2)
names(I) <- union(names(I1),names(I2))
results <- list(alpha = unname(alpha), se = drop(se), t = unname(tval),
pval = unname(pval), coefficients = unname(alpha), coefficient = unname(alpha),
coefficients.reg = coef(reg1), selection.index = I, residuals = res, call = match.call(),
samplesize = n)
}
class(results) <- "rlassoEffects"
return(results)
}
################################################################################################################################### Methods for rlassoEffects
#' Methods for S3 object \code{rlassoEffects}
#'
#' Objects of class \code{rlassoEffects} are constructed by \code{rlassoEffects}.
#' \code{print.rlassoEffects} prints and displays some information about fitted \code{rlassoEffect} objects.
#' summary.rlassoEffects summarizes information of a fitted \code{rlassoEffect} object and is described at \code{\link{summary.rlassoEffects}}.
#' \code{confint.rlassoEffects} extracts the confidence intervals.
#' \code{plot.rlassoEffects} plots the estimates with confidence intervals.
#'
#' @param x an object of class \code{rlassoEffects}
#' @param digits significant digits in printout
#' @param ... arguments passed to the print function and other methods.
#' @rdname methods.rlassoEffects
#' @aliases methods.rlassoEffects print.rlassoEffects confint.rlassoEffects plot.rlassoEffects
#' @export
print.rlassoEffects <- function(x, digits = max(3L, getOption("digits") -
3L), ...) {
cat("\nCall:\n", paste(deparse(x$call), sep = "\n", collapse = "\n"),
"\n\n", sep = "")
if (length(coef(x))) {
cat("Coefficients:\n")
print.default(format(coef(x), digits = digits), print.gap = 2L,
quote = FALSE)
} else cat("No coefficients\n")
cat("\n")
invisible(coef(x))
}
#' @rdname methods.rlassoEffects
#' @param object an object of class \code{rlassoEffects}
#' @param parm a specification of which parameters are to be given confidence intervals among the variables for which inference was done, either a vector of numbers or a vector of names. If missing, all parameters are considered.
#' @param level confidence level required
#' @param joint logical, if \code{TRUE} joint confidence intervals are calculated.
#' @export
# confint.rlassoEffects <- function(object, parm, level = 0.95, joint = FALSE,
# ...) {
# B <- 500 # number of bootstrap repitions
# n <- object$samplesize
# k <- p1 <- length(object$coefficients)
# cf <- coef(object)
# pnames <- names(cf)
# if (missing(parm))
# parm <- pnames else if (is.numeric(parm))
# parm <- pnames[parm]
# if (!joint) {
# a <- (1 - level)/2
# a <- c(a, 1 - a)
# # fac <- qt(a, n-k)
# fac <- qnorm(a)
# pct <- format_perc(a, 3)
# ci <- array(NA, dim = c(length(parm), 2L), dimnames = list(parm,
# pct))
# ses <- object$se[parm]
# ci[] <- cf[parm] + ses %o% fac
# }
#
# if (joint) {
# phi <- object$residuals$e * object$residuals$v
# m <- 1/sqrt(colMeans(phi^2))
# phi <- t(t(phi)/m)
# sigma <- sqrt(colMeans(phi^2))
# sim <- vector("numeric", length = B)
# for (i in 1:B) {
# xi <- rnorm(n)
# phi_temp <- phi * xi
# Nstar <- 1/sqrt(n) * colSums(phi_temp)
# sim[i] <- max(abs(Nstar))
# }
# a <- (1 - level)/2
# ab <- c(a, 1 - a)
# pct <- format_perc(ab, 3)
# ci <- array(NA, dim = c(length(parm), 2L), dimnames = list(parm,
# pct))
# hatc <- quantile(sim, probs = 1 - a)
# ci[, 1] <- cf[parm] - hatc * 1/sqrt(n) * sigma
# ci[, 2] <- cf[parm] + hatc * 1/sqrt(n) * sigma
# }
# return(ci)
# }
confint.rlassoEffects <- function(object, parm, level = 0.95, joint = FALSE,
...) {
B <- 500 # number of bootstrap repitions
n <- object$samplesize
k <- p1 <- length(object$coefficients)
cf <- coef(object)
pnames <- names(cf)
if (missing(parm))
parm <- pnames else if (is.numeric(parm))
parm <- pnames[parm]
if (!joint) {
a <- (1 - level)/2
a <- c(a, 1 - a)
# fac <- qt(a, n-k)
fac <- qnorm(a)
pct <- format_perc(a, 3)
ci <- array(NA, dim = c(length(parm), 2L), dimnames = list(parm,
pct))
ses <- object$se[parm]
ci[] <- cf[parm] + ses %o% fac
}
if (joint) {
e <- object$residuals$e
v <- object$residuals$v
ev <- e*v
Ev2 <- colMeans(v^2)
Omegahat <- matrix(NA, ncol=k, nrow=k)
for (j in 1:k) {
for (l in 1:k) {
Omegahat[j,l] = Omegahat[l,j] = 1/(Ev2[j]*Ev2[l]) * mean(ev[,j]*ev[,l])
}
}
var <- diag(Omegahat)
names(var) <- names(cf)
sim <- vector("numeric", length = B)
for (i in 1:B) {
beta_i <- MASS::mvrnorm(mu = rep(0,k), Sigma=Omegahat/n)
sim[i] <- max(abs(beta_i/sqrt(var)))
}
a <- (1 - level) #not dividing by 2!
ab <- c(a/2, 1 - a/2)
pct <- format_perc(ab, 3)
ci <- array(NA, dim = c(length(parm), 2L), dimnames = list(parm,
pct))
hatc <- quantile(sim, probs = 1 - a)
ci[, 1] <- cf[parm] - hatc * sqrt(var[parm])
ci[, 2] <- cf[parm] + hatc * sqrt(var[parm])
}
return(ci)
}
#' @rdname methods.rlassoEffects
#' @export
#' @param main an overall title for the plot
#' @param xlab a title for the x axis
#' @param ylab a title for the y axis
#' @param xlim vector of length two giving lower and upper bound of x axis
plot.rlassoEffects <- function(x, joint=FALSE, level= 0.95, main = "", xlab = "coef", ylab = "",
xlim = NULL, ...) {
# generate ordered KI-matrix
coefmatrix <- cbind(summary(x)$coef, confint(x, joint = joint, level=level))[, c(1, 5, 6)]
if (is.null(dim(coefmatrix))) {
vec <- coefmatrix
coefmatrix <- matrix(vec, ncol = 3)
colnames(coefmatrix) <- names(vec)
}
rownames(coefmatrix) <- names(x$coefficients)
coefmatrix <- as.data.frame(coefmatrix)
coefmatrix <- cbind(rownames(coefmatrix), coefmatrix)
colnames(coefmatrix) <- c("names", "coef", "lower", "upper")
coefmatrix <- coefmatrix[order(abs(coefmatrix[, 2])), ]
col <- "#000099"
# scale
if (missing(xlim)) {
low <- min(coefmatrix[, -1])
up <- max(coefmatrix[, -1])
} else {
low <- xlim[1]
up <- xlim[2]
}
# generate points
plotobject <- ggplot2::ggplot(coefmatrix, ggplot2::aes(y = coef, x = factor(names,
levels = names))) + ggplot2::geom_point(colour = col) +
ggplot2::geom_hline(colour = col, ggplot2::aes(width = 0.1, h = 0, yintercept=0))
# generate errorbars (KIs)
plotobject <- plotobject + ggplot2::geom_errorbar(ymin = coefmatrix$lower,
ymax = coefmatrix$upper, colour = col)
# further graphic parameter
plotobject <- plotobject + ggplot2::ggtitle(main) + ggplot2::ylim(low,
up) + ggplot2::xlab(ylab) + ggplot2::ylab(xlab)
## invert x and y axis
#plotobject <- plotobject + ggplot2::coord_flip()
# layout
plotobject <- plotobject + ggplot2::theme_bw() + ggplot2::geom_blank() +
ggplot2::theme(panel.grid.major.x = ggplot2::element_blank(), panel.grid.minor.x = ggplot2::element_blank())
plotobject <- plotobject + ggplot2::scale_x_discrete(labels = abbreviate)
# plot
plotobject
}
################ Methods: summary
#' Summarizing rlassoEffects fits
#'
#' Summary method for class \code{rlassoEffects}
#'
#' Summary of objects of class \code{rlassoEffects}
#'
#' @param object an object of class \code{rlassoEffects}, usually a result of a call to \code{rlassoEffects}
#' @param ... further arguments passed to or from other methods.
#' @rdname summary.rlassoEffects
#' @export
summary.rlassoEffects <- function(object, ...) {
ans <- NULL
k <- length(object$coefficients)
table <- matrix(NA, ncol = 4, nrow = k)
rownames(table) <- names(object$coefficient)
colnames(table) <- c("Estimate.", "Std. Error", "t value", "Pr(>|t|)")
table[, 1] <- object$coefficients
table[, 2] <- object$se
table[, 3] <- object$t
table[, 4] <- object$pval
ans$coefficients <- table
ans$object <- object
class(ans) <- "summary.rlassoEffects"
return(ans)
}
#' @param x an object of class \code{summary.rlassoEffects}, usually a result of a call or \code{summary.rlassoEffects}
#' @param digits the number of significant digits to use when printing.
#' @method print summary.rlassoEffects
#' @rdname summary.rlassoEffects
#' @export
print.summary.rlassoEffects <- function(x, digits = max(3L, getOption("digits") -
3L), ...) {
if (length(coef(x$object))) {
k <- dim(x$coefficients)[1]
table <- x$coefficients
print("Estimates and significance testing of the effect of target variables")
printCoefmat(table, digits = digits, P.values = TRUE, has.Pvalue = TRUE)
cat("\n")
} else {
cat("No coefficients\n")
}
cat("\n")
invisible(table)
}
#' Coefficients from S3 objects \code{rlassoEffects}
#'
#' Method to extract coefficients from objects of class \code{rlassoEffects}
#'
#' Printing coefficients and selection matrix for S3 object \code{rlassoEffects}. Interpretation of entries in the selection matrix
#' \itemize{
#' \item \code{"-"} indicates a target variable,
#' \item \code{"x"} indicates that a variable has been selected with rlassoEffects (coefficient is different from zero),
#' \item \code{"."} indicates that a variable has been de-selected with rlassoEffects (coefficient is zero).
#' }
#'
#' @param object an object of class \code{rlassoEffects}, usually a result of a call \code{rlassoEffect} or \code{rlassoEffects}.
#' @param selection.matrix if TRUE, a selection matrix is returned that indicates the selected variables from each auxiliary regression.
#' Default is set to FALSE.
#' @param include.targets if FALSE (by default) only the selected control variables are listed in the \code{selection.matrix}. If set to TRUE,
#' the selection matrix will also indicate the selection of the target coefficients that are specified in the \code{rlassoEffects} call.
#' @param complete general option of the function \code{coef}.
#' @param ... further arguments passed to functions coef or print.
#' @export
#' @rdname coef.rlassoEffects
#' @examples
#' library(hdm)
#' set.seed(1)
#' n = 100 #sample size
#' p = 100 # number of variables
#' s = 7 # number of non-zero variables
#' X = matrix(rnorm(n*p), ncol=p)
#' colnames(X) <- paste("X", 1:p, sep="")
#' beta = c(rep(3,s), rep(0,p-s))
#' y = 1 + X%*%beta + rnorm(n)
#' data = data.frame(cbind(y,X))
#' colnames(data)[1] <- "y"
#' lasso.effect = rlassoEffects(X, y, index=c(1,2,3,50),
#' method = "double selection")
#' coef(lasso.effect) # standard use of coef() - without selection matrix
#' # with selection matrix
#' coef(lasso.effect, selection.matrix = TRUE)
#' # prettier output with print_coef (identical options as coef())
#' print_coef(lasso.effect, selection.matrix = TRUE)
coef.rlassoEffects <- function(object, complete = TRUE, selection.matrix = FALSE, include.targets = FALSE, ...) {
cf <- object$coefficients
if (selection.matrix == TRUE) {
mat <- object$selection.matrix
if (is.null(mat)) {
mat <- cbind(object$selection.index)
dmat2 <- dim(mat)[2]
rnames <- rownames(mat)
targetindx <- stats::complete.cases(mat)
}
else {
dmat2 <- dim(mat)[2]
rnames <- rownames(mat)
targetindx <- stats::complete.cases(mat)
mat <- cbind(mat, as.logical(apply(mat, 1, sum)))
colnames(mat)[dim(mat)[2]] <- "global"
}
if (include.targets == FALSE) {
mat <- mat[targetindx, , drop = FALSE]
rnames <- rownames(mat)
}
else{
mat <- rbind(mat[targetindx == FALSE, , drop = FALSE], mat[targetindx, , drop = FALSE])
rnames <- rownames(mat)
}
mat <- rbind(mat, apply(mat, 2, sum, na.rm = TRUE))
mat <- apply(mat, 2, function(x) gsub(1, "x", x))
mat <- apply(mat, 2, function(x) gsub(0, ".", x))
mat[is.na(mat)] <- "-"
rownames(mat) <- c(rnames, "sum")
if (complete) {
coef <- list(cf = cf, selection.matrix = mat)
return(coef)
}
else {
coef <- list(cf = cf[!is.na(cf)], selection.matrix = mat)
return(coef)
}
}
else {
if (complete) {
return(cf)
}
else {
return(cf[!is.na(cf)])
}
}
}
#' Printing coefficients from S3 objects \code{rlassoEffects}
#'
#' Printing coefficients for class \code{rlassoEffects}
#'
#' Printing coefficients and selection matrix for S3 object \code{rlassoEffects}
#'
#' @param x an object of class \code{rlassoEffects}, usually a result of a call \code{rlassoEffect} or \code{rlassoEffects}.
#' @param selection.matrix if TRUE, a selection matrix is returned that indicates the selected variables from each auxiliary regression.
#' Default is set to FALSE.
#' @param include.targets if FALSE (by default) only the selected control variables are listed in the \code{selection.matrix}. If set to TRUE,
#' the selection matrix will also indicate the selection of the target coefficients that are specified in the \code{rlassoEffects} call.
#' @param complete general option of the function \code{coef}.
#' @param ... further arguments passed to functions coef or print.
#' @rdname print_coef
#' @aliases print_coef.rlassoEffects
#' @export
#' @examples
#' library(hdm)
#' set.seed(1)
#' n = 100 #sample size
#' p = 100 # number of variables
#' s = 7 # number of non-zero variables
#' X = matrix(rnorm(n*p), ncol=p)
#' colnames(X) <- paste("X", 1:p, sep="")
#' beta = c(rep(3,s), rep(0,p-s))
#' y = 1 + X%*%beta + rnorm(n)
#' data = data.frame(cbind(y,X))
#' colnames(data)[1] <- "y"
#' lasso.effect = rlassoEffects(X, y, index=c(1,2,3,50),
#' method = "double selection")
#' # without target coefficient estimates
#' print_coef(lasso.effect, selection.matrix = TRUE)
#' # with target coefficient estimates
#' print_coef(lasso.effect, selection.matrix = TRUE, targets = TRUE)
print_coef <- function(x, ...){
UseMethod("print_coef")
}
#' @rdname print_coef
#' @export
print_coef.rlassoEffects <- function(x, complete = TRUE, selection.matrix = FALSE, include.targets = TRUE, ...) {
checkmate::check_class(x, "rlassoEffects")
if (selection.matrix == FALSE) {
cat("\n")
print("Estimated target coefficients")
print(coef(x), complete = complete, ...)
cat("\n")
}
else {
sel.mat <- coef(x, selection.matrix = selection.matrix, include.targets = include.targets, complete = complete, ...)
cat("\n")
print("Estimated target coefficients")
print(sel.mat$cf)
cat("\n")
print("Matrix with selection index from auxiliary regressions")
cat("\n")
print(sel.mat$selection.matrix)
cat("_ _ _ \n")
print("'-' indicates a target variable; ")
print("'x' indicates that a variable has been selected with rlassoEffects (coefficient is different from zero);")
print("'.' indicates that a variable has been de-selected with rlassoEffects (coefficient is zero).")
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.