Learn R Programming

multiway (version 1.0-7)

mcr: Multiway Covariates Regression

Description

Fits Smilde and Kiers's Multiway Covariates Regression (MCR) model to connect a 3-way predictor array and a 2-way response array that share a common mode. Parameters are estimated via alternating least squares with optional constraints.

Usage

mcr(X, Y, nfac = 1, alpha = 0.5, nstart = 10, 
    model = c("parafac", "parafac2", "tucker"),
    const = NULL, control = NULL, weights = NULL,
    Afixed = NULL, Bfixed = NULL, Cfixed = NULL, Dfixed = NULL,
    Astart = NULL, Bstart = NULL, Cstart = NULL, Dstart = NULL,
    Astruc = NULL, Bstruc = NULL, Cstruc = NULL, Dstruc = NULL,
    Amodes = NULL, Bmodes = NULL, Cmodes = NULL, Dmodes = NULL,
    maxit = 500, ctol = 1e-4, parallel = FALSE, cl = NULL, 
    output = c("best", "all"), verbose = TRUE, backfit = FALSE)

Value

A

Predictor A weight matrix.

B

Predictor B weight matrix.

C

Common C weight matrix.

D

Response D weight matrix.

W

Coefficients. See Note.

LOSS

Value of LOSS function.

SSE

Sum of Squared Errors for X and Y.

Rsq

R-squared value for X and Y.

iter

Number of iterations.

cflag

Convergence flag. See Note.

model

See argument model.

const

See argument const.

control

See argument control.

weights

See argument weights.

alpha

See argument alpha.

fixed

Logical vector indicating whether 'fixed' weights were used for each matrix.

struc

Logical vector indicating whether 'struc' constraints were used for each matrix.

Phi

Mode A crossproduct matrix. Only if model = "parafac2".

G

Core array. Only if model = "tucker".

Arguments

X

Three-way predictor array with dim = c(I,J,K).

Y

Two-way response array with dim = c(K,L).

nfac

Number of factors.

alpha

Tuning parameter between 0 and 1.

nstart

Number of random starts.

model

Model for X. Defaults to "parafac".

const

Character vector of length 4 giving the constraints for A, B, C, and D (defaults to unconstrained). See const for the 24 available options. Ignored if model = "tucker".

control

List of parameters controlling options for smoothness constraints. This is passed to const.control, which describes the available options.

weights

Vector of length K giving non-negative weights for fitting via weighted least squares. Defaults to vector of ones.

Afixed

Used to fit model with fixed Mode A weights.

Bfixed

Used to fit model with fixed Mode B weights.

Cfixed

Used to fit model with fixed Mode C weights.

Dfixed

Used to fit model with fixed Mode D weights.

Astart

Starting Mode A weights. Default uses random weights.

Bstart

Starting Mode B weights. Default uses random weights.

Cstart

Starting Mode C weights. Default uses random weights.

Dstart

Starting Mode D weights. Default uses random weights.

Astruc

Structure constraints for Mode A weights. See Note.

Bstruc

Structure constraints for Mode B weights. See Note.

Cstruc

Structure constraints for Mode C weights. Ignored.

Dstruc

Structure constraints for Mode D weights. See Note.

Amodes

Mode ranges for Mode A weights (for unimodality constraints). See Note.

Bmodes

Mode ranges for Mode B weights (for unimodality constraints). See Note.

Cmodes

Mode ranges for Mode C weights (for unimodality constraints). Ignored.

Dmodes

Mode ranges for Mode D weights (for unimodality constraints). See Note.

maxit

Maximum number of iterations.

ctol

Convergence tolerance (R^2 change).

parallel

Logical indicating if parLapply should be used. See Examples.

cl

Cluster created by makeCluster. Only used when parallel=TRUE.

output

Output the best solution (default) or output all nstart solutions.

verbose

If TRUE, fitting progress is printed via txtProgressBar. Ignored if parallel=TRUE.

backfit

Should backfitting algorithm be used for cmls?

Author

Nathaniel E. Helwig <helwig@umn.edu>

Details

Given a predictor array X = array(x, dim=c(I,J,K)) and a response matrix Y = matrix(y, nrow=K, ncol=L), the multiway covariates regression (MCR) model assumes a tensor model for X and a bilinear model for Y, which are linked through a common C weight matrix. For example, using the Parafac model for X, the MCR model has the form

