R/matrixMethod.R

In stefvanbuuren/dsMice: DataSHIELD server-side functions for the mice package

#'
#' @title The Method of Matrices
#' @description The method of matrices consists in build the matrix \emph{A} and the matrix \emph{g}.
#' @importFrom Rdpack reprompt
#' @details If the experimenter has \emph{n} independent variables and \emph{m} observations,
#' each of which can be expressed by the equation:
#' \deqn{y[i] = \beta[0] + \beta[1]x[i] + \beta[2]x[2i] + ... + \beta[n]x[ni] + \epsilon[i]}
#' This model represents the equations system that describe the response values generated by fitting the regression.
#' Using the matrix notation, according to \insertCite{walpole1993probability}{dsMice},
#' the equations system can be written as
#' \deqn{y = X\beta + \epsilon}
#' In the \emph{X} matrix, apart from the initial element, the ith row represents the x-values that give
#' rise to the response \emph{y[i]}.
#' So \eqn{A = X'X}, \eqn{g = X'y} and the normal equations can be put in the matrix form \eqn{AB = g}.
#'
#' Requires the function \code{\link{getVarbyFormula}}.
#' @param formula a string character to be transformed as an object of class \code{formula}.
#' @param subset optional vector which specifies the subset of observations to be used in the analysis.
#' @param weight a list or numerical value for weigthed analysis.
#' @param family a string character with the name of the error distribution and link function to be used in the analysis.
#' @return A list with components:
#' \item{xtx}{a data matrix, the \emph{A} matrix.}
#' \item{xty}{a data matrix, the \emph{g} matrix.}
#' \item{sum.y}{a numerical value, the sum of elements for a given dependent variable \emph{y}.}
#' \item{n.rows}{a numerical value, the sample size.}
#'
#' @author Paula R. Costa e Silva
#' @section Dependencies:
#' \code{\link{getVarbyFormula}}
#' @references
#' \insertRef{walpole1993probability}{dsMice}
#'
#' @export
#'

matrixMethod <- function(formula, subset=NULL, weight=1, family=NULL) {
  bindxy <- getVarbyFormula(formula, subset, weight, family)
  
  bind.x <- data.matrix(bindxy$x)
  bind.y <- data.matrix(bindxy$y)
  
  xtx <- t(bind.x) %*% bind.x
  xty <- t(bind.x) %*% bind.y
  
  sum.y <- sum(bind.y)
  n.rows <- nrow(bind.y)
  
  return(list(xtx = xtx, xty = xty, sum.y = sum.y, n.rows = n.rows))
}