Thanks for following up. I’m not sure how to use the posterior_linpred() function to generate people’s factor scores. Below, I’ll provide a minimal reproducible example:

```
library("rstan")
library("brms")
library("mirt")
library("dplyr")
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
Sys.setenv(LOCAL_CPPFLAGS = '-march=native')
set.seed(1234)
a <- matrix(rlnorm(20,.2,.3))
diffs <- t(apply(matrix(runif(20*2, .3, 1), 20), 1, cumsum))
diffs <- -(diffs - rowMeans(diffs))
d <- diffs + rnorm(20)
itemData <- simdata(a, d, 2400, itemtype = 'graded')
participantData <- expand.grid(ID = 1:200, wave = 1:4, rater = 1:3)
myData <- cbind(participantData, itemData)
myData[which(myData$ID %% 2 == 0 & myData$wave == 2), 4:ncol(myData)] <- NA
myData[which(myData$ID %% 2 != 0 & myData$wave == 3), 4:ncol(myData)] <- NA
myData_long <- myData %>%
pivot_longer(-c("ID","wave","rater"), names_to = "item", values_to = "response")
myData_long$ID <- as.factor(myData_long$ID)
myData_long$wave <- as.factor(myData_long$wave)
myData_long$rater <- as.factor(myData_long$rater)
myData_long$item <- as.factor(myData_long$item)
myData_long$response <- factor(myData_long$response, levels = c(0, 1, 2), ordered = TRUE)
formula_va_ord_2pl <- bf(
response ~ 1 + rater + (1 |i| item) + (1 |w| wave) + (1 |p| ID),
disc ~ 1 + rater + (1 |i| item) + (1 |w| wave) + (1 |p| ID)
)
mixedEffectsGRM <- brm(
formula = formula_va_ord_2pl,
data = myData_long,
family = brmsfamily("cumulative", "logit"),
chains = 1,
iter = 100,
seed = 1234
)
mixedEffectsGRM_fitted <- fitted(mixedEffectsGRM)
mixedEffectsGRM_predict <- predict(mixedEffectsGRM)
mixedEffectsGRM_posterior1 <- posterior_linpred(mixedEffectsGRM)
mixedEffectsGRM_posterior2 <- posterior_linpred(mixedEffectsGRM, newdata = myData_long)
dim(mixedEffectsGRM_fitted) #36000 x 4 x 3
dim(mixedEffectsGRM_predict) #36000 x 3
dim(mixedEffectsGRM_posterior1) #50 x 36000
dim(mixedEffectsGRM_posterior2) #50 x 36000
nrow(unique(myData[,c("ID","wave")])) #800 ID-wave combinations
```

First, I simulate the data. Then I run the item response model in brms. Then I check the output for various functions, including fitted(), predict(), and posterior_linpred(). However, the output from all of these functions leads to output that is much larger than expected. I’m looking to obtain the factor scores for all 800 ID-wave combinations. Therefore, my final output should be a single vector of length 800. It seems that the fitted values are generated at the item-level, but I’m looking to output factor scores at the person-wave level. Paul showed how to do this at the person-level in his paper on item response modeling by using the ranef() function (the reference is below). However, I’m not sure how to do generate the factor scores at the person-wave level. Your help is greatly appreciated!

Bürkner, P.-C. (2019). Bayesian Item Response Modelling in R with brms and Stan. arXiv. doi: 1905.09501