Conditional intraclass correlation (ICC) from a linear mixed model as evidence for test-retest reliability? In my experiment with two conditions (between-subjects design), participants completed a single-item scale three times: (1) before the experimental manipulation, (2) after the experimental manipulation, and (3) at the end of the study. I want to test the test-retest reliability of my single-item scale over three measurements.
The issue here is that the experimental manipulation is known to affect the response in the measure. I'm assuming that I need to control for any effect of the manipulation, so my guess is to use a linear mixed model (hierarchical linear model) with measurements nested within participants and include an experimental condition as a factor. This way, I can get an intraclass correlation (ICC) from this conditional model.
Would it be possible to interpret this conditional ICC as "a measure of test-retest reliability after controlling for the experimental effect"? Are there any existing studies using a similar approach?
 A: This post really helped me and I wanted to thank you.
In case other users ran to the same issue I had - I am adding a slight change to the simulation above. The only thing here is that this shows that Pearson corr for two times measurements is exactly the same as $\rho$. Nothing special - only nice to see the numbers match :) Also, ever so slight correction in the participant vector to make this work.
Cheers
Nitzan
set.seed(807)

npart=1000 # number of particpants
ntime=2   # numer of observations (timepoints) per participant
mu=2.5     # mean value on the Likert item
sdp=1      # standard deviation of participant random effect (variance==1)
sd=.7071   # standard deviation of within participant (residual; variance = .5)

participant = rep(rep(1:npart, each = nobs),ntime)  # creating 1000 participants w/ 3 repeats
participant = participant[order(participant)]
time        = rep(rep(1:ntime, each=1),1000)        # creating a time variable

parteff = rnorm(npart, 0, sdp)     # drawing from normal for participant deviation
parteff = rep(parteff, each=ntime) # ensuring participant effect is same for three observations

timeeff = rnorm(npart*ntime, 0, sd) # drawing from normal for within-participant residual

dat=data.frame(participant, time, parteff, timeeff) # create data frame

dat$resp = with(dat, mu + parteff + timeeff ) # creating response for each individual

#Variance components model
library(lme4)

m1 <- lmer(resp ~ 1 + (1|participant), dat)
summary(m1) # estimates close to simulated values

#calculate pearson corr
library(reshape2)
df.wide   <-dcast(dat,participant~time,mean,value.var='resp')[,-1]
cor(df.wide)

#get the same from the HLM fit
print(VarCorr(m1))
.95478^2/(.95478^2+0.74685^2)
 ```

