Estimating people's factor scores (theta) from item response model

I assumed that you wanted factor scores for each ID, split by wave.

Here is what I meant in more detail.
(sorry, I am a data.table person and don’t use the tidyverse much, hope this can still be understood)

First transpose the results from posterior_linpred, make a data.table, and rename columns, and merge with myData_long.

library(data.table)
tmp = data.table(t(mixedEffectsGRM_posterior1))
setnames(tmp,names(tmp), paste0("iter.",1:50))
complete.myData_long = data.table(myData_long)[complete.cases(myData_long)]
plp_wide = cbind(complete.myData_long, tmp)

Next make really long format, that has also iterations in long format.

plp_long = 
melt(plp_wide,
     id.vars = names(myData_long),
     variable.name = "iter")

Now we can aggregate by ID, wave and iteration:

plp_long.ID.wave.iter = 
plp_long[,list(latent = mean(value)),
          by = c("ID","wave","iter")]

With this in hand, we can for example calculate means, and CIs by wave and ID:

plp.long.ID.wave = 
plp_long_ID_wave_iter[,list(m = mean(latent),
                            CI05 = quantile(latent,.05),
                            CI95 = quantile(latent,.95)),
                        by = c("ID","wave")]
1 Like