remlscore

0th

Percentile

REML for Heteroscedastic Regression

Fits a heteroscedastic regression model using residual maximum likelihood (REML).

Keywords
regression
Usage
remlscore(y, X, Z, trace=FALSE, tol=1e-5, maxit=40)
Arguments
y
numeric vector of responses
X
design matrix for predicting the mean
Z
design matrix for predicting the variance
trace
Logical variable. If true then output diagnostic information at each iteration.
tol
Convergence tolerance
maxit
Maximum number of iterations allowed
Details

Write $\mu_i=E(y_i)$ for the expectation of the $i$th response and $s_i=\var(y_i)$. We assume the heteroscedastic regression model $$\mu_i=\bold{x}_i^T\bold{\beta}$$ $$\log(\sigma^2_i)=\bold{z}_i^T\bold{\gamma},$$ where $x_i$ and $z_i$ are vectors of covariates, and $\beta$ and $\gamma$ are vectors of regression coefficients affecting the mean and variance respectively. Parameters are estimated by maximizing the REML likelihood using REML scoring as described in Smyth (2002).

Value

  • List with the following components:
  • betaVector of regression coefficients for predicting the mean
  • se.beta
  • gammaVector of regression coefficients for predicting the variance
  • se.gamStandard errors for gamma
  • muEstimated means
  • phiEstimated variances
  • devianceMinus twice the REML log-likelihood
  • hLeverages

References

Smyth, G. K. (2002). An efficient algorithm for REML in heteroscedastic regression. Journal of Computational and Graphical Statistics 11, 1-12.

Aliases
  • remlscore
Examples
data(welding)
attach(welding)
y <- Strength
# Reproduce results from Table 1 of Smyth (2002)
X <- cbind(1,(Drying+1)/2,(Material+1)/2)
colnames(X) <- c("1","B","C")
Z <- cbind(1,(Material+1)/2,(Method+1)/2,(Preheating+1)/2)
colnames(Z) <- c("1","C","H","I")
out <- remlscore(y,X,Z)
cbind(Estimate=out$gamma,SE=out$se.gam)
Documentation reproduced from package statmod, version 1.2.0, License: LGPL version 2 or newer

Community examples

Looks like there are no examples yet.