Load xxM and data

library(xxm)
data(faces.xxm, package = "xxm")
# Two observed responses at level 1: PA(physical attractiveness),
# SYM(symmetry) Level 1 responses are cross-classified within Rater(levek 2)
# and Target(level 3) faces.xxm contains 3 data frames: 'response' (level1),
# rater (level2), target (level3)

Define local variables for conveinence

# 

Construct R-matrices

For each parameter matrix, construct three related matrices:

  1. pattern matrix: A matrix indicating free or fixed parameters.
  2. value matrix: with start or fixed values for corresponding parameters.
  3. label matrix: with user friendly label for each parameter. label matrix is optional.
## WITHIN-LEVEL MATRICES Model Matrices for Level 1 (response) Response:
## Observed residual covariance matrix
l1_theta_pat <- matrix(1, 2, 2)
l1_theta_val <- diag(3, 2)
# Response: Observed variable intercepts
l1_nu_pat <- matrix(c(1, 1), 2, 1)
l1_nu_val <- matrix(c(4.2671, 5.2066), 2, 1)

# Model Matrices for Level 2 (rater) Rater: factor-covariance matrix
l2_psi_pat <- matrix(1, 2, 2)
l2_psi_val <- diag(0.1, 2)

# Model Matrices for Level 3 (target) Target: factor-covariance matrix
l3_psi_pat <- matrix(1, 2, 2)
l3_psi_val <- diag(0.1, 2)

## ACROSS-LEVEL MATRICES Rater -> Response factor-loading matrix
l1_l2_lambda_pat <- matrix(0, 2, 2)
l1_l2_lambda_val <- matrix(c(1, 0, 0, 1), 2, 2)
# Target -> Response factor-loading matrix
l1_l3_lambda_pat <- matrix(0, 2, 2)
l1_l3_lambda_val <- matrix(c(1, 0, 0, 1), 2, 2)

Construct main model object

xxmModel() is used to declare level names. The function returns a model object that is passed as a parameter to subsequent statements.

xm <- xxmModel(levels = c("response", "rater", "target"))
## A new model was created.

Add submodels to the model objects

For each declared level xxmSubmodel() is invoked to add corresponding submodel to the model object. The function adds three pieces of information: 1. parents declares a list of parents of the current level. 2. variables declares names of observed dependent (ys), observed independent (xs) and latent variables (etas) for the level. 3. data R data object for the current level.

# Level 1: Response submodel
xm <- xxmSubmodel(model = xm, level = "response", parents = c("rater", "target"), 
    ys = c("PA", "SYM"), xs = , etas = , data = faces.response, siblings = )
## Submodel for level `response` was created.
# Level 2: Rater submodel
xm <- xxmSubmodel(model = xm, level = "rater", parents = , ys = , xs = , etas = c("eta1", 
    "eta2"), data = faces.rater, siblings = )
## Submodel for level `rater` was created.
# Level 3: Target submodel
xm <- xxmSubmodel(model = xm, level = "target", parents = , ys = , xs = , etas = c("eta3", 
    "eta4"), data = faces.target, siblings = )
## Submodel for level `target` was created.

Add Within-level parameter matrices for each submodel

For each declared level xxmWithinMatrix() is used to add within-level parameter matrices. For each parameter matrix, the function adds the three matrices constructed earlier:

## Level 1 (Response) within matrices (theta and nu)
xm <- xxmWithinMatrix(model = xm, level = "response", type = "theta", pattern = l1_theta_pat, 
    value = l1_theta_val, label = , name = )
## 
## 'theta' matrix does not exist and will be added.
##  Added `theta` matrix.

xm <- xxmWithinMatrix(model = xm, level = "response", type = "nu", pattern = l1_nu_pat, 
    value = l1_nu_val, label = , name = )
## 
## 'nu' matrix does not exist and will be added.
##  Added `nu` matrix.
## Level 2 (Rater) within matrices (psi)
xm <- xxmWithinMatrix(model = xm, level = "rater", type = "psi", pattern = l2_psi_pat, 
    value = l2_psi_val, label = , name = )
## 
## 'psi' matrix does not exist and will be added.
##  Added `psi` matrix.
## Level 3 (Target) within matrices (psi)
xm <- xxmWithinMatrix(model = xm, level = "target", type = "psi", pattern = l3_psi_pat, 
    value = l3_psi_val, label = , name = )