X[i,j,k] = sum A[i,r]*B[j,r]*C[k,r] + Ex[i,j,k]
and
Y[k,l] = sum C[k,r]*D[l,r] + Ey[k,l]

Parameter matrices are estimated by minimizing the loss function

LOSS = alpha * (SSE.X / SSX) + (1 - alpha) * (SSE.Y / SSY)

where

SSE.X = sum((X - Xhat)^2)
SSE.Y = sum((Y - Yhat)^2)
SSX = sum(X^2)
SSY = sum(Y^2)

When weights are input, SSE.X, SSE.Y, SSX, and SSY are replaced by the corresponding weighted versions.

References

Smilde, A. K., & Kiers, H. A. L. (1999). Multiway covariates regression models, Journal of Chemometrics, 13, 31-48. tools:::Rd_expr_doi("10.1002/(SICI)1099-128X(199901/02)13:1<31::aid-cem528>3.0.CO;2-P")

See Also

The fitted.mcr function creates the model-implied fitted values from a fit "mcr" object.

The resign.mcr function can be used to resign factors from a fit "mcr" object.

The rescale.mcr function can be used to rescale factors from a fit "mcr" object.

The reorder.mcr function can be used to reorder factors from a fit "mcr" object.

The cmls function (from CMLS package) is called as a part of the alternating least squares algorithm.

See parafac, parafac2, and tucker for more information about the Parafac, Parafac2, and Tucker models.

Examples

Run this code

##########   multiway covariates regression   ##########

# create random data array with Parafac structure
set.seed(3)
mydim <- c(10, 20, 100)
nf <- 2
Amat <- matrix(rnorm(mydim[1]*nf), mydim[1], nf)
Bmat <- matrix(rnorm(mydim[2]*nf), mydim[2], nf)
Cmat <- matrix(rnorm(mydim[3]*nf), mydim[3], nf)
Xmat <- tcrossprod(Amat, krprod(Cmat, Bmat))
Xmat <- array(Xmat, dim = mydim)
EX <- array(rnorm(prod(mydim)), dim = mydim)
EX <- nscale(EX, 0, ssnew = sumsq(Xmat))   # SNR = 1
X <- Xmat + EX

# create response array
ydim <- c(mydim[3], 4)
Dmat <- matrix(rnorm(ydim[2]*nf), ydim[2], nf)
Ymat <- tcrossprod(Cmat, Dmat)
EY <- array(rnorm(prod(ydim)), dim = ydim)
EY <- nscale(EY, 0, ssnew = sumsq(Ymat))   # SNR = 1
Y <- Ymat + EY

# fit MCR model
mcr(X, Y, nfac = nf, nstart = 1)
mcr(X, Y, nfac = nf, nstart = 1, model = "parafac2")
mcr(X, Y, nfac = nf, nstart = 1, model = "tucker")



if (FALSE) {

##########   parallel computation   ##########

# create random data array with Parafac structure
set.seed(3)
mydim <- c(10, 20, 100)
nf <- 2
Amat <- matrix(rnorm(mydim[1]*nf), mydim[1], nf)
Bmat <- matrix(rnorm(mydim[2]*nf), mydim[2], nf)
Cmat <- matrix(rnorm(mydim[3]*nf), mydim[3], nf)
Xmat <- tcrossprod(Amat, krprod(Cmat, Bmat))
Xmat <- array(Xmat, dim = mydim)
EX <- array(rnorm(prod(mydim)), dim = mydim)
EX <- nscale(EX, 0, ssnew = sumsq(Xmat))   # SNR = 1
X <- Xmat + EX

# create response array
ydim <- c(mydim[3], 4)
Dmat <- matrix(rnorm(ydim[2]*nf), ydim[2], nf)
Ymat <- tcrossprod(Cmat, Dmat)
EY <- array(rnorm(prod(ydim)), dim = ydim)
EY <- nscale(EY, 0, ssnew = sumsq(Ymat))   # SNR = 1
Y <- Ymat + EY

# fit MCR-Parafac model (10 random starts -- sequential computation)
set.seed(1)
system.time({mod <- mcr(X, Y, nfac = nf)})
mod

# fit MCR-Parafac model (10 random starts -- parallel computation)
cl <- makeCluster(detectCores())
ce <- clusterEvalQ(cl, library(multiway))
clusterSetRNGStream(cl, 1)
system.time({mod <- mcr(X, Y, nfac = nf, parallel = TRUE, cl = cl)})
mod
stopCluster(cl)

}


Run the code above in your browser using DataLab