Summary: This paper is concerned with the analysis of a time series comprising the eruption inter-arrival times of the Old Faithful geyser in 2009. The series is much longer than other well-documented ones and thus gives a more comprehensive insight into the dynamics of the geyser. Basic hidden Markov models with gamma state-dependent distributions and several extensions are implemented. In order to better capture the stochastic dynamics exhibited by Old Faithful, the different non-standard models under consideration seek to increase the flexibility of the basic models in various ways: (i) by allowing non-geometric distributions for the times spent in the different states; (ii) by increasing the memory of the underlying Markov chain, with or without assuming additional structure implied by mixture transition distribution models; and (iii) by incorporating feedback from the observation process on the latent process. In each case it is shown how the likelihood can be formulated as a matrix product which can be conveniently maximized numerically.