## 
## 'psi' matrix does not exist and will be added.
##  Added `psi` matrix.

Add Across-level parameter matrices to the model

Pairs of levels that share parent-child relationship have regression relationships. xxmBetweenMatrix() is used to add corresponding parameter matrices connecting the two levels.

For each parameter matrix, the function adds the three matrices constructed earlier:

# Rater -> Response factor-loading matrix
xm <- xxmBetweenMatrix(model = xm, parent = "rater", child = "response", type = "lambda", 
    pattern = l1_l2_lambda_pat, value = l1_l2_lambda_val, label = , name = )
## 
## 'lambda' matrix does not exist and will be added.
##  Added `lambda` matrix.
# Target -> Response factor-loading matrix
xm <- xxmBetweenMatrix(model = xm, parent = "target", child = "response", type = "lambda", 
    pattern = l1_l3_lambda_pat, value = l1_l3_lambda_val, label = , name = )
## 
## 'lambda' matrix does not exist and will be added.
##  Added `lambda` matrix.

Estimate model parameters

Estimation process is initiated by xxmRun(). If all goes well, a quick printed summary of results is produced.

xm <- xxmRun(xm)
## Model specification seems reasonable.
## Validity: Checked
## Matrix Constructed
## ------------------------------------------------------------------------------
## Estimating model parameters
## ------------------------------------------------------------------------------
##               26431.8098017981 
##               26371.4492959767 
##               26072.0902076243 
##               25792.9179129407 
##               25773.1848346248 
##               25752.1331005728 
##               25748.4539827775 
##               25729.8795218688 
##               25707.2225123686 
##               25704.3121614445 
##               25702.6115933138 
##               25684.0215206468 
##               25683.5860533049 
##               25681.1605990959 
##               25677.4380503965 
##               25672.7694135890 
##               25671.5746114201 
##               25667.4435282819 
##               25659.4221562202 
##               25654.1276371467 
##               25648.8632487349 
##               25643.8122398250 
##               25635.0140061606 
##               25630.8073218753 
##               25626.7287888638 
##               25624.8326060910 
##               25622.9430818210 
##               25622.3331290390 
##               25621.6295114769 
##               25621.1621980843 
##               25620.8004991919 
##               25620.7407825865 
##               25620.7338131240 
##               25620.7332302781 
##               25620.7332222210 
##               25620.7332173246 
## Model converged normally
## ------------------------------------------------------------------------------
## *
##  1:                            response_theta_1_1 ::      1.772 [     0.000,      0.000]
## 
##  2:                            response_theta_1_2 ::      0.624 [     0.000,      0.000]
## 
##  3:                            response_theta_2_2 ::      2.288 [     0.000,      0.000]
## 
##  4:                               response_nu_1_1 ::      4.267 [     0.000,      0.000]
## 
##  5:                               response_nu_2_1 ::      5.207 [     0.000,      0.000]
## 
##  6:                                 rater_psi_1_1 ::      1.130 [     0.000,      0.000]
## 
##  7:                                 rater_psi_1_2 ::      0.462 [     0.000,      0.000]
## 
##  8:                                 rater_psi_2_2 ::      1.011 [     0.000,      0.000]
## 
##  9:                                target_psi_1_1 ::      0.916 [     0.000,      0.000]
## 
## 10:                                target_psi_1_2 ::      0.456 [     0.000,      0.000]
## 
## 11:                                target_psi_2_2 ::      0.327 [     0.000,      0.000]
## 
## ------------------------------------------------------------------------------

Estimate profile-likelihood confidence intervals

Once parameters are estimated, confidence inetrvals are estimated by invoking xxmCI() . Depending on the the number of observations and the complexity of the dependence structure xxmCI() may take very long. xxMCI() displays a table of parameter estimates and CIS.

View results

A summary of results may be retrived as an R list by a call to xxmSummary()

s <- xxmSummary(xm)

Free moodel object

xxM model object may hog a large amount of RAM outside of R's memory. This memory will automatically be released, when R's workspace is cleared by a call to rm(list=ls()) or at the end of the R session. Alternatively, xxmFree() may be called to release memory.

xxmFree(xm)
##          used (Mb) gc trigger (Mb) max used (Mb)
## Ncells 259095  7.0     407500 10.9   350000  9.4
## Vcells 247901  1.9     786432  6.0   786320  6.0
rm(list = ls())