#install polychor library
library(polycor)
#read in transition probabilites from 3 attribute example
tprobs = read.table('transprobs3.txt', header = T)
#read the reliability function into the environment
growth_reliability = function(transprobs, numtraits){
rel1 = matrix(NA, nrow = numtraits)
rel2 = matrix(NA, nrow = numtraits)
#loop over each trait
for(a in 1:numtraits){
probs00 = transprobs[,1+4*(a-1)]
probs01 = transprobs[,2+4*(a-1)]
probs10 = transprobs[,3+4*(a-1)]
probs11 = transprobs[,4+4*(a-1)]
all = as.matrix(cbind(probs00, probs01, probs10, probs11))
#create 4x4 contingency table of experiment-re-experiment probabilities
cell11 = probs00^2
cell22 = probs01^2
cell33 = probs10^2
cell44 = probs11^2
cell12 = probs00*probs01
cell13 = probs00*probs10
cell14 = probs00*probs11
cell23 = probs01*probs10
cell24 = probs01*probs11
cell34 = probs10*probs11
fourmat = matrix(NA, nrow = 4, ncol = 4)
fourmat[1,1] = mean(cell11)
fourmat[2,2] = mean(cell22)
fourmat[3,3] = mean(cell33)
fourmat[4,4] = mean(cell44)
fourmat[1,2] = fourmat[2,1] = mean(cell12)
fourmat[1,3] = fourmat[3,1] = mean(cell13)
fourmat[1,4] = fourmat[4,1] = mean(cell14)
fourmat[2,3] = fourmat[3,2] = mean(cell23)
fourmat[2,4] = fourmat[4,2] = mean(cell24)
fourmat[3,4] = fourmat[4,3] = mean(cell34)
#compute polychoric correlation of the aggregate table
rel1[a] = polychor(fourmat)
#obtain most likely transition for each examinee
pmax = matrix(NA, nrow = nrow(transprobs))
for(j in 1:nrow(transprobs)){
pmax[j] = max(all[j,])
}
#compute average of most likely class probabilities
rel2[a] = mean(pmax)
}
results = as.matrix(cbind(rel1, rel2))
#column 1 is metric 1, experiment-re-experiment reliability
#column 2 is metric 2, average most likely transition probability
colnames(results) = c("Exp-Reexp","Avg MLT")
return(results)
}
#run function on transtion probs file
growth_reliability(tprobs, 3)