From 795926bf83dfa00ec142500bc2968b5109f7f5cc Mon Sep 17 00:00:00 2001 From: wbonat <wbonat@gmail.com> Date: Fri, 23 Oct 2015 17:40:17 +0200 Subject: [PATCH] Fix small errors --- R/mc_bias_correct_std.R | 2 +- R/mc_influence.R | 27 +++++++++++++++------------ R/mc_initial_values.R | 6 +++--- man/mc_dfbetaOij.Rd | 31 +++++++++++++++++++++++++++++++ man/mc_influence.Rd | 20 ++++++++++++-------- man/mc_sic_covariance.Rd | 31 +++++++++++++++++++++++++++++++ 6 files changed, 93 insertions(+), 24 deletions(-) create mode 100644 man/mc_dfbetaOij.Rd create mode 100644 man/mc_sic_covariance.Rd diff --git a/R/mc_bias_correct_std.R b/R/mc_bias_correct_std.R index 718e9c1..c9c5c16 100644 --- a/R/mc_bias_correct_std.R +++ b/R/mc_bias_correct_std.R @@ -34,7 +34,7 @@ mc_bias_corrected_std <- function(object, id) { Hi[[i]] <- Di[[i]]%*%inv_M%*%t(Di[[i]])%*%inv_Ci[[i]] } H <- bdiag(Hi) - I <- Diagonal(dim(data)[1], 1) + I <- Diagonal(dim(temp_data)[1], 1) inv_IH <- solve(I - H) Vbeta = inv_M%*%(t(D)%*%object$inv_C%*%inv_IH%*%R%*%inv_IH%*% object$inv_C%*%D)%*%inv_M diff --git a/R/mc_influence.R b/R/mc_influence.R index 437c6b5..fc5e4c4 100644 --- a/R/mc_influence.R +++ b/R/mc_influence.R @@ -14,6 +14,7 @@ mc_influence <- function(object, id) { temp_data <- data.frame(object$residuals, id) + names(temp_data) <- c("res", "id") data_id <- split(temp_data, temp_data$id) n_id <- length(data_id) D <- bdiag(lapply(object$mu_list, function(x) x$D)) @@ -24,7 +25,7 @@ mc_influence <- function(object, id) { # Leverage for observations inv_M <- -object$inv_S_beta M <- solve(inv_M) - uni_id <- unique(id) + uni_id <- as.character(unique(id)) Hi <- list() Di <- list() inv_Ci <- list() @@ -42,7 +43,7 @@ mc_influence <- function(object, id) { # Leverage for observations Hobs <- diag(H) # Leverage for clusters - Hid <- tapply(Hobs, id, sum) + Hid <- tapply(Hobs, id, mean) # DFBETAs for clusters I <- Diagonal(dim(temp_data)[1], 1) inv_IH <- solve(I - H) @@ -60,7 +61,7 @@ mc_influence <- function(object, id) { } inv_IH_i = inv_IH[idTF,idTF] r_i <- object$residuals[idTF] - DFBETACi[[i]] <- t(inv_M%*%t(D_i)%*%inv_C_i%*%inv_IH_i%*%r_i) + DFBETACi[[i]] <- t(inv_M%*%t(D_i)%*%inv_C_i%*%inv_IH_i%*%r_i) DCLSi[i] <- as.numeric(DFBETACi[[i]]%*%M%*%t(DFBETACi[[i]]))/n_beta if(dim(inv_C_i)[1] == 1) {DFBETAOij[[i]] <- DFBETACi[[i]] } else { @@ -71,19 +72,21 @@ mc_influence <- function(object, id) { } DFBETACi <- do.call(rbind, DFBETACi) DFBETAOij <- do.call(rbind, DFBETAOij) - DCOij = do.call(rbind, apply(DFBETAOij, MARGIN = 1, - function(x,M,n_beta)t(x)%*%M%*%x/n_beta, - M = M, n_beta = n_beta)) - std.error <- coef(object, std.error = TRUE)$Std.error + DCOij <- apply(DFBETAOij, MARGIN = 1, + function(x,M,n_beta)as.numeric((t(x)%*%M%*%x)/n_beta), + M = M, n_beta = n_beta) + std.error <- coef(object, type = "beta", std.error = TRUE)$Std.error DFBETACi = data.frame(as.matrix(DFBETACi))/std.error DFBETAOij = data.frame(as.matrix(DFBETAOij))/std.error names(DFBETACi) <- do.call(c, object$beta_names) names(DFBETAOij) <- do.call(c, object$beta_names) - output_id <- data.frame(DFBETACi, "Leverage" = as.numeric(Hid), - "Cook" = as.numeric(DCLSi)) - output_obs <- data.frame(DFBETAOij, "Leverage" = as.numeric(Hobs), - "Cook" = as.numeric(DCOij)) - output <- list("Id" = output_id, "Obs" = output_obs) + DFBETACi$Leverage <- as.numeric(Hid) + DFBETACi$Cook <- as.numeric(DCLSi) + DFBETAOij$Leverage <- as.numeric(Hobs) + DFBETAOij$Cook <- as.numeric(DCOij) + DFBETAOij$Id <- rownames(object$data) + DFBETACi$Id <- uni_id + output <- list("Id" = DFBETACi, "Obs" = DFBETAOij) return(output) } diff --git a/R/mc_initial_values.R b/R/mc_initial_values.R index 53f2dce..bbbaa5c 100644 --- a/R/mc_initial_values.R +++ b/R/mc_initial_values.R @@ -84,13 +84,13 @@ mc_initial_values <- function(linear_pred, matrix_pred, link, variance, covarian tau_extra <- lapply(matrix_pred, length) list_initial$tau <- list() for (i in 1:n_resp) { - if (covariance == "identity") { + if (covariance[i] == "identity") { list_initial$tau[[i]] <- as.numeric(c(tau0_initial[[i]], rep(0, c(tau_extra[[i]] - 1)))) } - if (covariance == "inverse") { + if (covariance[i] == "inverse") { list_initial$tau[[i]] <- as.numeric(c(1/tau0_initial[[i]], rep(0, c(tau_extra[[i]] - 1)))) } - if (covariance == "expm") { + if (covariance[i] == "expm") { list_initial$tau[[i]] <- as.numeric(c(exp(tau0_initial[[i]]), rep(0.1, c(tau_extra[[i]] - 1)))) } } diff --git a/man/mc_dfbetaOij.Rd b/man/mc_dfbetaOij.Rd new file mode 100644 index 0000000..b7f16ea --- /dev/null +++ b/man/mc_dfbetaOij.Rd @@ -0,0 +1,31 @@ +% Generated by roxygen2 (4.1.1): do not edit by hand +% Please edit documentation in R/mc_influence.R +\name{mc_dfbetaOij} +\alias{mc_dfbetaOij} +\title{Influence measures for McGLMs} +\usage{ +mc_dfbetaOij(Di, Ci, inv_Ci, ri, inv_M) +} +\arguments{ +\item{Di}{D matrix for the cluster i.} + +\item{Ci}{C matrix for the cluster i.} + +\item{inv_Ci}{Inverse of C matrix for the cluster i.} + +\item{ri}{Residual vector for the cluster i.} + +\item{inv_M}{Inverse of variance/covariance of regression parameters.} +} +\value{ +Matrix with the DFBETA for observation in the cluster i. +} +\description{ +Compute influence measures for multivariate covariance +generalized linear models. Leverage, DFBETA and Cook's distance +for observations. Auxiliar function for \code{mc_influence}. +} +\author{ +Wagner Hugo Bonat, \email{wbonat@ufpr.br} +} + diff --git a/man/mc_influence.Rd b/man/mc_influence.Rd index 2886d73..6f80a60 100644 --- a/man/mc_influence.Rd +++ b/man/mc_influence.Rd @@ -2,22 +2,26 @@ % Please edit documentation in R/mc_influence.R \name{mc_influence} \alias{mc_influence} -\title{Influence measures} +\title{Influence measures for McGLMs} \usage{ mc_influence(object, id) } \arguments{ -\item{object}{An object of mcglm class.} +\item{object}{An object of \code{mcglm} class.} -\item{id}{a vector which identifies the clusters. The length and order of id should be the -same as the number of observations. Data are assumed to be sorted so that observations on a cluster -are contiguous rows for all entities in the formula.} +\item{id}{a vector which identifies the clusters. +The length and order of id should match with the number of +observations.} } \value{ -A matrix. Note that the function assumes that the data are in the correct order. +A list with influence measures for cluster and observations. } \description{ -Compute influence measures for multivariate covariance generalized linear models. -Leverage, DFBETA and Cook's distance for unit sample and observations. +Compute influence measures for multivariate covariance +generalized linear models. Leverage, DFBETA and Cook's distance +for unit sample and observations. +} +\author{ +Wagner Hugo Bonat, \email{wbonat@ufpr.br} } diff --git a/man/mc_sic_covariance.Rd b/man/mc_sic_covariance.Rd new file mode 100644 index 0000000..01fab7d --- /dev/null +++ b/man/mc_sic_covariance.Rd @@ -0,0 +1,31 @@ +% Generated by roxygen2 (4.1.1): do not edit by hand +% Please edit documentation in R/mc_sic_covariance.R +\name{mc_sic_covariance} +\alias{mc_sic_covariance} +\title{Compute the score information criterion (SIC) for multivariate +covariance generalized linear models.} +\usage{ +mc_sic_covariance(object, scope, idx, data, penalty = 2, response) +} +\arguments{ +\item{object}{an object representing a model of \code{mcglm} class.} + +\item{scope}{a list of matrices to be tested in the matrix linear +predictor.} + +\item{idx}{Indicator of matrices belong to the same effect.} + +\item{data}{data frame containing all variables envolved in the model.} + +\item{penalty}{penalty term (default = 2).} + +\item{response}{Indicate for which response variable SIC is computed.} +} +\value{ +A data frame with SIC values for each matrix in the scope +argument. +} +\description{ +Compute SIC for covariance parameters in McGLMS. +} + -- GitLab