Main Page | File List | Globals

ssclme.c File Reference

#include "ssclme.h"

Defines

#define slot_dup(dest, src, sym)   SET_SLOT(dest, sym, duplicate(GET_SLOT(src, sym)))

Functions

void ssclme_copy_ctab (int nf, const int nc[], SEXP ctab, SEXP ssc)
void ssclme_calc_maxod (int n, int Parent[])
SEXP ssclme_create (SEXP facs, SEXP ncv)
void bVj_to_A (int ncj, int Gpj, int Gpjp, const double bVj[], const int Ap[], const int Ai[], double Ax[])
SEXP ssclme_transfer_dimnames (SEXP x, SEXP facs, SEXP mmats)
SEXP ssclme_update_mm (SEXP x, SEXP facs, SEXP mmats)
SEXP ssclme_inflate_and_factor (SEXP x)
SEXP ssclme_factor (SEXP x)
int ldl_update_ind (int probe, int start, const int ind[])
SEXP ldl_inverse (SEXP x)
SEXP ssclme_invert (SEXP x)
SEXP ssclme_initial (SEXP x)
SEXP ssclme_fixef (SEXP x)
SEXP ssclme_ranef (SEXP x)
SEXP ssclme_sigma (SEXP x, SEXP REML)
int coef_length (int nf, const int nc[])
SEXP ssclme_coef (SEXP x)
SEXP ssclme_coefUnc (SEXP x)
SEXP ssclme_coefGetsUnc (SEXP x, SEXP coef)
SEXP ssclme_coefGets (SEXP x, SEXP coef)
SEXP ssclme_EMsteps (SEXP x, SEXP nsteps, SEXP REMLp, SEXP verb)
SEXP ssclme_gradient (SEXP x, SEXP REMLp, SEXP Uncp)
SEXP ssclme_fitted (SEXP x, SEXP facs, SEXP mmats, SEXP useRf)
SEXP ssclme_variances (SEXP x)
SEXP ssclme_collapse (SEXP x)
SEXP ssclme_to_lme (SEXP call, SEXP facs, SEXP x, SEXP model, SEXP REML, SEXP rep, SEXP fitted, SEXP residuals)

Define Documentation

#define slot_dup dest,
src,
sym   )     SET_SLOT(dest, sym, duplicate(GET_SLOT(src, sym)))
 


Function Documentation

void bVj_to_A int  ncj,
int  Gpj,
int  Gpjp,
const double  bVj[],
const int  Ap[],
const int  Ai[],
double  Ax[]
[static]
 

Copy information on Z'Z accumulated in the bVar array to Z'Z

Parameters:
ncj number of columns in this block
Gpj initial column for this group
Gpjp initial column for the next group
bVj pointer to the ncj x ncj x mj array to be filled
Ap column pointer array for Z'Z
Ai row indices for Z'Z
Ax elements of Z'Z

int coef_length int  nf,
const int  nc[]
[static]
 

Calculate the length of the parameter vector, which is called coef for historical reasons.

Parameters:
nf number of factors
nc number of columns in the model matrices for each factor
Returns:
total length of the coefficient vector

SEXP ldl_inverse SEXP  x  )  [static]
 

Update the diagonal blocks of the inverse of LDL' (=Z'Z+W). The lower Cholesky factors of the updated blocks are stored in the bVar slot.

Parameters:
x pointer to an ssclme object
Returns:
R_NilValue (x is updated in place)

int ldl_update_ind int  probe,
int  start,
const int  ind[]
[static]
 

Return the position of probe in the sorted index vector ind. It is known that the position is greater than or equal to start so a linear search from start is used.

Parameters:
probe value to be matched
start index at which to start
ind vector of indices
Returns:
index of the entry matching probe

void ssclme_calc_maxod int  n,
int  Parent[]
[static]
 

Calculate and store the maximum number of off-diagonal elements in the inverse of L, based on the elimination tree. The maximum is itself stored in the Parent array. (FIXME: come up with a better design.)

Parameters:
n number of columns in the matrix
Parent elimination tree for the matrix

SEXP ssclme_coef SEXP  x  ) 
 

Extract the upper triangles of the Omega matrices. These aren't "coefficients" but the extractor is called coef for historical reasons. Within each group these values are in the order of the diagonal entries first then the strict upper triangle in row order.

Parameters:
x pointer to an ssclme object
Returns:
numeric vector of the values in the upper triangles of the Omega matrices

SEXP ssclme_coefGets SEXP  x,
SEXP  coef
 

Assign the upper triangles of the Omega matrices. (Called coef for historical reasons.)

Parameters:
x pointer to an ssclme object
coef pointer to an numeric vector of appropriate length
Returns:
R_NilValue

SEXP ssclme_coefGetsUnc SEXP  x,
SEXP  coef
 

Assign the Omega matrices from the unconstrained parameterization.

Parameters:
x pointer to an ssclme object
coef pointer to an numeric vector of appropriate length
Returns:
R_NilValue

SEXP ssclme_coefUnc SEXP  x  ) 
 

Extract the unconstrained parameters that determine the Omega matrices. (Called coef for historical reasons.) The unconstrained parameters are derived from the LDL' decomposition of Omega_i. The first nc[i] entries in each group are the diagonals of log(D) followed by the strict lower triangle of L in column order.

Parameters:
x pointer to an ssclme object
Returns:
numeric vector of unconstrained parameters that determine the Omega matrices

SEXP ssclme_collapse SEXP  x  ) 
 

Copy an ssclme object collapsing the fixed effects slots to the response only.

Parameters:
x pointer to an ssclme object
Returns:
a duplicate of x with the fixed effects slots collapsed to the response only

void ssclme_copy_ctab int  nf,
const int  nc[],
SEXP  ctab,
SEXP  ssc
[static]
 

Using the sscCrosstab object from the grouping factors, generate the slots in an ssclme object related to the symmetric sparse matrix representation of Z'Z. If the model matrices for the grouping factors have only one column each then the structure can be copied, otherwise it must be generated from the sscCrosstab and the number of columns per grouping factor.

Parameters:
nf number of factors
nc vector of length nf+2 with number of columns in model matrices
ctab pointer to the sscCrosstab object
ssc pointer to an ssclme object to be filled out

SEXP ssclme_create SEXP  facs,
SEXP  ncv
 

Create an ssclme object from a list of grouping factors, sorted in order of non-increasing numbers of levels, and an integer vector of the number of columns in the model matrices. There is one more element in ncv than in facs. The last element is the number of columns in the model matrix for the fixed effects plus the response. (i.e. p+1)

Parameters:
facs pointer to a list of grouping factors
ncv pointer to an integer vector of number of columns per model matrix
Returns:
pointer to an ssclme object

SEXP ssclme_EMsteps SEXP  x,
SEXP  nsteps,
SEXP  REMLp,
SEXP  verb
 

Perform a number of ECME steps for the REML or ML criterion.

Parameters:
x pointer to an ssclme object
nsteps pointer to an integer scalar giving the number of ECME steps to perform
REMLp pointer to a logical scalar indicating if REML is to be used
verb pointer to a logical scalar indicating verbose mode
Returns:
NULL

SEXP ssclme_factor SEXP  x  ) 
 

If status[["factored"]] is FALSE, create and factor Z'Z+Omega, then create RZX and RXX, the deviance components, and the value of the deviance for both ML and REML.

Parameters:
x pointer to an ssclme object
Returns:
NULL

SEXP ssclme_fitted SEXP  x,
SEXP  facs,
SEXP  mmats,
SEXP  useRf
 

Calculate and return the fitted values.

Parameters:
x pointer to an ssclme object
facs list of grouping factors
mmats list of model matrices
useRf pointer to a logical scalar indicating if the random effects should be used
Returns:
pointer to a numeric array of fitted values

SEXP ssclme_fixef SEXP  x  ) 
 

Extract the conditional estimates of the fixed effects

Parameters:
x Pointer to an ssclme object
Returns:
a numeric vector containing the conditional estimates of the fixed effects

SEXP ssclme_gradient SEXP  x,
SEXP  REMLp,
SEXP  Uncp
 

Return the gradient of the ML or REML deviance.

Parameters:
x pointer to an ssclme object
REMLp pointer to a logical scalar indicating if REML is to be used
Uncp pointer to a logical scalar indicating if the unconstrained parameterization is to be used
Returns:
pointer to a numeric vector of the gradient.

SEXP ssclme_inflate_and_factor SEXP  x  ) 
 

Inflate Z'Z according to Omega and create the factorization LDL'

Parameters:
x pointer to an ssclme object
Returns:
NULL

SEXP ssclme_initial SEXP  x  ) 
 

Create and insert initial values for Omega_i.

Parameters:
x pointer to an ssclme object
Returns:
NULL

SEXP ssclme_invert SEXP  x  ) 
 

If necessary, factor Z'Z+Omega, ZtX, and XtX then, if necessary, form RZX, RXX, and bVar for the inverse of the Cholesky factor.

Parameters:
x pointer to an ssclme object
Returns:
NULL (x is updated in place)

SEXP ssclme_ranef SEXP  x  ) 
 

Extract the conditional modes of the random effects.

Parameters:
x Pointer to an ssclme object
Returns:
a vector containing the conditional modes of the random effects

SEXP ssclme_sigma SEXP  x,
SEXP  REML
 

Extract the ML or REML conditional estimate of sigma

Parameters:
x pointer to an ssclme object
REML logical scalar - TRUE if REML estimates are requested
Returns:
numeric scalar

SEXP ssclme_to_lme SEXP  call,
SEXP  facs,
SEXP  x,
SEXP  model,
SEXP  REML,
SEXP  rep,
SEXP  fitted,
SEXP  residuals
 

Create an lme object from its components. This is not done by new("lme", ...) at the R level because of the possibility of causing the copying of very large objects.

Parameters:
call Pointer to the original call
facs pointer to the list of grouping factors
x pointer to the model matrices (may be of length zero)
model pointer to the model frame
REML pointer to a logical scalar indicating if REML is used
rep pointer to the converged ssclme object
fitted pointer to the fitted values
residuals pointer to the residuals
Returns:
an lme object

SEXP ssclme_transfer_dimnames SEXP  x,
SEXP  facs,
SEXP  mmats
 

Copy the dimnames from the list of grouping factors and the model matrices for the grouping factors into the appropriate parts of the ssclme object.

Parameters:
x pointer to an ssclme object
facs pointer to a list of factors
mmats pointer to a list of model matrices
Returns:
NULL

SEXP ssclme_update_mm SEXP  x,
SEXP  facs,
SEXP  mmats
 

Update the numerical entries x, ZtX, and XtX in an ssclme object according to a set of model matrices.

Parameters:
x pointer to an ssclme object
facs pointer to a list of grouping factors
mmats pointer to a list of model matrices
Returns:
NULL

SEXP ssclme_variances SEXP  x  ) 
 

Return the unscaled variances

Parameters:
x pointer to an ssclme object
Returns:
a list similar to the Omega list with the unscaled variances


Generated on Mon May 31 14:04:37 2004 for Matrix by doxygen 1.3.6-20040222