11 The telescoping sampler
Every chapter so far has taken \(K\) as given. Chapter 8 explained why choosing it is awkward — the classical theory breaks down at exactly the boundary where a \(K\)-state model sits inside a \((K+1)\)-state one — and left the question open. This chapter closes it by treating \(K\) as a parameter with a prior and sampling it along with everything else.
The obstacle has always been computational. Changing \(K\) changes the dimension of \(\theta\), which normally means reversible jump MCMC (Green 1995) and the design of dimension-matching proposals that actually get accepted. Richardson and Green (1997) did this for univariate normal mixtures, with split-and-merge moves constructed by hand; the effort involved is a fair measure of the difficulty. For HMMs, where a new state needs both an emission distribution and a full row and column of the transition matrix, it is worse. Allocation-based samplers (Nobile 2004; Nobile and Fearnside 2007) avoid some of that work by operating on the partition, and the telescoping sampler can be seen as completing that programme.
The telescoping sampler of Frühwirth-Schnatter et al. (2021) removes the obstacle. It samples \(K\) explicitly, needs no reversible jump machinery, no Jacobians and no bespoke proposals, and uses only the standard updates of a fixed-\(K\) model. The trick is to marginalise the component weights before updating \(K\), so that the move is a draw from a closed-form discrete distribution rather than a trans-dimensional jump.
By the end of this chapter you will be able to
- explain why changing \(K\) normally requires reversible jump, and how marginalising the weights avoids it;
- derive the closed-form conditional posterior of \(K\) given the partition, and identify where the dimension change is absorbed;
- distinguish static from dynamic mixtures of finite mixtures, and say what each implies about \(K\) versus \(K_+\);
- write out the telescoping sweep for an HMM and say which steps differ from an ordinary Gibbs sweep.
The chapter turns on a distinction that is easy to miss and does most of the work.
\(K\) is the number of components in the model — how many states the sampler currently carries.
\(K_+\) is the number of components occupied by the data — how many states are actually visited by the current state sequence.
Always \(K_+ \le K\). The components that make up the difference are empty: they have parameters, drawn from the prior, but no observations. Both have posteriors, and they answer different questions. \(p(K_+ \mid y)\) is usually the one of scientific interest.
11.1 Mixtures of finite mixtures
We develop the method for finite mixtures first, following Frühwirth-Schnatter et al. (2021) exactly, and then transfer it to HMMs in Section 11.4.
11.1.1 The model
Definition 11.1 (Mixture of finite mixtures) Observations \(y_1,\dots,y_N\) are generated by \[ K \sim p(K), \qquad \boldsymbol{\eta} \mid K \sim \mathrm{Dir}_K(\gamma_K, \dots, \gamma_K), \qquad \theta_k \mid K \stackrel{\text{iid}}{\sim} G_0 \ \ (k = 1,\dots,K), \] \[ \mathbb{P}(S_i = k \mid \boldsymbol{\eta}) = \eta_k, \qquad y_i \mid S_i = k, \, \theta \sim f(\cdot \mid \theta_k), \] independently across \(i = 1,\dots,N\).
The only unusual feature is that the Dirichlet parameter \(\gamma_K\) is allowed to depend on \(K\). This apparently minor freedom generates two qualitatively different model classes.
Definition 11.2 (Static and dynamic MFMs) The MFM is static if \(\gamma_K = \alpha\) does not depend on \(K\), and dynamic if \(\gamma_K = \alpha / K\).
The distinction controls what happens to the weights as \(K\) grows. Under the static specification the total Dirichlet mass \(K\gamma_K = K\alpha\) grows with \(K\), so extra components are pushed towards having comparable weight; small \(\alpha\) makes them sparse, and this is the sparse finite mixture of Malsiner-Walli et al. (2016). Under the dynamic specification the total mass \(K\gamma_K = \alpha\) is constant, so the prior on the weights stabilises as \(K\) increases — the behaviour of a Dirichlet process, and the reason the dynamic class connects to Bayesian nonparametrics.
The consequence is a different relationship between \(K\) and \(K_+\). Under a static MFM the prior on \(K_+\) tracks \(K\) fairly closely when \(\alpha\) is not small, so the number of components and the number of clusters are close. Under a dynamic MFM the two separate: \(K\) can be large while \(K_+\) stays small and grows only slowly with \(N\). Neither is right in the abstract; the choice is a modelling statement about whether “component” and “cluster” are meant to be the same thing.
11.1.2 The prior on the number of components
Common choices for \(p(K)\) on \(\{1,2,3,\dots\}\) are a shifted Poisson, \(K - 1 \sim \mathrm{Poi}(\lambda)\); a shifted geometric or negative binomial, which has heavier tails and is less informative; or a uniform prior on \(\{1,\dots,K_{\max}\}\), which is only defensible if \(K_{\max}\) is genuinely a hard bound. Miller and Harrison (2018) use the shifted Poisson and give the resulting properties of \(K_+\); the beta-negative-binomial prior favoured in Frühwirth-Schnatter et al. (2021) has heavier tails and is the safer default when little is known.
The posterior of \(K\) is more sensitive to \(p(K)\) than one would like, and considerably more sensitive than the posterior of \(K_+\). This is an argument for reporting \(p(K_+ \mid y)\) as the primary summary, and for checking sensitivity to \(p(K)\) explicitly. Grün et al. (2022) work through the sensitivity in detail on a benchmark data set.
11.2 The key computation
Everything rests on being able to sample \(K\) from its full conditional, given the current partition of the data but with the weights and the empty components’ parameters integrated out.
11.2.1 Marginalising the weights
Proposition 11.1 (Allocation likelihood with weights integrated out) Let \(\boldsymbol{\eta} \sim \mathrm{Dir}_K(\gamma_K,\dots,\gamma_K)\) and let \(S_1,\dots,S_N\) be conditionally i.i.d. with \(\mathbb{P}(S_i = k) = \eta_k\). Writing \(N_k = \#\{i : S_i = k\}\), \[ p(S_{1:N} \mid K) \;=\; \frac{\Gamma(K\gamma_K)}{\Gamma(N + K\gamma_K)} \prod_{k=1}^{K} \frac{\Gamma(N_k + \gamma_K)}{\Gamma(\gamma_K)} . \]
Proof. Integrate the multinomial likelihood against the Dirichlet density: \[ p(S_{1:N} \mid K) = \int \prod_{k=1}^K \eta_k^{N_k} \cdot \frac{\Gamma(K\gamma_K)}{\Gamma(\gamma_K)^K} \prod_{k=1}^{K} \eta_k^{\gamma_K - 1} \, d\boldsymbol{\eta} , \] and recognise the remaining integral as the normalising constant of a \(\mathrm{Dir}_K(\gamma_K + N_1, \dots, \gamma_K + N_K)\) density, namely \(\prod_k \Gamma(N_k+\gamma_K) / \Gamma(N + K\gamma_K)\).
The expression involves \(K\) only through the two Gamma functions and the product, and empty components contribute \(\Gamma(0 + \gamma_K)/\Gamma(\gamma_K) = 1\). So the product effectively runs over the \(K_+\) occupied components only, and the dependence on \(K\) is entirely explicit and cheap to evaluate.
11.2.2 From allocations to partitions
The allocation vector \(S_{1:N}\) carries labels; the partition \(\mathcal{C}\) — the grouping of observations, without labels — does not. Since the model is label-invariant, we work with the partition.
Proposition 11.2 (Partition likelihood) Let \(\mathcal{C}\) be a partition of \(\{1,\dots,N\}\) into \(K_+\) non-empty blocks of sizes \(N_1,\dots,N_{K_+}\). Then for \(K \ge K_+\), \[ p(\mathcal{C} \mid K) \;=\; \frac{K!}{(K-K_+)!} \cdot \frac{\Gamma(K\gamma_K)}{\Gamma(N + K\gamma_K)} \prod_{k=1}^{K_+} \frac{\Gamma(N_k + \gamma_K)}{\Gamma(\gamma_K)}, \] and \(p(\mathcal{C} \mid K) = 0\) for \(K < K_+\).
Proof. Each labelled allocation consistent with \(\mathcal{C}\) assigns the \(K_+\) blocks to \(K_+\) distinct labels out of \(K\), in some order. There are \(K(K-1)\cdots(K-K_++1) = K!/(K-K_+)!\) such assignments, and by Proposition 11.1 each has the same probability, since the expression there depends on the block sizes but not on which labels they carry. Summing gives the falling factorial times the common value.
The falling factorial is where the dimension change is absorbed. It counts the ways of embedding \(K_+\) occupied components into a model of size \(K\), and it is a combinatorial factor, not a Jacobian: no continuous change of variables is involved, which is exactly why reversible jump is unnecessary.
11.2.3 The conditional posterior of \(K\)
Theorem 11.1 (Full conditional for the number of components) Under Definition 11.1, with the weights \(\boldsymbol{\eta}\) and the parameters of empty components integrated out, \[ p(K \mid \mathcal{C}, \gamma_K) \;\propto\; p(K) \cdot \frac{K!}{(K-K_+)!} \cdot \frac{\Gamma(K\gamma_K)}{\Gamma(N + K\gamma_K)} \cdot \frac{1}{\Gamma(\gamma_K)^{K_+}} \prod_{k=1}^{K_+} \Gamma(N_k + \gamma_K) \] for \(K \ge K_+\), and zero otherwise.
Proof. Bayes’ theorem applied to Proposition 11.2, noting that the parameters \(\theta_k\) of empty components integrate out to one against their prior \(G_0\) — they enter no likelihood term — and that the parameters of occupied components have the same posterior contribution whatever \(K\) is, so they cancel from the ratio.
The distribution is supported on \(\{K_+, K_+ + 1, \dots\}\) and, with a prior having light enough tails, can be normalised by direct evaluation over a truncation range. Sampling from it is a single categorical draw.
Note what has not happened. No proposal was designed, no acceptance probability computed, no Jacobian evaluated. The move is a Gibbs step in a model where the awkward continuous parameters have been integrated away, leaving a discrete quantity whose conditional distribution is available in closed form. The word “telescoping” refers to the model expanding and contracting through this step: components are created empty and destroyed when empty, and the number of them simply follows the draw.
11.3 The algorithm
Definition 11.3 (Telescoping sampler for an MFM) Given the current \((K, \boldsymbol{\eta}, \theta_{1:K}, S_{1:N}, \alpha)\), one sweep is:
- Allocate. For each \(i\), draw \(S_i\) from \(\ \mathbb{P}(S_i = k) \propto \eta_k \, f(y_i \mid \theta_k)\), \(k = 1,\dots,K\).
- Count. Determine the occupied components, set \(K_+\) and the sizes \(N_1,\dots,N_{K_+}\), and relabel so that the occupied components come first.
- Update occupied parameters. For \(k = 1,\dots,K_+\), draw \(\theta_k\) from its posterior given \(\{y_i : S_i = k\}\) and \(G_0\).
- Update the hyperparameter. Draw \(\alpha\) by a Metropolis–Hastings step targeting \(p(\alpha \mid K, \mathcal{C}) \propto p(\alpha) \, p(\mathcal{C} \mid K, \gamma_K)\), using Proposition 11.2.
- Sample \(K\). Draw \(K\) from Theorem 11.1, on \(\{K_+, K_+ + 1, \dots, K_{\max}\}\).
- Refill. For \(k = K_+ + 1, \dots, K\), draw \(\theta_k \sim G_0\) — fresh empty components.
- Weights. Draw \[ \boldsymbol{\eta} \sim \mathrm{Dir}_K\bigl( \gamma_K + N_1, \dots, \gamma_K + N_{K_+}, \, \gamma_K, \dots, \gamma_K \bigr). \]
Steps 1, 3 and 7 are the ordinary updates of a fixed-\(K\) finite mixture, unchanged. Steps 2 and 6 are bookkeeping. The whole of the trans-dimensional problem is contained in steps 4 and 5, and neither requires anything beyond a Metropolis step on a scalar and a categorical draw.
Step 4 deserves a comment. Learning \(\alpha\) rather than fixing it matters, because \(\alpha\) controls how readily the model creates small components, and fixing it silently fixes the answer to the question being asked. Frühwirth-Schnatter et al. (2021) use a random-walk proposal on \(\log \alpha\) with a moderately diffuse prior.
11.3.1 Reading the output
The sampler produces draws of both \(K\) and \(K_+\), and they should be reported separately.
\(p(K_+ \mid y)\) estimates how many clusters the data support. It is the summary to lead with, it is comparatively robust to \(p(K)\), and it is permutation-invariant, so none of the relabelling machinery of Chapter 9 is needed to compute it.
\(p(K \mid y)\) estimates the number of components in the generating mechanism, including any that happen to be unoccupied in this sample. It is more prior-sensitive, and under a dynamic MFM it can be substantially larger than \(K_+\) without that indicating anything is wrong.
For state-specific summaries, the label switching problem returns in full force, and relabelling is complicated by the fact that the number of components changes between draws. The usual approach is to condition on the modal value of \(K_+\), discard the other draws, and relabel within that subsample.
11.4 Telescoping for HMMs
The construction transfers to hidden Markov models. What changes is that the single weight vector \(\boldsymbol{\eta}\) becomes the \(K\) rows of the transition matrix, and the independent allocations become a Markov chain — so the counts are transition counts and the allocation step becomes FFBS.
Frühwirth-Schnatter et al. (2021) develop and validate the telescoping sampler for finite mixtures, and the telescope package (Malsiner-Walli et al. 2026) implements that case. The HMM version below follows the same construction and is derived here rather than quoted from a published HMM implementation. It is presented as the natural extension, and anyone using it should validate it on simulated data before trusting it. Established alternatives for HMMs with unknown \(K\) include reversible jump (Robert et al. 2000), overfitted models with sparse priors (Rousseau and Mengersen 2011; Havre et al. 2015), and the nonparametric route of Chapter 10.
11.4.1 The model
Take the finite HMM of Definition 3.2 with a prior on \(K\): \[ K \sim p(K), \qquad \boldsymbol{\Gamma}_{j,\cdot} \mid K \stackrel{\text{iid}}{\sim} \mathrm{Dir}_K(\gamma_K,\dots,\gamma_K), \qquad \phi_k \mid K \stackrel{\text{iid}}{\sim} G_0 , \] with \(\gamma_K = \alpha\) (static) or \(\gamma_K = \alpha/K\) (dynamic), and a Dirichlet prior on \(\boldsymbol{\delta}\) or the stationarity constraint.
Given a state sequence \(s_{1:T}\), define the transition counts and row totals \[ N_{jk} = \#\{ t \ge 2 : s_{t-1} = j, \; s_t = k \}, \qquad N_{j\cdot} = \sum_{k=1}^{K} N_{jk} . \] Call state \(j\) occupied if it is visited by \(s_{1:T}\), and let \(K_+\) be the number of occupied states.
11.4.2 Marginalising the transition matrix
Proposition 11.3 (State-sequence likelihood with \(\boldsymbol{\Gamma}\) integrated out) With independent symmetric Dirichlet rows, \[ p(s_{1:T} \mid K, \boldsymbol{\delta}) \;=\; \delta_{s_1} \prod_{j=1}^{K} \left[ \frac{\Gamma(K\gamma_K)}{\Gamma(N_{j\cdot} + K\gamma_K)} \prod_{k=1}^{K} \frac{\Gamma(N_{jk} + \gamma_K)}{\Gamma(\gamma_K)} \right] . \]
Proof. The complete-data likelihood factorises across rows, \(\prod_{j,k} \gamma_{jk}^{N_{jk}}\), and the rows are a priori independent, so the integral separates into \(K\) independent Dirichlet–multinomial integrals, each evaluated as in Proposition 11.1.
The structure that made the mixture case work survives. A row \(j\) belonging to an unoccupied state has \(N_{j\cdot} = 0\) and \(N_{jk} = 0\) for all \(k\), so its bracket equals \(\Gamma(K\gamma_K)/\Gamma(K\gamma_K) = 1\) and it drops out. Within an occupied row, columns corresponding to unoccupied states also contribute \(1\). So the product runs effectively over the \(K_+ \times K_+\) block of occupied states, and the dependence on \(K\) is once again explicit and cheap.
11.4.3 The conditional posterior of \(K\)
Theorem 11.2 (Full conditional for \(K\) in an HMM) With \(\boldsymbol{\Gamma}\) and the emission parameters of unoccupied states integrated out, and states relabelled so the occupied ones come first, \[ p(K \mid s_{1:T}, \alpha) \;\propto\; p(K) \cdot \frac{K!}{(K - K_+)!} \cdot \prod_{j=1}^{K_+} \left[ \frac{\Gamma(K\gamma_K)}{\Gamma(N_{j\cdot} + K\gamma_K)} \prod_{k=1}^{K_+} \frac{\Gamma(N_{jk} + \gamma_K)}{\Gamma(\gamma_K)} \right] \] for \(K \ge K_+\), and zero otherwise.
Proof. Apply Bayes’ theorem to Proposition 11.3, using the argument of Proposition 11.2 to convert from labelled state sequences to their unlabelled equivalence class: the \(K_+\) occupied states can be assigned to \(K\) available labels in \(K!/(K-K_+)!\) ways, each giving the same value of Proposition 11.3. Rows of unoccupied states contribute a factor of \(1\), and the emission parameters of unoccupied states integrate to \(1\) against \(G_0\).
Compare with Theorem 11.1. The single Dirichlet–multinomial factor has become a product of \(K_+\) of them, one per occupied row, because a mixture has one weight vector and an HMM has one per state. Everything else — the falling factorial, the support on \(K \ge K_+\), the closed form — is unchanged.
11.4.4 The sampler
Definition 11.4 (Telescoping sampler for a finite HMM) Given \((K, \boldsymbol{\delta}, \boldsymbol{\Gamma}, \phi_{1:K}, s_{1:T}, \alpha)\), one sweep is:
- States. Draw \(s_{1:T}\) by forward filtering, backward sampling (Theorem 9.1) using the current \(K\)-state model.
- Count. Identify the occupied states, set \(K_+\), form the transition counts \(N_{jk}\), and relabel so occupied states come first.
- Emissions. For each occupied \(k\), draw \(\phi_k\) from its conjugate posterior given \(\{y_t : s_t = k\}\).
- Hyperparameter. Update \(\alpha\) by a Metropolis step targeting \(p(\alpha) \, p(s_{1:T} \mid K, \alpha)\) from Proposition 11.3.
- Number of states. Draw \(K\) from Theorem 11.2.
- Refill. For \(k = K_+ + 1, \dots, K\), draw \(\phi_k \sim G_0\).
- Transitions. Draw the rows independently, \[ \boldsymbol{\Gamma}_{j,\cdot} \sim \mathrm{Dir}_K\bigl( \gamma_K + N_{j1}, \dots, \gamma_K + N_{jK} \bigr), \qquad j = 1,\dots,K, \] where \(N_{jk} = 0\) whenever \(j\) or \(k\) is unoccupied.
- Initial distribution. Draw \(\boldsymbol{\delta} \sim \mathrm{Dir}_K(d + \mathbf{1}\{s_1 = \cdot\})\), or recompute the stationary distribution of the new \(\boldsymbol{\Gamma}\).
Steps 1, 3, 7 and 8 are exactly the sampler of Definition 9.2. The only additions are the two steps that move \(K\), and they cost a scalar Metropolis update and one categorical draw per sweep — negligible beside the \(\mathcal{O}(K^2T)\) FFBS pass.
Unoccupied states are visitable. Because unoccupied states have parameters drawn from \(G_0\) and rows drawn from the prior, FFBS in step 1 can assign observations to them. This is the mechanism by which new states are born, and it is why \(G_0\) must be a sensible prior over plausible emission parameters rather than something vague: a wildly diffuse \(G_0\) makes empty states unattractive and the sampler will never populate them.
Mixing is harder than for mixtures. Persistence in the transition matrix makes it costly to reassign a run of consecutive observations to a new state, so births of genuinely useful states can be slow. Running several chains from different \(K\) and comparing the posteriors of \(K_+\) is the minimum diagnostic.
Truncate honestly. Step 5 requires evaluating Theorem 11.2 over a finite range. The upper limit \(K_{\max}\) must be large enough that the tail mass beyond it is negligible; check this rather than assuming it.
Report \(K_+\). As in the mixture case, the posterior of the number of occupied states is the more stable and more interpretable summary.
Check your understanding. Your sampler reports a posterior for \(K\) concentrated on \(\{6,7,8\}\) but a posterior for \(K_+\) concentrated on \(3\). Is the model overfitting?
No — that pattern is expected, particularly under a dynamic MFM. \(K\) counts components carried by the model, including empty ones created afresh from \(G_0\) at every sweep; \(K_+\) counts states the data actually occupy. Empty components cost nothing and contribute nothing to the fit, which is exactly what makes the sampler’s dimension changes cheap.
Report \(p(K_+ \mid y)\) as the substantive answer: the data support three states. A large \(K\) with small \(K_+\) says the prior on \(K\) has heavy tails, not that the fit is bloated. If \(K\) and \(K_+\) were both concentrated on 8 while a three-state fit looked adequate, that would warrant concern.
11.5 Where this fits
Four ways of handling an unknown number of states, and what distinguishes them:
| Approach | Mechanism | Cost |
|---|---|---|
| Information criteria (Chapter 8) | fit each \(K\), compare scores | regularity conditions fail; heuristic only |
| Reversible jump (Robert et al. 2000) | trans-dimensional MH moves | proposals must be designed per model; tuning is delicate |
| Overfitted / sparse priors (Rousseau and Mengersen 2011) | fix \(K\) large, let the prior empty states | \(K\) is never sampled; asymptotic justification only |
| Telescoping (Frühwirth-Schnatter et al. 2021) | marginalise weights, sample \(K\) directly | needs conjugate weights; otherwise standard steps |
The overfitted-mixture approach is the closest competitor and deserves a fair comparison. Setting \(K\) deliberately too large with \(e_0 < 1\) causes superfluous states to empty out, and Rousseau and Mengersen (2011) show that the posterior concentrates on the correct number asymptotically. It is simpler to implement than telescoping, since it is just a fixed-\(K\) sampler. What it does not give is a posterior distribution over \(K\): the result is one number, sensitive to the chosen \(K\) and to \(e_0\), with no accompanying statement of uncertainty. Telescoping gives the distribution, at a small additional cost.
The requirement that limits telescoping is conjugacy of the weight prior. The marginalisation in Proposition 11.1 is what makes the \(K\) move a closed-form draw, and it depends on the Dirichlet–multinomial pair. A logistic-normal prior on the rows — natural when transition probabilities depend on covariates, as in Chapter 10 — has no such conjugacy, and the method as stated does not apply.
11.6 Summary
- Treating \(K\) as unknown is the honest response to the failure of classical model selection theory at the boundary between \(K\) and \(K+1\) states.
- The telescoping sampler marginalises the component weights before updating \(K\), turning a trans-dimensional move into a categorical draw from a closed-form distribution (Theorem 11.1). No reversible jump, no Jacobians.
- The falling factorial \(K!/(K-K_+)!\) absorbs the change in dimension combinatorially, by counting the ways of embedding \(K_+\) occupied components into \(K\) available slots.
- Static (\(\gamma_K = \alpha\)) and dynamic (\(\gamma_K = \alpha/K\)) specifications differ in how \(K\) and \(K_+\) relate; the choice is a modelling statement about whether components and clusters are the same thing.
- The construction transfers to HMMs by replacing the single weight vector with the rows of \(\boldsymbol{\Gamma}\), giving Theorem 11.2; the sampler is the standard HMM Gibbs sweep plus two cheap steps.
- Report \(p(K_+ \mid y)\), the posterior on occupied states, as the primary summary.
11.7 Further reading
Frühwirth-Schnatter et al. (2021) is the source for everything in this chapter and should be read directly; the exposition is careful about exactly the points that matter. Malsiner-Walli et al. (2016) introduced sparse finite mixtures and is the natural precursor. Miller and Harrison (2018) develop the MFM framework and its relationship to Dirichlet process mixtures, and Argiento and De Iorio (2022) examine that relationship from the nonparametric side. Grün et al. (2022) walk through a full applied analysis with sensitivity checks and is the best model for how to report results. For HMM-specific alternatives, Robert et al. (2000) for reversible jump and Rousseau and Mengersen (2011) for the asymptotics of overfitted models. Software for the mixture case is in Malsiner-Walli et al. (2026).