SOP (version 1.0-1)

sop.fit: Fitting generalised linear mixed models with overlapping precision matrices.

Description

This is an internal function of package SOP. It is used to fit SOP models by specifying the design matrices for the fixed and random effects as well as the precision matrices for each variance component in the model.

Usage

sop.fit(y, X, Z, weights = NULL, G = NULL, vcstart = NULL, 
  etastart = NULL, mustart = NULL, offset = NULL, 
  family = gaussian(), control = sop.control())

Value

A list containing the following objects:

b.fixed

the estimated fixed effect coefficients.

b.random

the predicted random effect coefficients.

residuals

the (deviance) residuals.

fitted.values

the fitted values.

linear.predictor

the values of the linear predictor.

X

the fixed effect design matrix.

Z

the random effect design matrix.

y

the response.

weights

the prior weights.

family

the distribution family.

out

a list with i) tol.ol tolerance parameter (outer loop); ii) it.ol number of iteration (outer loop); iii) tol.il tolerance parameter (inner loop); it.il (number of iteration (inner loop)), iv) vc variance components estimates, v) edf effective degrees of freedom.

deviance

the deviance.

null.deviance

the null deviance.

Vp

Bayesian posterior covariance matrix for the coefficients.

Arguments

y

vector of observations of length n.

X

design matrix for the fixed effects (dimension n x p).

Z

design matrix for the random effects (of dimension n x q).

weights

an optional vector of 'prior weights' to be used in the fitting process. If NULL (default), the weights are considered to be one.

G

a list with the diagonal elements of the precision matrices for each variance component in the model. Each element of the list is a vector of the same length as the number of columns in Z (i.e. q). The vector can be padded out with zeroes to indicate random coefficient not `affected' by the variance component (see details).

vcstart

optional numeric vector. Initial values for the variance components (including the error variance as the first element of the vector). If NULL, all variance components are initialised to one.

etastart

initial values for the linear predictor.

mustart

initial values for the expected response.

offset

this can be used to specify an a priori known component to be included in the linear predictor during fitting. This should be NULL or a numeric vector of length equal to the number of observations.

family

object of class family specifying the distribution and link function.

control

a list of control values to replace the default values returned by the function sop.control.

Details

sop.fit is the workhorse function: it is typically not normally called directly but can be more efficient where the response vector 'y', design matrixs 'X' and 'Z', and precision matrices 'G' have already been calculated. Currently, the funcion only allows for diagonal precision matrices (possibly overlappping).

References

Rodriguez-Alvarez, M.X., Lee, D. J., Kneib, T., Durban, M., and Eilers, P. (2015). Fast smoothing parameter separation in multidimensional generalized P-splines: the SAP algorithm. Statistics and Computing, 25 (5), 941--957.

Rodriguez-Alvarez, M.X., Durban, M., Lee, D. J. and Eilers, P. (2019). On the estimation of variance parameters in non-standard generalised linear mixed models: application to penalised smoothing. Statistics and Computing, 29 (3), 483--500.

Examples

library(SOP)
# Simulate the data
set.seed(123)
n <- 1000
sigma <- 0.1
x1 <- runif(n, -1, 1)
x2 <- runif(n, -1, 1)
f0 <- function(x1, x2) cos(2*pi*sqrt((x1 - 0.5)^2 + (x2 - 0.5)^2))
f <- f0(x1, x2)
y <- f + rnorm(n, 0, sigma)
dat <- data.frame(x1 = x1, x2 = x2, y = y)

# Save but not fit the model
m0_nfit <- sop(formula = y ~ f(x1, x2, nseg = 10), data = dat, 
	    	fit = FALSE)

# Now fit using sop.fit()
m0 <- sop.fit(X = m0_nfit$X, Z = m0_nfit$Z, G = m0_nfit$G, 
	y = m0_nfit$y, weights = m0_nfit$weights, 
	   control = list(trace = FALSE))

names(m0)