R/model_stat.R

Defines functions model_stat

#' @export
model_stat <- function(fit, lag = FALSE, data, y, X) {
    if (any(class(fit$model) == "ARIMA")) {
        k <- length(fit$model$coef)
        n <- fit$model$nobs
    } else {
        k <- fit$model$rank
        n <- NROW(fit$model$x)
    }
    df <- n - k
    if (any(class(fit$model) == "ARIMA")) {
        rss <- sum((fit$fitted - mean(y))^2)
        ess <- sum(fit$residuals^2)
    } else {
        rss <- sum((fit$model$fitted.values - mean(y))^2)
        ess <- sum(fit$model$residuals^2)
    }
    tss <- rss + ess
    r_squared_adj <- 1 - (ess / (df)) / (tss / (n - 1))
    if (k != 0) {
        fisher_test <- (r_squared_adj / (k - 1)) / ((1 - r_squared_adj) / df)
        p_value <- pf(q = fisher_test, df1 = k - 1, df2 = df, lower.tail = FALSE)
        coefficients <- coeftest(x = fit$model)
        if (!is.null(X)) {
            tryCatch(
                expr = rownames(coefficients)[
                    (length(names(fit$model$coef)) - NCOL(X) + 1):
                        (length(names(fit$model$coef)))
                    ] <- if (lag) {
                        paste0(names(data[-1]), paste0(" lag ", x_lag))
                    } else {
                        colnames(data[-1])
                    },
                error = function(e) {
                    cat(
                        "Message: names of the factors for",
                        paste0("'", names(data)[1], "'"),
                        "were not replaced.\n"
                    )
                }
            )
        }
        coefficients <- matrix(
            data = round(x = coefficients, digits = 4),
            nrow = k, ncol = 4, byrow = FALSE,
            dimnames = list(
                rownames(coefficients),
                colnames(coefficients)
            )
        )
    } else {
        fisher_test <- NULL
        p_value <- NULL
        coefficients <- NULL
    }
    other_stat <- round(
        x = c(
            `N` = n, `df` = df,
            `R^2 adj` = r_squared_adj,
            `F value` = fisher_test, `Pr(>|F|)` = p_value,
            AIC = fit$model$aicc, `sigma^2` = fit$model$sigma2
        ),
        digits = 4
    )
    return(
        list(
            `Forecast method` = fit$method,
            Model = other_stat,
            Coefficients = coefficients
        )
    )
}
faganok/scenario documentation built on Nov. 28, 2017, 4:06 p.m.