Let $(X_{n})_{n \in \mathbb{N}}$ be a discrete-time markov chain with state space $S := \{0,...,r-1\}$. Let $p_{i,j} = \mathbb{P}(X_{n} = j | X_{n-1} = i)$ be the transition probability and $(p_{0},...,p_{r-1})$ the initial distribution vector. Let $K:= n_{0} + ... + n_{r-1}$ and define $q(n_{0},...,n_{r-1})$ to be the probability that up to time $K$ state $0$ is visited exactly $n_{0}$ times, state $1$ is visited exactly $n_{1}$ times,..., state $r-1$ is visited exactly $n_{r-1}$.
My question is how to formally write down the probability $q(n_{0},...,n_{r-1})$
Define $V$ as the multiset containing all the visits, ie. $$V := \left\{\underbrace{0, \dots, 0}_{n_{0}}, \underbrace{1, \dots, 1}_{n_{1}}, \dots, \underbrace{r - 1, \dots, r - 1}_{n_{r-1}}\right\}$$ Then one crude way to express $q\left(n_{0}, \dots, n_{r - 1}\right)$ is by using the Law of Total Probability: $$ q\left(n_{0}, \dots, n_{r - 1}\right) = \sum_{\left\{i_{1}: i_{1} \in V\right\}}\sum_{\left\{i_{2}: i_{2} \in V\setminus\left\{i_{1}\right\}\right\}}\sum_{\left\{i_{3}: i_{3} \in V\setminus\left\{i_{1}, i_{2}\right\}\right\}}\dots\sum_{\left\{i_{K}: i_{K} \in V\setminus\left\{i_{1}, \dots, i_{K - 1}\right\}\right\}}p_{i_{1}}p_{i_{1}, i_{2}}\dots p_{i_{K - 1}, i_{K}}$$ Note that I have counted the initial state as a visit. There is possibly a more elegant/efficient method of calculating the probability which I do not know of. But it seems possible at the very least to write an algorithm which computes the probability this way for small enough $r$ and $K$.