Skip to content
Snippets Groups Projects
Commit 01595a21 authored by Wagner Hugo Bonat's avatar Wagner Hugo Bonat
Browse files

Merge branch 'issue#6' into 'devel'

Issue#6: Review and organize S3 methods

This branch:
  - Review all generic S3 methods, those functions that has `.mcglm` sulfix.
  - Modifications predominantly were apply 72 column rule, add `@author`, `@method`, `@title` and `@name` in the documentation, add `\code{}` to occurrences of `mcglm`.
  - Ocurrences of `paste( , sep = "")` where replaced by `paste0( )`.
  - All the S3 generic functions were placed in one file named `mc_S3_methods.R`. Methods in general are small functions and they aren't the core of the package. Keep them in the same file reduces the number of files, is easier to find. Modifications in the function that produce a `mcglm` will, potentially, affect all generic methods, so the corrections must be applied in just one file.
  - Individual files containing one S3 methods were deleted.
  - The `qic.mcglm` where not placed in `mc_S3_methods.R`. Also, no `@method` were included in documentation. Its was treated as a regular function, not a method as the others. See CMT 
e406660a for details.
  - Put `README.Rmd` and `Examples/` in `.Rbuildignore`.

In the building process, there was a warning/note about the absence of the `ahs` object.

See merge request !4
parents a66e44be 94c33bed
Branches
No related tags found
No related merge requests found
......@@ -6,3 +6,5 @@
.Rd2pdf5516
.Rd2pdf*
buildPkg.R
README.Rmd
Examples/
......@@ -17,7 +17,7 @@ Suggests:
plyr
Imports:
Matrix,
asserthat
assertthat
License: GPL-2
LazyData: TRUE
URL: http://git.leg.ufpr.br/wbonat/mcglm
......
#' @title ANOVA method for McGLMs.
#' @name anova.mcglm
#'
#' @description ANOVA method for object of class McGLMS.
#'
#' @param object an object of class \code{mcglm}, usually, a result of a
#' call to \code{mcglm()}.
#' @param ... additional arguments affecting the summary produced. Note
#' that there is no extra options for mcglm object class.
#'
#' @return A \code{data.frame} with Chi-square statistic to test the
#' null hypothesis of a parameter, or a set of parameters, be
#' zero. The Wald test based on the observed covariance matrix of
#' the parameters is used.
#'
#' @method anova mcglm
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#' @export
anova.mcglm <- function(object, ...) {
n_resp <- length(object$mu_list)
n_beta <- lapply(object$list_X, ncol)
idx.list <- list()
for (i in 1:n_resp) {
idx.list[[i]] <- rep(i, n_beta[i])
}
vv <- vcov(object)
n_par <- dim(vv)[1]
idx.vec <- do.call(c, idx.list)
n_cov <- n_par - length(idx.vec)
idx.vec <- c(idx.vec, rep(0, n_cov))
temp.vcov <- list()
temp.beta <- list()
for (i in 1:n_resp) {
idx.id = idx.vec == i
temp.vcov[[i]] <- vv[idx.id, idx.id]
temp.beta[[i]] <-
coef(object, type = "beta", response = i)$Estimates
}
saida <- list()
for (i in 1:n_resp) {
idx <- attr(object$list_X[[i]], "assign")
names <- colnames(object$list_X[[i]])
if (names[1] == "(Intercept)") {
idx <- idx[-1]
names <- names[-1]
temp.beta[[i]] <- temp.beta[[i]][-1]
temp.vcov[[i]] <- temp.vcov[[i]][-1, -1]
}
n_terms <- length(unique(idx))
X2.resp <- list()
for (j in 1:n_terms) {
idx.TF <- idx == j
temp <- as.numeric(
t(temp.beta[[i]][idx.TF]) %*%
solve(as.matrix(temp.vcov[[i]])[idx.TF, idx.TF]) %*%
temp.beta[[i]][idx.TF])
nbeta.test <- length(temp.beta[[i]][idx.TF])
X2.resp[[j]] <- data.frame(
Covariate = names[idx.TF][1],
Chi.Square = round(temp, 4),
Df = nbeta.test,
p.value = round(pchisq(temp, nbeta.test,
lower.tail = FALSE), 4))
}
saida[[i]] <- do.call(rbind, X2.resp)
}
cat("Wald test for fixed effects\n")
return(saida)
}
#' @title Extract model coefficients for mcglm class
#' @name coef.mcglm
#'
#' @description \code{coef.mcglm} is a function which extracts model
#' coefficients from objects of \code{mcglm} class.
#'
#' @param object An object of \code{mcglm} class.
#' @param std.error Logical. If \code{TRUE} returns the standard errors
#' of the estimates. Default is \code{FALSE}.
#' @param response A numeric vector specyfing for which response
#' variables the coefficients should be returned.
#' @param type A string vector (can be 1 element length) specyfing which
#' coefficients should be returned. Options are \code{"beta"},
#' \code{"tau"}, \code{"power"}, \code{"tau"} and
#' \code{"correlation"}.
#' @param ... additional arguments affecting the summary produced. Note
#' that there is no extra options for mcglm object class.
#'
#' @return A \code{data.frame} with parameters names, estimates,
#' response number and parameters type.
#'
#' @method coef mcglm
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#' @export
coef.mcglm <- function(object, std.error = FALSE,
response = c(NA, 1:length(object$beta_names)),
type = c("beta", "tau", "power", "correlation"),
...) {
n_resp <- length(object$beta_names)
cod_beta <- list()
cod_power <- list()
cod_tau <- list()
type_beta <- list()
type_power <- list()
type_tau <- list()
resp_beta <- list()
resp_power <- list()
resp_tau <- list()
response_for <- 1:n_resp
for (i in response_for) {
cod_beta[[i]] <- paste0(
paste0("beta", i), 0:c(object$Information$n_betas[[i]] - 1))
type_beta[[i]] <- rep("beta", length(cod_beta[[i]]))
resp_beta[[i]] <- rep(response_for[i], length(cod_beta[[i]]))
if (object$Information$n_power[[i]] != 0 |
object$power_fixed[[i]] == FALSE) {
cod_power[[i]] <- paste0(
paste0("power", i), 1:object$Information$n_power[[i]])
type_power[[i]] <- rep("power",
length(cod_power[[i]]))
resp_power[[i]] <- rep(response_for[i],
length(cod_power[[i]]))
}
if (object$Information$n_power[[i]] == 0) {
cod_power[[i]] <- rep(1, 0)
type_power[[i]] <- rep(1, 0)
resp_power[[i]] <- rep(1, 0)
}
cod_tau[[i]] <- paste0(
paste0("tau", i), 1:object$Information$n_tau[[i]])
type_tau[[i]] <- rep("tau", length(cod_tau[[i]]))
resp_tau[[i]] <- rep(response_for[i], length(cod_tau[[i]]))
}
rho_names <- c()
if (n_resp != 1) {
combination <- combn(n_resp, 2)
for (i in 1:dim(combination)[2]) {
rho_names[i] <- paste0(
paste0("rho", combination[1, i]), combination[2, i])
}
}
type_rho <- rep("correlation", length(rho_names))
resp_rho <- rep(NA, length(rho_names))
cod <- c(do.call(c, cod_beta), rho_names,
do.call(c, Map(c, cod_tau)))
type_cod <- c(do.call(c, type_beta), type_rho,
do.call(c, Map(c, type_tau)))
response_cod <- c(do.call(c, resp_beta), resp_rho,
do.call(c, Map(c, resp_tau)))
if (length(cod_power) != 0) {
cod <- c(do.call(c, cod_beta), rho_names,
do.call(c, Map(c, cod_power, cod_tau)))
type_cod <- c(do.call(c, type_beta), type_rho,
do.call(c, Map(c, type_power, type_tau)))
response_cod <- c(do.call(c, resp_beta), resp_rho,
do.call(c, Map(c, resp_power, resp_tau)))
}
Estimates <- c(object$Regression, object$Covariance)
coef_temp <- data.frame(
Estimates = Estimates,
Parameters = cod,
Type = type_cod,
Response = response_cod)
if (std.error == TRUE) {
coef_temp <- data.frame(
Estimates = Estimates,
Std.error = sqrt(diag(object$vcov)),
Parameters = cod, Type = type_cod,
Response = response_cod)
}
output <- coef_temp[
which(coef_temp$Response %in% response &
coef_temp$Type %in% type), ]
return(output)
}
#' @title Confidence Intervals for mcglm
#' @name confint.mcglm
#'
#' @description Computes confidence intervals for parameters in a fitted
#' \code{mcglm} model.
#'
#' @param object a fitted \code{mcglm} object.
#' @param parm a specification of which parameters are to be given
#' confidence intervals, either a vector of number or a vector of
#' strings. If missing, all parameters are considered.
#' @param level the nominal confidence level.
#' @param ... additional arguments affecting the confidence interval
#' produced. Note that there is no extra options for \code{mcglm}
#' object class.
#'
#' @return A \code{data.frame} with confidence intervals, parameters
#' names, response number and parameters type.
#'
#' @method confint mcglm
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#' @export
confint.mcglm <- function(object, parm, level = 0.95, ...) {
temp <- coef(object, std.error = TRUE)
if (missing(parm)) {
parm <- 1:length(temp$Estimates)
}
a <- (1 - level)/2
a <- c(a, 1 - a)
fac <- qnorm(a)
ci <- temp$Estimates + temp$Std.error %o% fac
colnames(ci) <- paste0(format(a, 2), "%")
rownames(ci) <- temp$Parameters
return(ci[parm])
}
#' @title Extract Model Fitted Values of McGLM
#' @name fitted.mcglm
#'
#' @description Extract fitted values for objects of \code{mcglm} class.
#'
#' @param object An object of \code{mcglm} class.
#' @param ... additional arguments affecting the summary produced. Note
#' that there is no extra options for \code{mcglm} object class.
#'
#' @return Depending on the number of response variable, the function
#' \code{fitted.mcglm} returns a vector (univariate models) or a
#' matrix (multivariate models) of fitted values.
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#'
#' @method fitted mcglm
#' @export
fitted.mcglm <- function(object, ...) {
n_resp <- length(object$beta_names)
output <- Matrix(object$fitted, ncol = n_resp, nrow = object$n_obs)
return(output)
}
#' @title Default Multivariate Covariance Generalized Linear models plotting
#' @name plot.mcglm
#'
#' @description takes a fitted \code{mcglm} object and do plots based on
#' residuals, influence diagnostic measures and algorithm check.
#'
#' @param x a fitted \code{mcglm} object.
#' @param type Specify which graphical analysis will be performed.
#' Options are: \code{"residuals"}, \code{"influence"} and
#' \code{"algorithm"}.
#' @param ... additional arguments affecting the plot produced. Note
#' that there is no extra options for mcglm object class.
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#'
#' @method plot mcglm
#' @export
plot.mcglm <- function(x, type = "residuals", ...) {
object = x
n_resp <- length(object$beta_names)
if (type == "residuals") {
par(mar = c(2.6, 2.5, 0.1, 0.1),
mgp = c(1.6, 0.6, 0),
mfrow = c(2, n_resp))
for (i in 1:n_resp) {
res <- residuals(object, type = "pearson")[, i]
fit_values <- fitted(object)[, i]
plot(res ~ fit_values,
ylab = "Pearson residuals",
xlab = "Fitted values")
temp <- loess.smooth(fitted(object)[, i],
residuals(object, type = "pearson")[, i])
lines(temp$x, temp$y)
qqnorm(res)
qqline(res)
}
}
if (type == "algorithm") {
n_iter <- length(na.exclude(object$IterationCovariance[, 1]))
par(mar = c(2.6, 2.5, 0.1, 0.1),
mgp = c(1.6, 0.6, 0),
mfrow = c(2, 2))
matplot(object$IterationRegression[1:c(n_iter + 5), ],
type = "l", lty = 2,
ylab = "Regression", xlab = "Iterations")
matplot(object$IterationCovariance[1:c(n_iter + 5), ],
type = "l", lty = 2,
ylab = "Covariance", xlab = "Iterations")
matplot(object$ScoreRegression[1:c(n_iter + 5), ],
type = "l", lty = 2,
ylab = "Quasi-score Regression", xlab = "Iterations")
matplot(object$ScoreCovariance[1:c(n_iter + 5), ],
type = "l", lty = 2,
ylab = "Quasi-score Covariance", xlab = "Iterations")
}
if (type == "partial_residuals") {
list_beta <- mc_updateBeta(list_initial = object$list_initial,
betas = object$Regression,
n_resp = n_resp,
information = object$Information)
comp_X <- list()
for (i in 1:n_resp) {
comp_X[[i]] <- as.matrix(object$list_X[[i]]) *
as.numeric(list_beta$regression[[i]])
}
for (i in 1:n_resp) {
res <- residuals(object, type = "pearson")[, i]
dev.new()
n_cov <- dim(comp_X[[i]])[2]
par(mar = c(2.6, 2.5, 0.5, 0.5),
mgp = c(1.6, 0.6, 0),
mfrow = c(1, c(n_cov - 1)))
for (j in 2:n_cov) {
p1 <- comp_X[[i]][, j] + res
plot(p1 ~ object$list_X[[i]][, j],
xlab = object$beta_names[[i]][j],
ylab = "Partial residuals ")
}
}
}
}
#' @title Print method for Multivariate Covariance Generalized Linear
#' Model
#' @name print.mcglm
#'
#' @description The default print method for a \code{mcglm} object.
#'
#' @param x fitted model objects of class mcglm as produced by mcglm().
#' @param ... further arguments passed to or from other methods.
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#'
#' @rdname print.mcglm
#' @method print mcglm
#' @export
print.mcglm <- function(x, ...) {
object <- x
n_resp <- length(object$beta_names)
regression <- mc_updateBeta(list_initial = list(),
betas = object$Regression,
information = object$Information,
n_resp = n_resp)
for (i in 1:n_resp) {
cat("Call: ")
print(object$linear_pred[[i]])
cat("\n")
cat("Link function:", object$link[[i]])
cat("\n")
cat("Variance function:", object$variance[[i]])
cat("\n")
cat("Covariance function:", object$covariance[[i]])
cat("\n")
names(regression[[1]][[i]]) <- object$beta_names[[i]]
cat("Regression:\n")
print(regression[[1]][[i]])
cat("\n")
cat("Dispersion:\n")
tau_temp <- coef(object, response = i, type = "tau")$Estimate
names(tau_temp) <- rep("", length(tau_temp))
print(tau_temp)
cat("\n")
power_temp <- coef(object, response = i, type = "power")$Estimate
if (length(power_temp) != 0) {
names(power_temp) <- ""
cat("Power:\n")
print(power_temp)
cat("\n")
}
}
}
#' @title Residuals for Multivariate Covariance Generalized Linear
#' Models (McGLM)
#' @name residuals.mcglm
#'
#' @description Compute residuals based on fitting \code{mcglm} models.
#'
#' @param object An of class \code{mcglm}, typically the result of a
#' call to \code{mcglm}.
#' @param type the type of residuals which should be returned. The
#' alternatives are: \code{"raw"} (default), \code{"pearson"} and
#' \code{"standardized"}.
#' @param ... additional arguments affecting the residuals
#' produced. Note that there is no extra options for mcglm object
#' class.
#'
#' @return Depending on the number of response variable the function
#' \code{residuals.mcglm} returns a vector (univariate models) or a
#' matrix (multivariate models) of residuals values.
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#'
#' @method residuals mcglm
#' @export
residuals.mcglm <- function(object, type = "raw", ...) {
n_resp <- length(object$beta_names)
output <- Matrix(object$residuals,
ncol = n_resp, nrow = object$n_obs)
if (type == "standardized") {
output <- Matrix(
as.numeric(object$residuals %*% chol(object$inv_C)),
ncol = n_resp, nrow = object$n_obs)
}
if (type == "pearson") {
output <- Matrix(
as.numeric(object$residuals/sqrt(diag(object$C))),
ncol = n_resp, nrow = object$n_obs)
}
return(output)
}
#' @title Summarizing Multivariate Covariance Generalized Linear Models
#' fits.
#' @name summary.mcglm
#'
#' @description Summary for McGLMs objects.
#'
#' @param object an object of class \code{mcglm}, usually, a result of a
#' call to \code{mcglm}.
#' @param ... additional arguments affecting the summary produced. Note
#' the there is no extra options for mcglm object class.
#'
#' @return Print an \code{mcglm} object.
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#'
#' @method summary mcglm
#' @export
summary.mcglm <- function(object, ...) {
n_resp <- length(object$beta_names)
output <- list()
for (i in 1:n_resp) {
cat("Call: ")
print(object$linear_pred[[i]])
cat("\n")
cat("Link function:", object$link[[i]])
cat("\n")
cat("Variance function:", object$variance[[i]])
cat("\n")
cat("Covariance function:", object$covariance[[i]])
cat("\n")
cat("Regression:\n")
tab_beta <- coef(object, std.error = TRUE, response = i, type = "beta")[, 1:2]
tab_beta$"Z value" <- tab_beta[, 1]/tab_beta[, 2]
rownames(tab_beta) <- object$beta_names[[i]]
output[i][[1]]$Regression <- tab_beta
print(tab_beta)
cat("\n")
tab_power <- coef(object, std.error = TRUE, response = i, type = "power")[, 1:2]
tab_power$"Z value" <- tab_power[, 1]/tab_power[, 2]
rownames(tab_power) <- NULL
if (dim(tab_power)[1] != 0) {
cat("Power:\n")
print(tab_power)
output[i][[1]]$Power <- tab_power
cat("\n")
}
cat("Dispersion:\n")
tab_tau <- coef(object, std.error = TRUE, response = i, type = "tau")[, 1:2]
tab_tau$"Z value" <- tab_tau[, 1]/tab_tau[, 2]
rownames(tab_tau) <- NULL
output[i][[1]]$tau <- tab_tau
print(tab_tau)
cat("\n")
}
tab_rho <- coef(object, std.error = TRUE, response = NA, type = "correlation")[, c(3, 1, 2)]
tab_rho$"Z value" <- tab_rho[, 2]/tab_rho[, 3]
if (dim(tab_rho)[1] != 0) {
cat("Correlation matrix:\n")
print(tab_rho)
cat("\n")
}
names(object$con$correct) <- ""
iteration_cov <- length(na.exclude(object$IterationCovariance[, 1]))
names(iteration_cov) <- ""
names(object$con$method) <- ""
cat("Algorithm:", object$con$method)
cat("\n")
cat("Correction:", object$con$correct)
cat("\n")
cat("Number iterations:", iteration_cov)
}
#' @title Calculate Variance-Covariance matrix for a fitted McGLM
#' object.
#' @name vcov.mcglm
#'
#' @description Returns the variance-covariance matrix for all
#' parameters of a \code{mcglm} fitted model object.
#'
#' @param object a fitted model \code{mcglm} object.
#' @param ... additional arguments affecting the summary produced. Note
#' that there is no extra options for mcglm object class.
#'
#' @return A variance-covariance matrix.
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#'
#' @method vcov mcglm
#' @export
vcov.mcglm <- function(object, ...) {
cod <- coef(object)$Parameters
colnames(object$vcov) <- cod
rownames(object$vcov) <- cod
return(object$vcov)
}
#' ANOVA method for McGLMs.
#'
#' @description ANOVA method for McGLMS.
#'
#' @param object an object of class mcglm, usually, a result of a call to \code{mcglm}.
#' @param ... additional arguments affecting the summary produced. Note that there is no extra options for
#' mcglm object class.
#' @export
anova.mcglm <- function(object, ...) {
n_resp <- length(object$mu_list)
n_beta <- lapply(object$list_X, ncol)
idx.list <- list()
for (i in 1:n_resp) {
idx.list[[i]] <- rep(i, n_beta[i])
}
vv <- vcov(object)
n_par <- dim(vv)[1]
idx.vec <- do.call(c, idx.list)
n_cov <- n_par - length(idx.vec)
idx.vec <- c(idx.vec, rep(0, n_cov))
temp.vcov <- list()
temp.beta <- list()
for (i in 1:n_resp) {
idx.id = idx.vec == i
temp.vcov[[i]] <- vv[idx.id, idx.id]
temp.beta[[i]] <- coef(object, type = "beta", response = i)$Estimates
}
saida <- list()
for (i in 1:n_resp) {
idx <- attr(object$list_X[[i]], "assign")
names <- colnames(object$list_X[[i]])
if (names[1] == "(Intercept)") {
idx <- idx[-1]
names <- names[-1]
temp.beta[[i]] <- temp.beta[[i]][-1]
temp.vcov[[i]] <- temp.vcov[[i]][-1, -1]
}
n_terms <- length(unique(idx))
X2.resp <- list()
for (j in 1:n_terms) {
idx.TF <- idx == j
temp <- as.numeric(t(temp.beta[[i]][idx.TF]) %*% solve(as.matrix(temp.vcov[[i]])[idx.TF, idx.TF]) %*% temp.beta[[i]][idx.TF])
nbeta.test <- length(temp.beta[[i]][idx.TF])
X2.resp[[j]] <- data.frame(Covariate = names[idx.TF][1], Chi.Square = round(temp, 4), Df = nbeta.test, p.value = round(pchisq(temp,
nbeta.test, lower.tail = FALSE), 4))
}
saida[[i]] <- do.call(rbind, X2.resp)
}
cat("Wald test for fixed effects", "\n")
return(saida)
}
#' Extract model coefficients for mcglm class
#'
#' coef.mcglm is a function which extracts model coefficients from objects of mcglm class.
#' @param object An object of mcglm class.
#' @param std.error Logical. Returns or not the standard errors.
#' @param response A numeric or vector specyfing for which response variables the coefficients
#' should be returned.
#' @param type A string or string vector specyfing which coefficients should be returned.
#' Options are 'beta', 'tau', 'power', 'tau' and 'correlation'.
#' @param ... additional arguments affecting the summary produced. Note that there is no extra options for
#' mcglm object class.
#' @return A data.frame with estimates, parameters names, response number and parameters type.
#' @export
coef.mcglm <- function(object, std.error = FALSE, response = c(NA, 1:length(object$beta_names)),
type = c("beta", "tau", "power", "correlation"), ...) {
n_resp <- length(object$beta_names)
cod_beta <- list()
cod_power <- list()
cod_tau <- list()
type_beta <- list()
type_power <- list()
type_tau <- list()
resp_beta <- list()
resp_power <- list()
resp_tau <- list()
response_for <- 1:n_resp
for (i in response_for) {
cod_beta[[i]] <- paste(paste("beta", i, sep = ""), 0:c(object$Information$n_betas[[i]] - 1), sep = "")
type_beta[[i]] <- rep("beta", length(cod_beta[[i]]))
resp_beta[[i]] <- rep(response_for[i], length(cod_beta[[i]]))
if (object$Information$n_power[[i]] != 0 | object$power_fixed[[i]] == FALSE) {
cod_power[[i]] <- paste(paste("power", i, sep = ""), 1:object$Information$n_power[[i]], sep = "")
type_power[[i]] <- rep("power", length(cod_power[[i]]))
resp_power[[i]] <- rep(response_for[i], length(cod_power[[i]]))
}
if (object$Information$n_power[[i]] == 0) {
cod_power[[i]] <- rep(1, 0)
type_power[[i]] <- rep(1, 0)
resp_power[[i]] <- rep(1, 0)
}
cod_tau[[i]] <- paste(paste("tau", i, sep = ""), 1:object$Information$n_tau[[i]], sep = "")
type_tau[[i]] <- rep("tau", length(cod_tau[[i]]))
resp_tau[[i]] <- rep(response_for[i], length(cod_tau[[i]]))
}
rho_names <- c()
if (n_resp != 1) {
combination <- combn(n_resp, 2)
for (i in 1:dim(combination)[2]) {
rho_names[i] <- paste(paste("rho", combination[1, i], sep = ""), combination[2, i], sep = "")
}
}
type_rho <- rep("correlation", length(rho_names))
resp_rho <- rep(NA, length(rho_names))
cod <- c(do.call(c, cod_beta), rho_names, do.call(c, Map(c, cod_tau)))
type_cod <- c(do.call(c, type_beta), type_rho, do.call(c, Map(c, type_tau)))
response_cod <- c(do.call(c, resp_beta), resp_rho, do.call(c, Map(c, resp_tau)))
if (length(cod_power) != 0) {
cod <- c(do.call(c, cod_beta), rho_names, do.call(c, Map(c, cod_power, cod_tau)))
type_cod <- c(do.call(c, type_beta), type_rho, do.call(c, Map(c, type_power, type_tau)))
response_cod <- c(do.call(c, resp_beta), resp_rho, do.call(c, Map(c, resp_power, resp_tau)))
}
Estimates <- c(object$Regression, object$Covariance)
coef_temp <- data.frame(Estimates = Estimates, Parameters = cod, Type = type_cod, Response = response_cod)
if (std.error == TRUE) {
coef_temp <- data.frame(Estimates = Estimates, Std.error = sqrt(diag(object$vcov)), Parameters = cod, Type = type_cod,
Response = response_cod)
}
output <- coef_temp[which(coef_temp$Response %in% response & coef_temp$Type %in% type), ]
return(output)
}
#' Confidence Intervals for mcglm
#'
#' @description Computes confidence intervals for parameters in a fitted mcglm.
#'
#' @param object a fitted mcglm object.
#' @param parm a specification of which parameters are to be given confidence intervals, either a vector of
#' number or a vector of names. If missing, all parameters are considered.
#' @param level the confidence level required.
#' @param ... additional arguments affecting the confidence interval produced.
#' Note that there is no extra options for
#' mcglm object class.
#' @return A data.frame with confidence intervals, parameters names,
#' response number and parameters type.
#' @export
confint.mcglm <- function(object, parm, level = 0.95, ...) {
temp <- coef(object, std.error = TRUE)
if(missing(parm)){parm <- 1:length(temp$Estimates)}
a <- (1 - level)/2
a <- c(a, 1 - a)
fac <- qnorm(a)
ci <- temp$Estimates + temp$Std.error %o% fac
colnames(ci) <- paste(format(a, 2), "%", sep = "")
rownames(ci) <- temp$Parameters
return(ci[parm])
}
#' Extract Model Fitted Values of McGLM
#'
#' @description Extract fitted values for objects of mcglm class.
#'
#' @param object An object of mcglm class.
#' @param ... additional arguments affecting the summary produced. Note that there is no extra options for
#' mcglm object class.
#' @return Depending on the number of response variable the function \code{fitted.mcglm} returns
#' a vector (univariate models) or a matrix (multivariate models) of fitted values.
#' @export
fitted.mcglm <- function(object, ...) {
n_resp <- length(object$beta_names)
output <- Matrix(object$fitted, ncol = n_resp, nrow = object$n_obs)
return(output)
}
#' Default Multivariate Covariance Generalized Linear models plotting
#'
#' @description takes a fitted mcglm object by \code{mcglm()} and plots, residuals,
#' influence diagnostic measures and algorithm check.
#'
#' @param x a fitted mcglm object as produced by \code{mcglm()}.
#' @param type Specify which graphical analysis will be performed, options are: residuals, influence
#' and algorithm.
#' @param ... additional arguments affecting the plot produced. Note that there is no extra options for
#' mcglm object class.
#' @export
plot.mcglm <- function(x, type = "residuals", ...) {
object = x
n_resp <- length(object$beta_names)
if (type == "residuals") {
par(mar = c(2.6, 2.5, 0.1, 0.1), mgp = c(1.6, 0.6, 0), mfrow = c(2, n_resp))
for (i in 1:n_resp) {
res <- residuals(object, type = "pearson")[, i]
fit_values <- fitted(object)[, i]
plot(res ~ fit_values, ylab = "Pearson residuals", xlab = "Fitted values")
temp <- loess.smooth(fitted(object)[, i], residuals(object, type = "pearson")[, i])
lines(temp$x, temp$y)
qqnorm(res)
qqline(res)
}
}
if (type == "algorithm") {
n_iter <- length(na.exclude(object$IterationCovariance[, 1]))
par(mar = c(2.6, 2.5, 0.1, 0.1), mgp = c(1.6, 0.6, 0), mfrow = c(2, 2))
matplot(object$IterationRegression[1:c(n_iter + 5), ], type = "l", lty = 2, ylab = "Regression", xlab = "Iterations")
matplot(object$IterationCovariance[1:c(n_iter + 5), ], type = "l", lty = 2, ylab = "Covariance", xlab = "Iterations")
matplot(object$ScoreRegression[1:c(n_iter + 5), ], type = "l", lty = 2, ylab = "Quasi-score Regression", xlab = "Iterations")
matplot(object$ScoreCovariance[1:c(n_iter + 5), ], type = "l", lty = 2, ylab = "Quasi-score Covariance", xlab = "Iterations")
}
if (type == "partial_residuals") {
list_beta <- mc_updateBeta(list_initial = object$list_initial, betas = object$Regression, n_resp = n_resp, information = object$Information)
comp_X <- list()
for (i in 1:n_resp) {
comp_X[[i]] <- as.matrix(object$list_X[[i]]) * as.numeric(list_beta$regression[[i]])
}
for (i in 1:n_resp) {
res <- residuals(object, type = "pearson")[, i]
dev.new()
n_cov <- dim(comp_X[[i]])[2]
par(mar = c(2.6, 2.5, 0.5, 0.5), mgp = c(1.6, 0.6, 0), mfrow = c(1, c(n_cov - 1)))
for (j in 2:n_cov) {
p1 <- comp_X[[i]][, j] + res
plot(p1 ~ object$list_X[[i]][, j], xlab = object$beta_names[[i]][j], ylab = "Partial residuals ")
}
}
}
}
#' Print method for Multivariate Covariance Generalized Linear Model
#'
#' @description The default print method for a mcglm object.
#' @method print mcglm
#' @param x fitted model objects of class mcglm as produced by mcglm().
#' @param ... further arguments passed to or from other methods.
#' @export
#' @rdname print.mcglm
print.mcglm <- function(x, ...) {
object <- x
n_resp <- length(object$beta_names)
regression <- mc_updateBeta(list_initial = list(), betas = object$Regression, information = object$Information, n_resp = n_resp)
for (i in 1:n_resp) {
cat("Call: ")
print(object$linear_pred[[i]])
cat("\n")
cat("Link function:", object$link[[i]])
cat("\n")
cat("Variance function:", object$variance[[i]])
cat("\n")
cat("Covariance function:", object$covariance[[i]])
cat("\n")
names(regression[[1]][[i]]) <- object$beta_names[[i]]
cat("Regression:\n")
print(regression[[1]][[i]])
cat("\n")
cat("Dispersion:\n")
tau_temp <- coef(object, response = i, type = "tau")$Estimate
names(tau_temp) <- rep("", length(tau_temp))
print(tau_temp)
cat("\n")
power_temp <- coef(object, response = i, type = "power")$Estimate
if (length(power_temp) != 0) {
names(power_temp) <- ""
cat("Power:\n")
print(power_temp)
cat("\n")
}
}
}
#' Compute Quasi Information Criterion (QIC) for McGLMs.
#' @title Compute Quasi Information Criterion (QIC) for McGLMs.
#' @name qic.mcglm
#'
#' @description \code{qic.mcglm} is a function which computes the QIC
#' for McGLMs.
#'
#' @param object An object of \code{mcglm} class.
#' @param object.iid An object of \code{mcglm} class contained the model
#' fitted using independent covariance structure.
#'
#' qic.mcglm is a function which computes the QIC for McGLMs.
#' @param object An object of mcglm class.
#' @param object.iid An object of mcglm class contained the model fitted using independent
#' covariance structure.
#' @return The QIC value.
#'
#' @author Wagner Hugo Bonat, \email{wbonat@@ufpr.br}
#'
#' @export
qic.mcglm <- function(object, object.iid) {
......@@ -19,7 +26,10 @@ qic.mcglm <- function(object, object.iid) {
if (object$power_fixed[[i]] == TRUE) {
power = object$list_initial$power[[i]]
}
Q[, i] <- mc_qll(y = obs[, i], mu = mu[, i], variance = object$variance[[i]], power = power)
Q[, i] <- mc_qll(y = obs[, i],
mu = mu[, i],
variance = object$variance[[i]],
power = power)
}
Vbeta <- -object$inv_S_beta
Vnull <- solve(-object.iid$inv_S_beta)
......
#' Residuals for Multivariate Covariance Generalized Linear Models (McGLM)
#'
#' @description Compute residuals based on fitting mcglm models.
#'
#' @param object An of class mcglm, typically the result of a call to \code{mcglm}.
#' @param type the type of residuals which should be returned. The alternatives are: 'raw'
#' (default), 'pearson' and 'standardized'.
#' @param ... additional arguments affecting the residuals produced. Note that there is no extra options for
#' mcglm object class.
#' @return Depending on the number of response variable the function \code{residuals.mcglm} returns
#' a vector (univariate models) or a matrix (multivariate models) of residuals values.
#' @export
residuals.mcglm <- function(object, type = "raw", ...) {
n_resp <- length(object$beta_names)
output <- Matrix(object$residuals, ncol = n_resp, nrow = object$n_obs)
if (type == "standardized") {
output <- Matrix(as.numeric(object$residuals %*% chol(object$inv_C)), ncol = n_resp, nrow = object$n_obs)
}
if (type == "pearson") {
output <- Matrix(as.numeric(object$residuals/sqrt(diag(object$C))), ncol = n_resp, nrow = object$n_obs)
}
return(output)
}
#' Summarizing Multivariate Covariance Generalized Linear Models fits.
#'
#' @description Summary for McGLMs objects.
#'
#' @param object an object of class mcglm, usually, a result of a call to \code{mcglm}.
#' @param ... additional arguments affecting the summary produced. Note the there is no extra options for
#' mcglm object class.
#' @return Print an mcglm object.
#' @method summary mcglm
#' @export
summary.mcglm <- function(object, ...) {
n_resp <- length(object$beta_names)
output <- list()
for (i in 1:n_resp) {
cat("Call: ")
print(object$linear_pred[[i]])
cat("\n")
cat("Link function:", object$link[[i]])
cat("\n")
cat("Variance function:", object$variance[[i]])
cat("\n")
cat("Covariance function:", object$covariance[[i]])
cat("\n")
cat("Regression:\n")
tab_beta <- coef(object, std.error = TRUE, response = i, type = "beta")[, 1:2]
tab_beta$"Z value" <- tab_beta[, 1]/tab_beta[, 2]
rownames(tab_beta) <- object$beta_names[[i]]
output[i][[1]]$Regression <- tab_beta
print(tab_beta)
cat("\n")
tab_power <- coef(object, std.error = TRUE, response = i, type = "power")[, 1:2]
tab_power$"Z value" <- tab_power[, 1]/tab_power[, 2]
rownames(tab_power) <- NULL
if (dim(tab_power)[1] != 0) {
cat("Power:\n")
print(tab_power)
output[i][[1]]$Power <- tab_power
cat("\n")
}
cat("Dispersion:\n")
tab_tau <- coef(object, std.error = TRUE, response = i, type = "tau")[, 1:2]
tab_tau$"Z value" <- tab_tau[, 1]/tab_tau[, 2]
rownames(tab_tau) <- NULL
output[i][[1]]$tau <- tab_tau
print(tab_tau)
cat("\n")
}
tab_rho <- coef(object, std.error = TRUE, response = NA, type = "correlation")[, c(3, 1, 2)]
tab_rho$"Z value" <- tab_rho[, 2]/tab_rho[, 3]
if (dim(tab_rho)[1] != 0) {
cat("Correlation matrix:\n")
print(tab_rho)
cat("\n")
}
names(object$con$correct) <- ""
iteration_cov <- length(na.exclude(object$IterationCovariance[, 1]))
names(iteration_cov) <- ""
names(object$con$method) <- ""
cat("Algorithm:", object$con$method)
cat("\n")
cat("Correction:", object$con$correct)
cat("\n")
cat("Number iterations:", iteration_cov)
}
#' Calculate Variance-Covariance matrix for a fitted McGLM object.
#'
#' @description Returns the variance-covariance matrix for all parameters of a mcglm fitted model object.
#'
#' @param object a fitted model mcglm object.
#' @param ... additional arguments affecting the summary produced. Note that there is no extra options for
#' mcglm object class.
#' @return A variance-covariance matrix.
#' @export
vcov.mcglm <- function(object, ...) {
cod <- coef(object)$Parameters
colnames(object$vcov) <- cod
rownames(object$vcov) <- cod
return(object$vcov)
}
##' @title Multivariate covariance generalized linear models (McGLMs)
##'
##'
##' @description Fits a multivariate covariance generalized linear models (McGLMs) to data.
##' McGLM is a general framework for non-normal multivariate data analysis, designed to handle
##' multivariate response variables, along with a wide range of temporal and spatial correlation
##' structures defined in terms of a covariance link function combined with a matrix linear predictor
##' involving known matrices. The models take non-normality into account in the conventional way by means
##' of a variance function, and the mean structure is modelled by means of a link function and a linear predictor.
##' The models are fitted using an efficient Newton scoring algorithm based on quasi-likelihood and
##' Pearson estimating functions, using only second-moment assumptions.
##' This provides a unified approach to a wide variety of different types of response variables and
##' covariance structures, including multivariate extensions of repeated measures, time series, longitudinal,
##' spatial and spatio-temporal structures.
##'
##' @description Fits a multivariate covariance generalized linear
##' models (McGLMs) to data. McGLM is a general framework for
##' non-normal multivariate data analysis, designed to handle
##' multivariate response variables, along with a wide range of
##' temporal and spatial correlation structures defined in terms of
##' a covariance link function combined with a matrix linear
##' predictor involving known matrices. The models take
##' non-normality into account in the conventional way by means of a
##' variance function, and the mean structure is modelled by means
##' of a link function and a linear predictor. The models are
##' fitted using an efficient Newton scoring algorithm based on
##' quasi-likelihood and Pearson estimating functions, using only
##' second-moment assumptions. This provides a unified approach to
##' a wide variety of different types of response variables and
##' covariance structures, including multivariate extensions of
##' repeated measures, time series, longitudinal, spatial and
##' spatio-temporal structures.
##'
##' @docType package
##' @name mcglm
NULL
#' @name ahs
#'
#' @title Australian health survey
#' @name ahs
#'
#' @description The Australian health survey was used by Bonat and Jorgensen (2015) as an example of multivariate
#' count regression model. The data consists of five count response variables concerning health system access
#' measures and nine covariates concerning social conditions in Australian for 1987-88.
#' @description The Australian health survey was used by Bonat and
#' Jorgensen (2015) as an example of multivariate count regression
#' model. The data consists of five count response variables
#' concerning health system access measures and nine covariates
#' concerning social conditions in Australian for 1987-88.
#'
#' \itemize{
#'
#' \item \code{sex} - Factor, two levels (0-Male; 1-Female).
#'
#' \item \code{age} - Respondent's age in years divided by 100.
#' \item \code{income} - Respondent's annual income in Australian dollars divided by 1000.
#' \item \code{levyplus} - Factor, two levels (1- if respondent is covered by private health
#' insurance fund for private patients in public hospital (with doctor of choice); 0 - otherwise).
#' \item \code{freepoor} - Factor, two levels (1 - if respondent is covered by government because low income,
#' recent immigrant, unemployed; 0 - otherwise).
#' \item \code{freerepa} - Factor, two levels (1 - if respondent is covered free by government because of
#' old-age or disability pension, or because invalid veteran or family of deceased veteran; 0 - otherwise).
#' \item \code{illnes} - Number of illnesses in past 2 weeks, with 5 or more weeks coded as 5.
#' \item \code{actdays} - Number of days of reduced activity in the past two weeks due to illness or injury.
#' \item \code{hscore} - Respondent's general health questionnaire score using Goldberg's method;
#' high score indicates poor health.
#' \item \code{chcond1} - Factor, two levels (1 - if respondent has chronic condition(s) but is not limited
#' in activity; 0 - otherwise).
#' \item \code{chcond2} - Factor, two levels (1 if respondent has chronic condition(s) and is limited in
#' activity; 0 - otherwise).
#' \item \code{Ndoc} - Number of consultations with a doctor or specialist (response variable).
#' \item \code{Nndoc} - Number of consultations with health professionals (response variable).
#' \item \code{Nadm} - Number of admissions to a hospital, psychiatric hospital, nursing or
#' convalescence home in the past 12 months (response variable).
#' \item \code{Nhosp} - Number of nights in a hospital during the most recent admission.
#' \item \code{Nmed} - Total number of prescribed and non prescribed medications used in the past two days.
#'
#' \item \code{income} - Respondent's annual income in Australian
#' dollars divided by 1000.
#'
#' \item \code{levyplus} - Factor, two levels (1- if respondent is
#' covered by private health insurance fund for private patients in
#' public hospital (with doctor of choice); 0 - otherwise).
#'
#' \item \code{freepoor} - Factor, two levels (1 - if respondent is
#' covered by government because low income, recent immigrant,
#' unemployed; 0 - otherwise).
#'
#' \item \code{freerepa} - Factor, two levels (1 - if respondent is
#' covered free by government because of old-age or disability
#' pension, or because invalid veteran or family of deceased
#' veteran; 0 - otherwise).
#'
#' \item \code{illnes} - Number of illnesses in past 2 weeks, with 5 or
#' more weeks coded as 5.
#'
#' \item \code{actdays} - Number of days of reduced activity in the past
#' two weeks due to illness or injury.
#'
#' \item \code{hscore} - Respondent's general health questionnaire score
#' using Goldberg's method; high score indicates poor health.
#'
#' \item \code{chcond1} - Factor, two levels (1 - if respondent has
#' chronic condition(s) but is not limited in activity; 0 -
#' otherwise).
#'
#' \item \code{chcond2} - Factor, two levels (1 if respondent has
#' chronic condition(s) and is limited in activity; 0 - otherwise).
#'
#' \item \code{Ndoc} - Number of consultations with a doctor or
#' specialist (response variable).
#'
#' \item \code{Nndoc} - Number of consultations with health
#' professionals (response variable).
#'
#' \item \code{Nadm} - Number of admissions to a hospital, psychiatric
#' hospital, nursing or convalescence home in the past 12 months
#' (response variable).
#'
#' \item \code{Nhosp} - Number of nights in a hospital during the most
#' recent admission.
#'
#' \item \code{Nmed} - Total number of prescribed and non prescribed
#' medications used in the past two days.
#'
#' \item \code{id} - Respondent's index.
#'
#' }
#'
#' @docType data
......@@ -61,7 +97,8 @@ NULL
#'
#' @format a \code{data.frame} with 5190 records and 17 variables.
#'
#' @source Deb, P. and Trivedi, P. K. (1997). Demand for medical care by the elderly: A finite mixture approach,
#' Journal of Applied Econometrics 12(3):313--336.
#' @source Deb, P. and Trivedi, P. K. (1997). Demand for medical care by
#' the elderly: A finite mixture approach, Journal of Applied
#' Econometrics 12(3):313--336.
#'
NULL
......@@ -6,42 +6,76 @@
\title{Australian health survey}
\format{a \code{data.frame} with 5190 records and 17 variables.}
\source{
Deb, P. and Trivedi, P. K. (1997). Demand for medical care by the elderly: A finite mixture approach,
Journal of Applied Econometrics 12(3):313--336.
Deb, P. and Trivedi, P. K. (1997). Demand for medical care by
the elderly: A finite mixture approach, Journal of Applied
Econometrics 12(3):313--336.
}
\usage{
data(ahs)
}
\description{
The Australian health survey was used by Bonat and Jorgensen (2015) as an example of multivariate
count regression model. The data consists of five count response variables concerning health system access
measures and nine covariates concerning social conditions in Australian for 1987-88.
The Australian health survey was used by Bonat and
Jorgensen (2015) as an example of multivariate count regression
model. The data consists of five count response variables
concerning health system access measures and nine covariates
concerning social conditions in Australian for 1987-88.
\itemize{
\item \code{sex} - Factor, two levels (0-Male; 1-Female).
\item \code{age} - Respondent's age in years divided by 100.
\item \code{income} - Respondent's annual income in Australian dollars divided by 1000.
\item \code{levyplus} - Factor, two levels (1- if respondent is covered by private health
insurance fund for private patients in public hospital (with doctor of choice); 0 - otherwise).
\item \code{freepoor} - Factor, two levels (1 - if respondent is covered by government because low income,
recent immigrant, unemployed; 0 - otherwise).
\item \code{freerepa} - Factor, two levels (1 - if respondent is covered free by government because of
old-age or disability pension, or because invalid veteran or family of deceased veteran; 0 - otherwise).
\item \code{illnes} - Number of illnesses in past 2 weeks, with 5 or more weeks coded as 5.
\item \code{actdays} - Number of days of reduced activity in the past two weeks due to illness or injury.
\item \code{hscore} - Respondent's general health questionnaire score using Goldberg's method;
high score indicates poor health.
\item \code{chcond1} - Factor, two levels (1 - if respondent has chronic condition(s) but is not limited
in activity; 0 - otherwise).
\item \code{chcond2} - Factor, two levels (1 if respondent has chronic condition(s) and is limited in
activity; 0 - otherwise).
\item \code{Ndoc} - Number of consultations with a doctor or specialist (response variable).
\item \code{Nndoc} - Number of consultations with health professionals (response variable).
\item \code{Nadm} - Number of admissions to a hospital, psychiatric hospital, nursing or
convalescence home in the past 12 months (response variable).
\item \code{Nhosp} - Number of nights in a hospital during the most recent admission.
\item \code{Nmed} - Total number of prescribed and non prescribed medications used in the past two days.
\item \code{income} - Respondent's annual income in Australian
dollars divided by 1000.
\item \code{levyplus} - Factor, two levels (1- if respondent is
covered by private health insurance fund for private patients in
public hospital (with doctor of choice); 0 - otherwise).
\item \code{freepoor} - Factor, two levels (1 - if respondent is
covered by government because low income, recent immigrant,
unemployed; 0 - otherwise).
\item \code{freerepa} - Factor, two levels (1 - if respondent is
covered free by government because of old-age or disability
pension, or because invalid veteran or family of deceased
veteran; 0 - otherwise).
\item \code{illnes} - Number of illnesses in past 2 weeks, with 5 or
more weeks coded as 5.
\item \code{actdays} - Number of days of reduced activity in the past
two weeks due to illness or injury.
\item \code{hscore} - Respondent's general health questionnaire score
using Goldberg's method; high score indicates poor health.
\item \code{chcond1} - Factor, two levels (1 - if respondent has
chronic condition(s) but is not limited in activity; 0 -
otherwise).
\item \code{chcond2} - Factor, two levels (1 if respondent has
chronic condition(s) and is limited in activity; 0 - otherwise).
\item \code{Ndoc} - Number of consultations with a doctor or
specialist (response variable).
\item \code{Nndoc} - Number of consultations with health
professionals (response variable).
\item \code{Nadm} - Number of admissions to a hospital, psychiatric
hospital, nursing or convalescence home in the past 12 months
(response variable).
\item \code{Nhosp} - Number of nights in a hospital during the most
recent admission.
\item \code{Nmed} - Total number of prescribed and non prescribed
medications used in the past two days.
\item \code{id} - Respondent's index.
}
}
\keyword{datasets}
......
% Generated by roxygen2 (4.1.1): do not edit by hand
% Please edit documentation in R/mc_anova.R
% Please edit documentation in R/mc_S3_methods.R
\name{anova.mcglm}
\alias{anova.mcglm}
\title{ANOVA method for McGLMs.}
......@@ -7,12 +7,22 @@
\method{anova}{mcglm}(object, ...)
}
\arguments{
\item{object}{an object of class mcglm, usually, a result of a call to \code{mcglm}.}
\item{object}{an object of class \code{mcglm}, usually, a result of a
call to \code{mcglm()}.}
\item{...}{additional arguments affecting the summary produced. Note that there is no extra options for
mcglm object class.}
\item{...}{additional arguments affecting the summary produced. Note
that there is no extra options for mcglm object class.}
}
\value{
A \code{data.frame} with Chi-square statistic to test the
null hypothesis of a parameter, or a set of parameters, be
zero. The Wald test based on the observed covariance matrix of
the parameters is used.
}
\description{
ANOVA method for McGLMS.
ANOVA method for object of class McGLMS.
}
\author{
Wagner Hugo Bonat, \email{wbonat@ufpr.br}
}
% Generated by roxygen2 (4.1.1): do not edit by hand
% Please edit documentation in R/mc_coef.R
% Please edit documentation in R/mc_S3_methods.R
\name{coef.mcglm}
\alias{coef.mcglm}
\title{Extract model coefficients for mcglm class}
......@@ -9,23 +9,31 @@
"correlation"), ...)
}
\arguments{
\item{object}{An object of mcglm class.}
\item{object}{An object of \code{mcglm} class.}
\item{std.error}{Logical. Returns or not the standard errors.}
\item{std.error}{Logical. If \code{TRUE} returns the standard errors
of the estimates. Default is \code{FALSE}.}
\item{response}{A numeric or vector specyfing for which response variables the coefficients
should be returned.}
\item{response}{A numeric vector specyfing for which response
variables the coefficients should be returned.}
\item{type}{A string or string vector specyfing which coefficients should be returned.
Options are 'beta', 'tau', 'power', 'tau' and 'correlation'.}
\item{type}{A string vector (can be 1 element length) specyfing which
coefficients should be returned. Options are \code{"beta"},
\code{"tau"}, \code{"power"}, \code{"tau"} and
\code{"correlation"}.}
\item{...}{additional arguments affecting the summary produced. Note that there is no extra options for
mcglm object class.}
\item{...}{additional arguments affecting the summary produced. Note
that there is no extra options for mcglm object class.}
}
\value{
A data.frame with estimates, parameters names, response number and parameters type.
A \code{data.frame} with parameters names, estimates,
response number and parameters type.
}
\description{
coef.mcglm is a function which extracts model coefficients from objects of mcglm class.
\code{coef.mcglm} is a function which extracts model
coefficients from objects of \code{mcglm} class.
}
\author{
Wagner Hugo Bonat, \email{wbonat@ufpr.br}
}
% Generated by roxygen2 (4.1.1): do not edit by hand
% Please edit documentation in R/mc_confint.mcglm.R
% Please edit documentation in R/mc_S3_methods.R
\name{confint.mcglm}
\alias{confint.mcglm}
\title{Confidence Intervals for mcglm}
......@@ -7,22 +7,27 @@
\method{confint}{mcglm}(object, parm, level = 0.95, ...)
}
\arguments{
\item{object}{a fitted mcglm object.}
\item{object}{a fitted \code{mcglm} object.}
\item{parm}{a specification of which parameters are to be given confidence intervals, either a vector of
number or a vector of names. If missing, all parameters are considered.}
\item{parm}{a specification of which parameters are to be given
confidence intervals, either a vector of number or a vector of
strings. If missing, all parameters are considered.}
\item{level}{the confidence level required.}
\item{level}{the nominal confidence level.}
\item{...}{additional arguments affecting the confidence interval produced.
Note that there is no extra options for
mcglm object class.}
\item{...}{additional arguments affecting the confidence interval
produced. Note that there is no extra options for \code{mcglm}
object class.}
}
\value{
A data.frame with confidence intervals, parameters names,
response number and parameters type.
A \code{data.frame} with confidence intervals, parameters
names, response number and parameters type.
}
\description{
Computes confidence intervals for parameters in a fitted mcglm.
Computes confidence intervals for parameters in a fitted
\code{mcglm} model.
}
\author{
Wagner Hugo Bonat, \email{wbonat@ufpr.br}
}
% Generated by roxygen2 (4.1.1): do not edit by hand
% Please edit documentation in R/mc_fitted.mcglm.R
% Please edit documentation in R/mc_S3_methods.R
\name{fitted.mcglm}
\alias{fitted.mcglm}
\title{Extract Model Fitted Values of McGLM}
......@@ -7,16 +7,20 @@
\method{fitted}{mcglm}(object, ...)
}
\arguments{
\item{object}{An object of mcglm class.}
\item{object}{An object of \code{mcglm} class.}
\item{...}{additional arguments affecting the summary produced. Note that there is no extra options for
mcglm object class.}
\item{...}{additional arguments affecting the summary produced. Note
that there is no extra options for \code{mcglm} object class.}
}
\value{
Depending on the number of response variable the function \code{fitted.mcglm} returns
a vector (univariate models) or a matrix (multivariate models) of fitted values.
Depending on the number of response variable, the function
\code{fitted.mcglm} returns a vector (univariate models) or a
matrix (multivariate models) of fitted values.
}
\description{
Extract fitted values for objects of mcglm class.
Extract fitted values for objects of \code{mcglm} class.
}
\author{
Wagner Hugo Bonat, \email{wbonat@ufpr.br}
}
......@@ -52,16 +52,22 @@ fitted using an estimating function approach, combining quasi-score functions fo
parameters and Pearson estimating function for covariance parameters. For details see Bonat and
Jorgensen (2015).
Fits a multivariate covariance generalized linear models (McGLMs) to data.
McGLM is a general framework for non-normal multivariate data analysis, designed to handle
multivariate response variables, along with a wide range of temporal and spatial correlation
structures defined in terms of a covariance link function combined with a matrix linear predictor
involving known matrices. The models take non-normality into account in the conventional way by means
of a variance function, and the mean structure is modelled by means of a link function and a linear predictor.
The models are fitted using an efficient Newton scoring algorithm based on quasi-likelihood and
Pearson estimating functions, using only second-moment assumptions.
This provides a unified approach to a wide variety of different types of response variables and
covariance structures, including multivariate extensions of repeated measures, time series, longitudinal,
spatial and spatio-temporal structures.
Fits a multivariate covariance generalized linear
models (McGLMs) to data. McGLM is a general framework for
non-normal multivariate data analysis, designed to handle
multivariate response variables, along with a wide range of
temporal and spatial correlation structures defined in terms of
a covariance link function combined with a matrix linear
predictor involving known matrices. The models take
non-normality into account in the conventional way by means of a
variance function, and the mean structure is modelled by means
of a link function and a linear predictor. The models are
fitted using an efficient Newton scoring algorithm based on
quasi-likelihood and Pearson estimating functions, using only
second-moment assumptions. This provides a unified approach to
a wide variety of different types of response variables and
covariance structures, including multivariate extensions of
repeated measures, time series, longitudinal, spatial and
spatio-temporal structures.
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment