Smooth PRC using LMMsolver (Boer, 2023)

Smooth PRC using LMMsolver (Boer, 2023)

smoothPRC(
  Y,
  lmm_model,
  options_iter = list(b_init = NULL, k = 3, mA = 4, tol = 1e-08, maxiter = 50),
  n_axes = 2
)

Arguments

Y

matrix, data frame of the abundance data (dimension sites x species)

lmm_model

The result of set_smoothPRC_model.

options_iter

list of five elements

• b_init: initial values of the species loadings. Default NULL.

• k: number of cycles in block-power updating to better separate the first eigen vector from the second. Default NULL.

• mA: number of vectors in Anderson acceleration. Default 3.

• tol: stopping criterion, tolarance. Default 1e-8.

• maxiter:maximum number of iterations. Default 50.

n_axes

the number of axes to extract. Default 2.

Value

A list. The eigenvalues are in result$eig.

Details

An atypical aspect of the function is that the time and treatment factors need to have the names Time and Treatment in data.

Classical PRC is a reduced rank model of rank 1 (or 2) with model Treatment*Time with the Time factor is a concomitant variable. In smoothPRC the factor Time is taken quantitatively and the main effect of time is not an un-ordered set of values but a smooth P-spline model. Also, the time-dependent effects of the treatments are smooth P-splines. As the main and interaction effects are splines, the main effect of time cannot be adjusted for in the usual way. It is taken out via a separate spline fit using the formulas in splineH0 and fixedH0 of set_smoothPRC_model.

The algorithm solves a penalized reduced rank model for the dominant eigen vector (i.e. rank 1). The model is a P-spline model and the penalty parameters are obtained by REML using LMMsolve (Boer 2023). The model is linear in, possibly tensor, B-splines, but through the B-splines non-linear from the perspective of the user. With quadratic ridge-type difference penalties on the coefficients the B-splines becomeP-splines (penalized B-splines). For an initial b of species scores, the unconstrainted sites Yb are calculated, which are then used as response vector in LMMsolve to give a constrained site scores x by the P-spline model specified in lmm_model. This is part of a power algorithm with block and Anderson acceleration to solve the dominant eigen vector of the underlying eigen problem while optimizing simultaneously the penalties using LMMsolve.

If lmm_model$scaling =="ms", then in the result B in the result is multiplied by mult = sqrt(ncol(Y)) and all other items are are divided by mult = sqrt(ncol(Y)), except Yb to be inline with fitted(result$obj). Note that Yb, a column in YB, is the response in LMMsolve model.

References

Boer, Martin P. 2023. Tensor Product P-Splines Using a Sparse Mixed Model Formulation. Statistical Modelling 23 (5-6): 465–79. doi:10.1177/1471082X231178591