1
$\begingroup$

I have a time series of an ordinal variable that I wish to model as a first-order Markov chain and estimate the matrix of transition probabilities. (I'm assuming the chain meets all the conditions to have a unique stationary distribution). As my ordinal variable has several categories and I don't have a vast amount of data, I'd like to make use of the ordering of the categories to reduce the number of model parameters by making some modelling assumption, analogous to an ordered logit or ordered probit model. I'd like to do this in R or python.

I've found a paper that proposes a model to do exactly what I want:
Varin, C. & Czado, C. (2010). "A mixed autoregressive probit model for ordinal longitudinal data." Biostatistics, 11, 127-138. doi: 10.1093/biostatistics/kxp042
Varin, C. & Vidoni, P. (2006). "Pairwise likelihood inference for ordinal categorical time series." Computational Statistics & Data Analysis, 51, 2365-2373. doi: 10.1016/j.csda.2006.09.009

The link in the paper to the authors' R code is dead however. It seems quite likely to me that this model or a similar one can be fitted in an R package published in the intervening ten years, maybe mvord (A flexible framework for fitting multivariate ordinal regression models with composite likelihood methods), or LMest (Latent Markov models for longitudinal categorical data). These packages are quite sophisticated and this isn't an area of statistics I'm familiar with, so I'm having trouble working out if or how I can do this. I'd be most grateful if anyone can show me a way, in one of these packages or somehow else.

The motivating example used in Varin & Czado's Vidoni's paper is a daily rainfall series for Alofi Island. This is available both as data(rain) in the markovchain R package and data(alofi) in the SMPracticals package. The three states are 1 (no rain), 2 (up to 5mm rain), 3 (over 5mm). The observed counts of one-step transitions are: \begin{pmatrix} 362 & 126 & 60 \\ 136 & 89 & 68 \\ 50 & 78 & 124 \end{pmatrix}

Modelling this as a first-order Markov chain ignoring the ordering of the states requires 6 parameters. More generally, a chain with $K$ states requires $K(K-1)$ parameters. Varin & Czado's Vidoni's model uses only $K$ parameters. I have up to 10 states, and would prefer a model with 10 parameters rather than 90 ! (I'll need to check it's a reasonable fit, of course).

$\endgroup$
3
  • $\begingroup$ The supplementary data with R code can still be found using a web search. Currently it gives 4 sources online: search term: "Supplementary material of the article A mixed autoregressive probit model for ordinal longitudinal data" $\endgroup$ Feb 13, 2020 at 16:29
  • 2
    $\begingroup$ I have several difficulties with reading this question. 1) What is your question? 2) Your example with rain transitions does not occur in the linked article from Varin and Czado, they had an example with reports about migraine in 6 levels/categories. 3) Is the question about simulating the $K(K-1)$ coefficients in the Markov chain (which Varin & Czado did not do), or is it about simulating the probabilities for the $K$ categories (which Varin & Czado did with a probit model with mixed effects and correlated/autoregressive errors)? $\endgroup$ Feb 13, 2020 at 16:36
  • $\begingroup$ @SextusEmpiricus Gah, awfully sorry, I somehow cited the wrong article! Now corrected in the question. I want to fit Varin & Vidoni's model in R (probably with mvord rather than LMest, having looked at both packages a bit more). I want fitted values for the K(K-1) transition probabilities from this model, which has K parameters. $\endgroup$
    – onestop
    Feb 15, 2020 at 20:07

1 Answer 1

2
+175
$\begingroup$

I believe one of the packages above can probably be used, but if you're willing to do a little work to do it easier (possible oxymoron) the following is probably the better way.

The lost supplementary material is here http://www.utstat.utoronto.ca/reid/sta4508s/VarinCzadoSupplement.pdf, and the original package description is here http://www.utstat.utoronto.ca/reid/sta4508s/maop/html/maop.html.

In the above two files you have code to replicate the Migraine experiment Varin and Czado carry out in the paper. However, the package is not available in a nice form anywhere, but I managed to bundle it together as a zip file, available here. To install, find your R package library and place the maop folder there. The only requirement is you use an R 2.x.x. installation rather than the current Version 3.x.x, as the package is rather old.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.