diff --git a/R/mc_fast_forward.R b/R/mc_fast_forward.R new file mode 100644 index 0000000000000000000000000000000000000000..d8a92452887df8ddae071aef0e55df903ddf31bd --- /dev/null +++ b/R/mc_fast_forward.R @@ -0,0 +1,45 @@ +#' Fast forward selection for multivariate covariance generalized linear +#' models. +#' +#' @description Perform fast forward model selection using the score +#' information criterion. This function works only for univariate months. +#' +#' @param object an object representing a model of \code{mcglm} class. +#' @param scope a vector specyfing the covariate to be tested. +#' @param interaction Maximum number of covariates interacting. +#' @param penalty penalty term (default = 2). +#' @return The selected model. +#' @export + +mc_fast_forward <- function(object, scope, interaction = 1, + penalty =2, n_max = 10) { + if (interaction > 1) { + int_terms <- list() + for (i in 2:interaction) { + int_terms[[c(i-1)]] <- combn(length(scope), i) + } + } + fun_temp <- function(int_terms) { + output <- c() + for(i in 1:dim(int_terms)[2]) { + output[i] <- paste(scope[int_terms[,i]],collapse = "*") + } + return(output) + } + scope <- c(scope, do.call(c, lapply(int_terms, fun_temp))) + for (i in 1:n_max) { + sic <- mc_sic(object = object, scope = scope, + data = data, response = 1) + if (all(sic$SIC < 0)) break + cov_new <- as.numeric(rownames(sic[which.max(sic$SIC),])) + cov_enter <- scope[as.numeric(rownames(sic[which.max(sic$SIC),]))] + next_cov <- paste("~. +", cov_enter) + new_formula <- update.formula(object$linear_pred[[1]], next_cov) + temp_models <- mcglm(c(new_formula), matrix_pred = matrix_pred, + link = object$link, variance = object$variance, + covariance = object$covariance, data = data) + object <- temp_models + scope <- scope[-cov_new] + } + return(object) +} diff --git a/R/mc_sic.R b/R/mc_sic.R new file mode 100644 index 0000000000000000000000000000000000000000..4d91101bb9317aedc9ab16d9537cc5a30529ab38 --- /dev/null +++ b/R/mc_sic.R @@ -0,0 +1,49 @@ +#' Compute the score information criterion (SIC) for multivariate +#' covariance generalized linear models. +#' +#' @description Compute the SIC for McGLMS. +#' @param object an object representing a model of \code{mcglm} class. +#' @param scope a vector containing all covariate names to be tested. +#' @param data data frame containing the all variables envolved +#' @param penalty penalty term (default = 2). +#' @param response Indicate for which response variable SIC is computed. +#' @return A data frame with SIC values for each covariate in the scope +#' argument. +#' @export + +mc_sic <- function (object, scope, data, response, penalty = 2) { + SIC <- c() + df <- c() + for(i in 1:length(scope)){ + ini_formula <- object$linear_pred[[response]] + ext_formula <- as.formula(paste("~", paste(ini_formula[3], + scope[i], sep = "+"))) + md <- model.frame(object$linear_pred[[response]], data = data) + Y = model.response(md) + ini_beta <- coef(object, type = "beta", response = response)$Estimates + ext_X <- model.matrix(ext_formula, data = data) + n_beta <- dim(ext_X)[2] + n_ini_beta <- length(ini_beta) + ext_beta <- c(ini_beta, rep(0, n_beta - n_ini_beta)) + n_total_beta <- length(ext_beta) + mu_temp <- mc_link_function(beta = ext_beta, X = ext_X, offset = NULL, + link = object$link[[response]]) + score_temp <- mc_quasi_score(D = mu_temp$D, inv_C = object$inv_C, + y_vec = Y, mu_vec = mu_temp$mu) + S11 <- score_temp$Variability[1:n_ini_beta,1:n_ini_beta] + S22 <- score_temp$Variability[c(n_ini_beta+1):n_total_beta, + c(n_ini_beta+1):n_total_beta] + S12 <- score_temp$Variability[1:n_ini_beta, + c(n_ini_beta+1):n_total_beta] + S21 <- score_temp$Variability[c(n_ini_beta+1):n_total_beta, + 1:n_ini_beta] + VB <- S22 - S21 %*% solve(S11) %*% S12 + Tu <- t(score_temp$Score[c(n_ini_beta+1):n_total_beta])%*% + solve(VB)%*%score_temp$Score[c(n_ini_beta+1):n_total_beta] + df[i] <- n_beta - n_ini_beta + SIC[i] <- as.numeric(sqrt(Tu)) - penalty*df[i] + } + output <- data.frame("SIC" = SIC, "Covariance" = scope, "df" = df) + return(output) +} + diff --git a/man/mc_fast_forward.Rd b/man/mc_fast_forward.Rd new file mode 100644 index 0000000000000000000000000000000000000000..0e8f2cc46d705647daa8d89c9cdbae02d81d5e92 --- /dev/null +++ b/man/mc_fast_forward.Rd @@ -0,0 +1,26 @@ +% Generated by roxygen2 (4.1.1): do not edit by hand +% Please edit documentation in R/mc_fast_forward.R +\name{mc_fast_forward} +\alias{mc_fast_forward} +\title{Fast forward selection for multivariate covariance generalized linear +models.} +\usage{ +mc_fast_forward(object, scope, interaction = 1, penalty = 2, n_max = 10) +} +\arguments{ +\item{object}{an object representing a model of \code{mcglm} class.} + +\item{scope}{a vector specyfing the covariate to be tested.} + +\item{interaction}{Maximum number of covariates interacting.} + +\item{penalty}{penalty term (default = 2).} +} +\value{ +The selected model. +} +\description{ +Perform fast forward model selection using the score +information criterion. This function works only for univariate months. +} + diff --git a/man/mc_sic.Rd b/man/mc_sic.Rd new file mode 100644 index 0000000000000000000000000000000000000000..e88fc662558ee3878db5ef2caf0fbc8831cacb60 --- /dev/null +++ b/man/mc_sic.Rd @@ -0,0 +1,28 @@ +% Generated by roxygen2 (4.1.1): do not edit by hand +% Please edit documentation in R/mc_sic.R +\name{mc_sic} +\alias{mc_sic} +\title{Compute the score information criterion (SIC) for multivariate +covariance generalized linear models.} +\usage{ +mc_sic(object, scope, data, response, penalty = 2) +} +\arguments{ +\item{object}{an object representing a model of \code{mcglm} class.} + +\item{scope}{a vector containing all covariate names to be tested.} + +\item{data}{data frame containing the all variables envolved} + +\item{response}{Indicate for which response variable SIC is computed.} + +\item{penalty}{penalty term (default = 2).} +} +\value{ +A data frame with SIC values for each covariate in the scope +argument. +} +\description{ +Compute the SIC for McGLMS. +} +