Autocorrelation for a categorical time series

Question

I am having trouble to compute the autocorrelation for different lags in a categorical time series. For instance, consider 3 possible categories: classA, classB and classC, and a vector x representing the series:

x <- c('classA', 'classB', 'classA', 'classC', 'classC', ...)

R contains the beautiful acf function, that computes the autocorrelation for different lags in a vector. Nevertheless the acf computes the Pearson correlation, which is (as far as I know) only suitable for continuous variables.

My first idea is to use the Cramer's V measure (available here) to compute the autocorrelation using something like cv.test(x, lag(x,1)), where the lag(x,1) function lags the x vector by 1, and then repeat the process for different lags. The problem is that I am not sure if this is the best approach, or even if this idea is correct.

Could you folks give-me some guidance about how to compute this autocorrelation?

Thanks

Answer 1

After reading the book (An Introduction to Discrete-Valued Time Series) suggested by @Ahmed Ali in the comments, I believe that the Cramer's V is the best answer to the problem.

For those who are interested in a R function to compute the autocorrelation using Cramer's V, here follow an adaptation of the R code that i found in the book.

data <- series#load your categorical series here. e.g. scan(file.choose(), what="character")
Tlen <- length(data)

states <- unique(data)
nostates <- length(states)
datanum <- match(data, states)

#Binarization
bincodes <- diag(1,nostates)
databin <- bincodes[datanum,]

#relative frequencies
hatpi <- colMeans(databin)

maxlag <- 100
hatbivprob <- array(0,c(nostates,nostates,maxlag))

for(k in c(1:maxlag)){#for each lag
    for(i in c(1:nostates)){#for each lagged vector representing a category
        for(j in c(1:nostates)){#for each vector representing a category
            hatbivprob[i,j,k] <- mean(databin[(k+1):Tlen,i]*databin[1:(Tlen-k),j])
        }
    }
}
#compare
indprob <- hatpi %*% t(hatpi)

#Cramer's v
cramer <- rep(0,maxlag)
for(k in c(1:maxlag)){
    cramer[k] <- sqrt(sum((hatbivprob[,,k]-indprob)^2/indprob)/(nostates-1))
}
plot(cramer, type="h", xlab = "k", ylab = "Cramer's   v(k)", lwd=4, ylim=c(-1,1))