## ----setup, include = FALSE-----------------------------------------------------------------------
knitr::opts_chunk$set(echo = TRUE, fig.align = "center", out.width = "100%")
knitr::opts_chunk$set(comment = "#>")
old_options <- options(width = 100)

## -------------------------------------------------------------------------------------------------
library(spifa)

data(ipixuna)
nfactors <- 3
nitems <- ncol(ipixuna$items)

## -------------------------------------------------------------------------------------------------
A <- matrix(1, nitems, nfactors)
A[c(4, 8), 1] <- 0
A[c(4, 5, 6, 7, 8, 10), 2] <- 0
A[c(5, 6), 3] <- 0
A

## -------------------------------------------------------------------------------------------------
ngp <- 2
T_loading <- matrix(c(1, 1, 0, 0, 0, 1), nfactors, ngp)
T_loading

## -------------------------------------------------------------------------------------------------
sd_fixed <- rep(0.2, nfactors)
sd_fixed

## -------------------------------------------------------------------------------------------------
# easiness hyperparameters
c_mean <- seq(-0.7, 0.9, length.out = nitems)
c_sd <- rep(0.3, nitems)

# discrimination hyperparameters
A_mean <- matrix(0, nitems, nfactors)
A_mean[c(5, 6), 1] <- 1
A_mean[c(4, 8), 3] <- -1
A_sd <- matrix(1, nitems, nfactors)
A_sd[A_mean != 0] <- 0.45

# predictor effects hyperparameters
B_mean <- matrix(c(-0.5, 0.5, 0.2), nrow = 1)
B_sd <- matrix(0.5, nrow = 1, ncol = nfactors)

# correlation initial value
corr_initial <- diag(nfactors)

# range parameter hyperparameters
phi_mean <- 150
phi_sd <- 0.4

## -------------------------------------------------------------------------------------------------
set.seed(1)
fit_informed <- spifa(
  items ~ wealth, data = ipixuna, nfactors = nfactors, ngp = 3,
  burnin = 200, niter = 200, thin = 2,
  constraints = list(discrimination = A, loading = diag(nfactors), sd = sd_fixed),
  priors = list(
    easiness = list(initial = c_mean, mean = c_mean, sd = c_sd),
    discrimination = list(initial = A_mean, mean = A_mean, sd = A_sd),
    effect = list(initial = B_mean, mean = B_mean, sd = B_sd),
    corr = list(initial = corr_initial, eta = 2),
    loading = list(initial = 0.6, mean = 0.6, sd = 0.3),
    range = list(initial = phi_mean, mean = phi_mean, sd = phi_sd)))
fit_informed

## ----warning = FALSE------------------------------------------------------------------------------
fit1 <- update(fit_informed, niter = 1000)
fit1

## ----warning = FALSE------------------------------------------------------------------------------
fit2 <- spifa(items ~ wealth, data = ipixuna, nfactors = 2, niter = 1000, burnin = 1000)
fit2

## -------------------------------------------------------------------------------------------------
# M1: 3 factors
dic(fit1, burnin = 100, thin = 2)

## -------------------------------------------------------------------------------------------------
# M2: 2 factors
dic(fit2, burnin = 100, thin = 2)

## ----include = FALSE----------------------------------------------------------
options(old_options)