R/lm.R

#' WLS as Feasible GLS
#' @author Kazumasa Hirata
#' @description Estimate with Weighted Least Squares as feasible Generalized Least Squares
#' @param lm lm
#' @param logarithmic logical
#' @param aux.x character
#' @return list
#' @importFrom dplyr %>% count mutate select
#' @export
#'
lm.WLS.as.FGLS <- function(lm, logarithmic, aux.x = c("x", "pred+predsq")) {
    # get residuals from lm
    residsq <- residuals(lm)^2

    # decide regressors for auxiliary regression
    aux.x <- match.arg(aux.x)
    aux.data <- NULL
    if (aux.x == "x") {
        # regressor
        aux.x <- attr(terms(lm), "term.labels")
        aux.data <- dplyr::select(model.frame(lm), aux.x)
    } else if (aux.x == "pred+predsq") {
        # regressor
        pred <- predict(lm)
        aux.data <- data.frame(pred = pred, predsq = pred^2)
    }

    aux.lm <-NULL
    h <- NULL
    if (logarithmic) {
        # apply logarithmic transformation to resid
        aux.data %>% dplyr::mutate(lresidsq = log(residsq)) -> aux.data
        aux.lm <- lm(lresidsq ~ ., data = aux.data)

        # get predictions
        aux.pred <- predict(aux.lm)

        # apply exponential transformation
        h <- exp(aux.pred)
    } else {
        aux.data %>% dplyr::mutate(residsq = residsq) -> aux.data
        aux.lm <- lm(residsq ~ ., data = aux.data)

        # get predictions
        h <- predict(aux.lm)
    }

    if (length(h[h < 0]) > 0) {
        print("The calculated weights contain negative values and cannot execute WLS.")
        return(list(
            weight = h,
            weight.neg.n = length(h[h < 0])
        ))
    } else {
        return(list(
            weight = h,
            fit = lm(formula = eval(lm$call[[2]]), weights = 1 / h, data = model.frame(lm))
        ))
    }
}

#' Two Stage Least Squares (TSLS)
#' @author Kazumasa Hirata
#' @description Estimate model with Two Stage Least Squares with Instrumental Variables.
#' @param lm lm
#' @return list
#' @importFrom dplyr %>% count mutate select
#' @export
#'
lm.TSLS <- function() {

}
sheeputech/econometrics documentation built on June 18, 2019, 7:33 a.m.