We consider a bivariate, possibly non-homogeneous, finite-state Markov chain $(X,U)=\{(X_t,U_t)\}_{t=1}^n$. We are interested in the marginal process $X$, which typically is not a Markov chain. The goal is to find a realization (path) $x=(x_1,\ldots,x_n)$ with maximal probability $P(X=x)$. If $X$ is Markov chain, then such path can be efficiently found using the celebrated Viterbi algorithm. However, when $X$ is not Markovian, identifying the most probable path -- hereafter referred to as the Viterbi path -- becomes computationally expensive. In this paper, we explore the branch-and-bound method for finding Viterbi paths. The method is based on the lower and upper bounds on maximum probability $\max_x P(X=x)$, and the objective of the paper is to exploit the joint Markov property of $(X,Y)$ to calculate possibly good bounds in possibly cheap way. This research is motivated by decoding or segmentation problem in triplet Markov models. A triplet Markov model is trivariate homogeneous Markov process $(X,U,Y)$. In decoding, a realization of one marginal process $Y$ is observed (representing the data), while $X$ and $U$ are latent processes. The process $U$ serves as a nuisance variable, whereas $X$ is the process of primary interest. Decoding refers to estimating the hidden sequence $X$ based solely on the observation $Y$. Conditional on $Y$, the latent processes $(X, U)$ form a non-homogeneous Markov chain. In this context, the Viterbi path corresponds to the maximum a posteriori (MAP) estimate of $X$, making it a natural choice for signal reconstruction.