Inference for Stochastically Contaminated Variable Length Markov Chains

TL;DR

This paper introduces a methodology to estimate parameters of stochastically contaminated variable length Markov chains, using a two-step process involving an adapted Baum-Welch algorithm and a bootstrap BIC, with proven consistency and validated through simulations.

Contribution

The paper develops a novel two-step estimation method for contaminated variable length Markov chains, including an adaptation of Baum-Welch and a bootstrap BIC for structure identification.

Findings

Method accurately recovers parameters in simulations

Estimates transition probabilities and noise parameters reliably

Effective for both additive and multiplicative noise regimes

Abstract

In this paper, we present a methodology to estimate the parameters of stochastically contaminated models under two contamination regimes. In both regimes, we assume that the original process is a variable length Markov chain that is contaminated by a random noise. In the first regime we consider that the random noise is added to the original source and in the second regime, the random noise is multiplied by the original source. Given a contaminated sample of these models, the original process is hidden. Then we propose a two steps estimator for the parameters of these models, that is, the probability transitions and the noise parameter, and prove its consistency. The first step is an adaptation of the Baum-Welch algorithm for Hidden Markov Models. This step provides an estimate of a complete order $k$ Markov chain, where $k$ is bigger than the order of the variable length Markov chain if it has finite order and is a constant depending on the sample size if the hidden process has infinite order. In the second estimation step, we propose a bootstrap Bayesian Information Criterion, given a sample of the Markov chain estimated in the first step, to obtain the variable length time dependence structure associated with the hidden process. We present a simulation study showing that our methodology is able to accurately recover the parameters of the models for a reasonable interval of random noises.

Figures

Tables

Table 1. Table 5.1: Transition probability matrix of the VLMC 𝑿 𝑿 \boldsymbol{X} .

$\omega$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$
010	0.05	0.95
110	0.87	0.13
00	0.27	0.73
1	0.38	0.62

Table 2. Table 5.2: Estimates of some noise parameter for TSCM and TPCM.

N=10.000			N=30.000
Noise	Estimate		Estimate
Real	TSCM	TPCM	TSCM	TPCM
0.01	0.019$\pm$ 0.011	0.020$\pm$ 0.013	0.015$\pm$ 0.008	0.017$\pm$ 0.009
0.05	0.055$\pm$ 0.012	0.058$\pm$ 0.013	0.046$\pm$ 0.008	0.054$\pm$ 0.009
0.25	0.256$\pm$ 0.013	0.253$\pm$ 0.012	0.245$\pm$ 0.007	0.246$\pm$ 0.008
0.45	0.457$\pm$ 0.012	0.443$\pm$ 0.011	0.454$\pm$ 0.008	0.455$\pm$ 0.009
0.55	0.557$\pm$ 0.011	0.544$\pm$ 0.012	0.553$\pm$ 0.006	0.556$\pm$ 0.007
0.75	0.742$\pm$ 0.013	0.758$\pm$ 0.014	0.753$\pm$ 0.006	0.746$\pm$ 0.007
0.95	0.954$\pm$ 0.012	-	0.947$\pm$ 0.007	-
0.99	0.986$\pm$ 0.011	-	0.992$\pm$ 0.006	-

Table 3. Table 5.3: Estimate of the transition matrix with TSCM regime for ϵ = 0.01 italic-ϵ 0.01 \epsilon=0.01

N=10.000			N=30.000
$\omega$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$
010	0.060$\pm$ 0.016	0.940$\pm$ 0.016	0.046$\pm$ 0.010	0.954$\pm$ 0.010
110	0.880$\pm$ 0.018	0.120$\pm$ 0.018	0.874$\pm$ 0.009	0.126$\pm$ 0.009
00	0.261$\pm$ 0.019	0.739$\pm$ 0.019	0.274$\pm$ 0.009	0.726$\pm$ 0.009
1	0.369$\pm$ 0.018	0.631$\pm$ 0.018	0.374$\pm$ 0.011	0.626$\pm$ 0.011

Table 4. Table 5.4: Estimate of the Transition Matrix with TSCM regime for ϵ = 0.05 italic-ϵ 0.05 \epsilon=0.05

N=10.000			N=30.000
$\omega$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$
010	0.076$\pm$ 0.020	0.924$\pm$ 0.020	0.068$\pm$ 0.015	0.932$\pm$ 0.015
110	0.885$\pm$ 0.021	0.115$\pm$ 0.021	0.862$\pm$ 0.014	0.138$\pm$ 0.014
00	0.279$\pm$ 0.022	0.731$\pm$ 0.022	0.275$\pm$ 0.013	0.725$\pm$ 0.013
1	0.350$\pm$ 0.021	0.650$\pm$ 0.021	0.362$\pm$ 0.013	0.638$\pm$ 0.013

Table 5. Table 5.5: Estimated transition matrix for TPCM with ϵ = 0.01 italic-ϵ 0.01 \epsilon=0.01 .

N=10.000			N=30.000
$\omega$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$
010	0.062$\pm$ 0.015	0.980$\pm$ 0.015	0.055$\pm$ 0.010	0.945$\pm$ 0.010
110	0.882$\pm$ 0.018	0.128$\pm$ 0.018	0.871$\pm$ 0.008	0.129$\pm$ 0.008
00	0.264$\pm$ 0.019	0.737$\pm$ 0.019	0.277$\pm$ 0.008	0.723$\pm$ 0.008
1	0.371$\pm$ 0.018	0.629$\pm$ 0.018	0.376$\pm$ 0.011	0.628$\pm$ 0.011

Table 6. Table 5.6: Transition probability matrix associated to 𝑿 𝑿 \boldsymbol{X}

$\omega$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$
0000	0.10	0.90
1000	0.50	0.50
100	0.83	0.17
10	0.25	0.75
1	0.25	0.75

Table 7. Table 5.7: Estimates of the transition probabilities with TSCM for ϵ = 0.01 italic-ϵ 0.01 \epsilon=0.01 .

N=10.000			N=30.000
$\omega$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$	$P\left(0|\omega \right)$	$P\left(1|\omega \right)$
0000	0.132$\pm$ 0.019	0.868$\pm$ 0.019	0.112$\pm$ 0.012	0.888$\pm$ 0.012
1000	0.532$\pm$ 0.018	0.468$\pm$ 0.018	0.515$\pm$ 0.011	0.485$\pm$ 0.011
100	0.838$\pm$ 0.015	0.162$\pm$ 0.015	0.825$\pm$ 0.009	0.175$\pm$ 0.009
10	0.258$\pm$ 0.016	0.742$\pm$ 0.016	0.246$\pm$ 0.011	0.754$\pm$ 0.011
1	0.243$\pm$ 0.018	0.757$\pm$ 0.018	0.253$\pm$ 0.011	0.747$\pm$ 0.011

Equations

$$p(a|\omega ):=\mathbb{P}(X_0=a|X_{-1}={\omega }_{-1},X_{-2}={\omega }_{-2},...)$$

$$p(\omega ):=\mathbb{P}(X_{-j}^{-1}=\omega )$$

$$\mathbb{P}\left(X_0=a|X_{-\infty }^{-1}=x_{-\infty }^{-1}\right)=p(a|\omega ),$$

$$\bullet \text{ignorespaces}\text{ignorespaces}\text{ignorespaces}\text{ignorespaces}\text{ignorespaces}\text{ignorespaces}\text{ignorespaces}\text{ignorespaces}$$

$$Z_t=X_t\oplus {\xi }_t,$$

$$\mathbf{A}_{\boldsymbol{S}}=\left\{p(a|\omega)\right\}=\mathbb{P}\left(X_{0}=a\bigg{|}X_{-l(\omega)}^{-1}=\omega\right),\forall\ a\in E,\forall\ \omega\in\cal{T}$$

$$\mathbf{B}_{\boldsymbol{S}}=\left\{b_{\omega}(z)\right\}=\mathbb{P}\left(Z_{t}=z\bigg{|}X_{t-l(\omega)+1}^{t}=\omega\right),\forall\ \omega\in\ E^{l(\omega)},\forall z\in E$$

$${\mathbf{\pi }}_{\mathbf{S}}=\left\{{\pi }_{\omega }\right\}=\mathbb{P}\left(X_{-j}^{-1}=\omega \right),\forall \omega \in \mathcal{T}$$

$$\mathbb{L}({\mathbf{\lambda }}_{\mathbf{S}}|Z)=\mathbb{P}({\mathbf{\lambda }}_{\mathbf{S}}|Z).$$

$$\mathbb{P}\left(Z_t=z_t|X_{-l(\omega )+1}^t=\omega \right)=\mathbb{P}\left(Z_t=z_t|X_t=a_t\right)=\mathbb{P}\left({\xi }_t=b_t\right)I_{\left\{z_t=a_t\oplus b_t\right\}},$$

$$\sum _{a_t,b_t\in E}\prod _{t=1}^T\left[\mathbb{P}({\xi }_t=b_t)\right]\left[\mathbb{P}\left(X_1^k=a_1^k\right)\prod _{t=k+1}^T\mathbb{P}\left(X_t=a_t|X_{t-l(\omega )}^{t-1}=a_{t-l(\omega )}^{t-1}\right)\right]\prod _{t=k+1}^T{I}_{\left\{z_t=a_t\oplus b_t\right\}}.$$

$$Z_t=X_t\cdot {\xi }_t.$$

$$\mathbf{A}_{\boldsymbol{P}}=\left\{p(a|\omega)\right\}=\mathbb{P}\left(X_{t}=a\bigg{|}X_{-l(\omega)}^{-1}=\omega\right),\forall\ a\in E,\forall\ \omega\in\cal{T}$$

$$\mathbf{B}_{\boldsymbol{P}}=\left\{b_{\omega}(z)\right\}=\mathbb{P}\left(Z_{t}=z\bigg{|}X_{t-l(\omega)+1}^{t}=\omega\right),\forall\ \omega\in\ E^{l(\omega)},\forall z\in\ E$$

$${\mathbf{\pi }}_{\mathbf{P}}=\left\{{\pi }_{\omega }\right\}=\mathbb{P}\left(X_{-j}^{-1}=\omega \right),\forall \omega \in \mathcal{T}$$

$$\mathbb{P}\left(Z_t=z_t|X_{-l(\omega )+1}^t=\omega \right)=\mathbb{P}\left(Z_t=z_t|X_t=a_t\right)=\mathbb{P}\left({\xi }_t=b_t\right)I_{\left\{z_t=a_t.b_t\right\}}.$$

$$\sum _{a_t,b_t\in E}\prod _{t=1}^T\left[\mathbb{P}({\xi }_t=b_t)\right]\left[\mathbb{P}\left(X_1^k=a_1^k\right)\prod _{t=k+1}^T\mathbb{P}\left(X_t=a_t|X_{t-l(\omega )}^{t-1}=a_{t-l(\omega )}^{t-1}\right)\right]\prod _{t=1}^T{I}_{\left\{z_t=a_t.b_t\right\}}.$$

$$X_r^{\ast }:=X_{(r+k)-1}^r,r=1,...,(T-k)+1,T\in \mathbb{N}.$$

$$Z_r^{\ast }=Z_{(r+k)-1}^r,r=1,...,(T-k)+1.$$

$$X^{\ast }=x_1^{r\ast }=\left\{00,01,10,01,11,10,01,...,01\right\}$$

$$Z^{\ast }=z_1^{r\ast }=\left\{00,01,11,10,00,01,...,10\right\}$$

$${\beta }_r\left(\omega \right)=\sum _{\nu \in E^k}{p}^{\ast }\left(\omega |\nu \right)b_{\omega }\left(z_{r+1}^{\ast }\right){\beta }_{r+1}\left(\omega \right),\forall \omega \in E^k,r=(T-k)+2,(T-k)+3,...,1.$$

$${\delta }_r\left(\omega ,\nu \right)=\mathbb{P}\left(X_r^{\ast }=\omega ,X_{r+1}^{\ast }=\nu |Z^{\ast },{\mathbf{\lambda }}^{\ast }\right).$$

$${\gamma }_r\left(\omega \right)=\frac{{\alpha }_r\left(\omega \right){\beta }_r\left(\omega \right)}{\sum _{\omega \in E^k}{\alpha }_r\left(\omega \right){\beta }_r\left(\omega \right)},$$

$${\delta }_r\left(\omega ,\nu \right)=\frac{{\alpha }_r\left(\omega \right)p^{\ast }\left(\omega |\nu \right)b_{\omega }\left(z_{r+1}^{\ast }\right){\beta }_{r+1}\left(\omega \right)}{\sum _{\omega \in E^k}\sum _{\nu \in E^k}{\alpha }_r\left(\omega \right)p^{\ast }\left(\omega |\nu \right)b_{\omega }\left(z_{r+1}^{\ast }\right){\beta }_{r+1}\left(\omega \right)},\forall \omega ,\nu \in E^k.$$

$${\hat{\mathbf{\lambda }}}^{\ast }=\arg \underset{{\tilde{\mathbf{\lambda }}}^{\ast }\in \mathbf{\Lambda }}{\max }\mathbb{L}\left({\tilde{\mathbf{\lambda }}}^{\ast }|Z^{\ast }\right).$$

$$\widehat{ϵ}=\arg \underset{{ϵ}_j}{\max }\mathbb{L}\left({\tilde{\mathbf{\lambda }}}^{\ast }({ϵ}_j)|Z^{\ast }\right).$$

$${\hat{N}}_m(\omega ,a)=\left|\left\{i:D(m)<i\le m,{\hat{x}}_{i-l(\omega )}^{i-1}=\omega ,{\hat{x}}_i=a\right\}\right|$$

$${\hat{N}}_m(\omega )=\left|\left\{i:D(m)<i\le m,{\hat{x}}_{i-l(\omega )}^{i-1}=\omega \right\}\right|.$$

$$BI{C}_{\mathcal{T}}({\hat{x}}_1^m)=-\log M{L}_{\mathcal{T}}({\hat{x}}_1^m)+\frac{(|E|-1)|\mathcal{T}|}{2}\log m,$$

Full text

Keywords: Contaminated Process; Variable Length Markov Chain; Bayesian Information Criterion; Bootstrap; EM algorithm

1 Introduction

In this paper, we present a new methodology to estimate the parameters of some stochastically contaminated processes. We assume that a hidden, original, process is contaminated by some noise and only the contaminated process is observable. In the case where the original hidden process is a Markov chain, this model is known in the literature as Hidden Markov Model (HMM) introduced in 1966 by [1]. This model has a large amount of work devoted to it due its importance and applications in subjects such as bioinformatics, communications engineering, finance and many others. A comprehensive treatment of inference for hidden Markov models can be found in [3].

We analyze this problem considering that the hidden original source belongs to a larger class of process where the order of dependence in the past is not fixed as in a Markov chain. These models are known in the literature as Variable Length Hidden Markov Models (VLHMM). As far as we know they appeared in the first time in a paper about the analysis of human movement [15], [16]. In [16], the author analyzes the 3D motion by rotating 19 major joints of the human body. The authors claim that VLHMM is superior to other models studied in relation to their efficiency and accuracy in modeling multivariate time series.

There are some previous works that analyze these class of model from a theoretical point of view, which we consider as a starting point. In [4] the authors assume that the original source is a chain with an infinite order with a binary alphabet which is contaminated by adding to each symbol a random Bernoulli noise, independent of the original source. In [9] the authors also consider a model where each symbol is multiplied by a Bernoulli random noise. In both articles, the authors showed that the difference between the transition probabilities of the contaminated process and the original process is limited by a constant $c$, where $c$ is a linear nondecreasing function of the random Bernoulli noise (more details in [4, 9]). Henceforth if the random Bernoulli noise is small enough then the contaminated sample can be used to estimate the transition probability matrix of the original hidden process. However, if the random noise is not small enough, the approximation of the hidden transition probabilities by the estimated transition probabilities of the contaminated process is not satisfactory. Then it is crucial to estimate this noise parameter in order to know if this approximation can be applied or not. But they do not address the problem of model parameters estimation. This estimation is the main goal of the present work.

In [8] it is presented an important result in parameter estimation for a class of contaminated models similar to those analyzed in this paper. The class of models discussed in [8] is, on one hand larger than the one discussed here, since it allows random noise with a greater variety of distributions, but it is more restrictive on the other hand since only the last symbol seen in the past is considered in the conditional distributions. The author proposes an estimator based on a penalized likelihood function but, according to the author, the penalty proposed in his paper is worse than the Bayesian Information Criterion penalty, which is used in our methodology.

In this paper, we present consistent estimators for hidden parameters of the contaminated models. The simplicity of the models considered here allows us to propose an inferential methodology based on an EM algorithm for Hidden Markov Models and in the Bayesian Information Criterion (BIC).

Besides, we present a sensitivity study of the estimators when we let the random noise to increase. Our goal with this study is to know for which interval of contamination noise the estimation procedure provides accurate estimates.

This paper is organized as follows: Section 2 presents the basic notations and some preliminary definitions. Section 3 presents the models discussed in this paper and some theoretical results. Section 4 presents an inferential methodology and the main results for model parameter estimation. In Section 5 we perform a sensitivity simulation study concerning the influence of the random noise in the estimation procedure. Section 6 presents some conclusions. Finally, in Appendix, we present the proofs of the results presented in Sections 3 and 4.

2 Basic Notation and Definitions

Let us consider a finite discrete alphabet $E=\left\{0,1,...,N-1\right\}$, with cardinality $|E|=N$. Given two integers $m,n\in\mathbb{Z}$, with $m\leq n$, we shall use the short notation $\omega_{m}^{n}$ to denote the string $(\omega_{m},...,\omega_{n})$ of symbols in $E$, and let $E^{l(\omega_{m}^{n})}$ denote the set containing such strings, where $l(\omega_{m}^{n})=n-m+1$ is the length of the string $\omega_{m}^{n}$. An empty string is denoted by $\emptyset$ and $l(\emptyset)=0$.

Given two strings $\omega$ and $\upsilon$, such that $l(\omega)<\infty$, we denote by $\upsilon\omega$ the string with length $l(\upsilon)+l(\omega)$ obtained by concatenation of this two strings. The concatenation can be extended to the case when the strings are semi-infinite $\upsilon=...\omega_{-2}\omega_{-1}$.

We shall say that a string $\nu$ is a suffix of the string $\omega$ if there exists a substring $\eta$ such that $\omega=\eta\nu$. If $1\leq l(\eta)<l(\omega)$, $\nu$ is a proper suffix of $\omega$ and we write $\nu\prec\omega$. When $\nu=\omega$ we denote $\nu\preceq\omega$.

In this work we consider $\boldsymbol{X}=\left\{X_{t}\right\}_{t\in\mathbb{Z}}$ as an ergodic stochastic process on the discrete alphabet $E$. Given an infinite string $\omega_{-\infty}^{-1}$ and $a\in E$, we denote by

[TABLE]

the transition probabilities of the process $\boldsymbol{X}$ and for a finite string $\omega\in E^{j}$, we denote by

[TABLE]

the initial probability distribution.

Definition 2.1**.**

A finite string $\omega\in\cup_{j=1}^{\infty}E^{j}$ is a context for $\boldsymbol{X}$ if it satisfies:

(i) For every semi-infinite string $x_{-\infty}^{-1}$ with $\omega$ as a suffix,

[TABLE]

for every $a\in E$.

(ii) No proper suffix of $\omega$ satisfies (1).

An infinite context is a semi-infinite string $\omega_{-\infty}^{-1}$ such that no suffix $\omega_{-j}^{-1},\;j\in\mathbb{N}$, is a context.

Definition 2.2**.**

A set $\cal{T}$ of contexts is called Context Tree, associated to the process $\boldsymbol{X}$, if no $\omega_{1}\in\cal{T}$ is a proper suffix any other $\omega_{2}\in\cal{T}$. A context tree satisfying condition (ii): is called irreducible.

Each context $\omega\in\cal{T}$ can be viewed as a path from a leave to a root (see Figure 2.1). The branches of the tree $\cal{T}$ are identified with the context (finite or infinite) $\omega\in\cal{T}$ in the past, the root represents the present time and it is represented by the empty context $\emptyset$.

Figure 2.1 shows an order 3 context tree $\cal{T}$ taking values in $E=\left\{0,1\right\}$.

Definition 2.3**.**

A tree $\cal{T}$ is complete if each node has $|E|$ branches.

Definition 2.4**.**

If $\cal{T}$ is a tree such that $l(\omega)=L,\forall\omega\in\cal{T}$ then it is a order $L$ Markov chain denoted by $L_{full}$ .

We denote d(\cal{T}):=\max$$\left\{l\left(\omega\right):\omega\in\cal{T}\right\} as the depth, or order, of the tree.

A stationary stochastic process $\boldsymbol{X}$ in $E$ is a Variable Length Markov Chain (VLMC) compatible with a pair $(\cal{T}$, $p(a|\omega))$ if it satisfies Definition 2.1.

Definition 2.5**.**

Given a positive integer $k$, define the truncated tree \cal{T}\big{|_{k}} of order $k$ as

\cal{T}\big{|_{k}}:=\left\{\omega\in\cal{T}\right.: $\left.l(\omega)\leq k\right\}\bigcup\left\{\omega:l(\omega)=k\right.$ and $\omega\prec\upsilon$, for some $\left.\upsilon\in\cal{T}\right\}$.

Definition 2.6**.**

A Variable Length Hidden Markov Model (VLHMM) is a bivariate stochastic process $(\boldsymbol{X},\boldsymbol{Z})$ characterized by:

(i) $\boldsymbol{X}$, the hidden VLMC, with tree $\cal{T}$, taking values in the alphabet $E$;

(ii) $\boldsymbol{Z}$, the observable process assuming values in a set $O$;

(iii) $\bf{A}$, the transition probability matrix of the hidden process $\boldsymbol{X}$ given by $p(a|\omega)$, $\forall a\in E,\forall\omega\in\cal{T}$;

(iv) $\bf{B}$ (emission distribution), the conditional probability distribution for a symbol in the observable process, given the context $\omega$ in the hidden process defined by $\mathbb{P}(Z_{t}=z|X_{(t-l(\omega))+1}^{t}=\omega)$, $\forall\omega\in\cal{T}$, $\forall z\in O$;

(v) $\boldsymbol{\pi}$, the initial distribution of the hidden process defined by $\mathbb{P}(X^{l(\omega)}_{1}=\omega),\ \forall\omega\in\cal{T}$.

Remark 2.1**.**

If the hidden process $\boldsymbol{X}$ is markovian and $\mathbb{P}(Z_{t}=k|X_{(t-l(\omega))+1}^{t}=\omega)=\mathbb{P}(Z_{t}=k|X_{t}=j)$, that is, if the emission distribution loses memory of the whole context then this process is an HMM, a particular case of a VLHMM.

3 Stochastic Contaminated Models

Let us consider $\boldsymbol{X}$ a VLMC as in definition 2.1 taking values in $E=\left\{0,1,\ldots,N-1\right\},N\in\mathbb{N}$ and let $\boldsymbol{\xi}=\left\{\xi_{t}\right\}_{t\in\mathbb{Z}}$ be a sequence of independent random variables with ${\mathbb{P}}(\xi_{t}=i)=\epsilon_{i}$ and such that $\displaystyle\sum_{i=0}^{N-1}\epsilon_{i}=1$, independently of $\boldsymbol{X}$. Closely following the models presented in [4] and [9], we propose two stochastic contaminated models as follows.

3.1 Type Sum Contamination Model

A Type Sum Contaminated Model (TSCM) is a bivariate process $(\boldsymbol{Z},\boldsymbol{X})$ where $\boldsymbol{X}$ is a VLMC, with associated tree $\cal{T}$ and the contaminated process $\boldsymbol{Z}$ is defined by

[TABLE]

where for $a,b\in E$ we define $a\oplus b=a+b\;mod(|E|)$. The vector parameter of this model is $\boldsymbol{\lambda}_{\boldsymbol{S}}=(\mathbf{A}_{\boldsymbol{S}},\mathbf{B}_{\boldsymbol{S}},\boldsymbol{\pi}_{\boldsymbol{S}})$, where

[TABLE]

is the transition probability matrix of the hidden process $\boldsymbol{X}$,

[TABLE]

is the probability distribution of the observed symbol given the hidden string of the original process (emission distribution),

[TABLE]

is the initial distribution of the original process $\boldsymbol{X}$.

Let $\boldsymbol{\lambda_{S}}$ be the set of parameters of the bivariate process $(\boldsymbol{X},\boldsymbol{Z})$. Given an observable sample, $Z=z_{1}^{T}$, $T\in\mathbb{N}$, the likelihood function $\mathbb{L}(\boldsymbol{\lambda_{S}}|Z)$ is defined by

[TABLE]

Remark 3.1.1**.**

We observe that in the model proposed in [4] $\boldsymbol{X}$ is a chain with infinite order and $E=\{0,1\}$.

Remark 3.1.2**.**

We observe that the TSCM is a VLHMM. In consequence of the model, the emission distribution depends only on the last symbol of the context instead of the whole context. This fact allows us to propose some adaptations in the Expectation-Maximization algorithm [7], originally for HMM and known in the literature as the Baum-Welch algorithm, to estimate $\boldsymbol{\lambda_{S}}$, the parameters of the VLHMM.

The following proposition shows that: item $\bf{(i)}$ the emission distribution considering a TSCM loses memory of the past symbols of the context; Item $\bf{(ii)}$ the number of computations needed to calculate the likelihood for a sample of a TSCM is of order $|E|^{T}$, where $T$ is the sample size, and hence a direct computation is not viable.

Proposition 3.1.1**.**

Let $\boldsymbol{Z}$ be a contaminated process as in TSCM.

i) For every $z_{t},a_{t},b_{t}\in E$ and every $\omega\in\cal{T}$, where $\cal{T}$ is the context tree of the process $\boldsymbol{X}$, the Emission distribution is

[TABLE]

where

I_{\left\{z_{t}=a_{t}\oplus b_{t}\right\}}=\left\{\begin{array}[]{rl}1&\mbox{ if z_{t}=a_{t}\oplus b_{t}}\\ 0&\mbox{ otherwise.}\end{array}\right.**

**ii) Let us consider a sample $Z=z_{1}^{T}$, $T\in\mathbb{N}$ of the contaminated process $\boldsymbol{Z}$, such that $l(\omega)\leq T,\forall\omega\in\cal{T}$, and $k=\max\left\{l(\omega):\omega\in\cal{T}\right\}$. Then the likelihood function $\mathbb{L}(\boldsymbol{\lambda_{S}}|Z)$ for the contaminated process $\boldsymbol{Z}$ can be written as: **

[TABLE]

Proof See Appendix.

Remark 3.1.3**.**

If the VLMC $\boldsymbol{X}$ has infinite order, we can obtain, in an analogous way, a truncated version of the likelihood function $\mathbb{L}(\boldsymbol{\lambda_{S}}|Z)$ at some finite order $L<T$.

Remark 3.1.4**.**

Due to item $\bf{i)}$ the emission distribution $\mathbf{B}_{\boldsymbol{S}}=\mathbf{B}_{\boldsymbol{S}}(\epsilon)$ depends only on the random noise.

3.2 Type Product Contaminated Model

We define a Type Product Contaminated Model (TPCM) also as a bivariate process ($\boldsymbol{Z},\boldsymbol{X}$) where $\boldsymbol{X}$ is a VLMC and the contaminated process $\boldsymbol{Z}$ is obtained by

[TABLE]

We denote the vector of parameters of this model by $\boldsymbol{\lambda}_{\boldsymbol{P}}=(\mathbf{A}_{\boldsymbol{P}},\mathbf{B}_{\boldsymbol{P}},\boldsymbol{\pi}_{\boldsymbol{P}})$, where

[TABLE]

is the transition probability matrix of the hidden process $\boldsymbol{X}$,

[TABLE]

is the probability distribution of the observed symbol given the hidden string of the original process (emission distribution),

[TABLE]

is the initial distribution of the original process $\boldsymbol{X}$.

Remark 3.2.1**.**

The model proposed in [9] is a TPCM with $E=\{0,1\}$.

Remark 3.2.2**.**

Observe that the TPCM is also a VLHMM and, also as consequence of the model, the emission distribution depends only on the last symbols of the context and not on the whole context.

In the same way as in TSCM, the following proposition shows that: Item $\bf{(i)}$, the emission distribution considering a TPCM also loses memory of the past symbols of the context; Item $\bf{(ii)}$, the number of computations needed to calculate the likelihood for a sample of a TSCM is also of order $|E|^{T}$, where $T$ is the sample size, and hence a direct computation is not viable.

Proposition 3.2.1**.**

*Let $\boldsymbol{Z}$ be a contaminated process as in a TPCM.

i) For every $z_{t},a_{t},b_{t}\in E$ and every $\omega\in\cal{T}$, where $\cal{T}$ is the context tree of the process $\boldsymbol{X}$, the emission distribution is*

[TABLE]

**ii) Let us consider a sample $Z=z_{1}^{T}$, $T\in\mathbb{N}$ of the contaminated process $\boldsymbol{Z}$, such that $l(\omega)\leq T,\forall\omega\in\cal{T}$, and $k=\max\left\{l(\omega):\omega\in\cal{T}\right\}$. Then the likelihood function $\mathbb{L}(\boldsymbol{\lambda_{P}}|Z)$ for the contaminated process $\boldsymbol{Z}$ can be written as: **

[TABLE]

Proof See Appendix.

Remark 3.2.3**.**

If the VLMC $\boldsymbol{X}$ has infinite order, we can obtain, in an analogous way, a truncated version of the likelihood function $\mathbb{L}(\boldsymbol{\lambda_{P}}|Z)$ at some finite order $L<T$.

Remark 3.2.4**.**

Due to item $\bf{i)}$, $\mathbf{B}_{\boldsymbol{P}}=\mathbf{B}_{\boldsymbol{P}}(\epsilon)$, the emission distribution, depends only on the random noise.

4 Inference for TSCM and TPCM

Considering that a direct computation of the likelihood is intractable, we propose an EM algorithm, based on the Baum-Welch algorithm for HMM, to iteratively compute the likelihood for VLHMM. To this end, we propose a new parameterization of the models as follows.

Let us first consider a finite VLHMM $(\boldsymbol{X},\boldsymbol{Z})$, where $\boldsymbol{X}$ has a finite tree $\cal{T}$ and let $k$ be the length of the biggest context, $k=\max\left\{l(\omega):\omega\in\cal{T}\right\}$. We can rewrite the order $k$ VLHMM as an order $k$ Markov chain as $\boldsymbol{X^{*}}=\left\{X^{*}_{r}\right\}_{r\in\mathbb{N}}$, which is a $k_{full}$ tree, assuming values in $E^{k}$ where

[TABLE]

The transitions probabilities of $\boldsymbol{X^{*}}$ are defined by $\boldsymbol{A^{*}}=\left\{p^{*}(\omega|\nu)\right\},\forall\omega,\nu\in\ E^{k}$ and with initial distribution $\boldsymbol{\pi^{*}}=\left\{\mathbb{P}(X^{*}_{1}=\omega)\right\},\forall\omega\in\ E^{k}$.

Similarly we define a new observable process $\boldsymbol{Z^{*}}=\left\{Z^{*}_{r}\right\}_{r\in\mathbb{N}}$, assuming values in $E^{k}$, where

[TABLE]

In this way, the VLHMM $(\boldsymbol{X},\boldsymbol{Z})$ can be viewed as an HMM $(\boldsymbol{X^{*}},\boldsymbol{Z^{*}})$ with set of parameters $\boldsymbol{\lambda^{*}}=(\boldsymbol{A^{*}},\boldsymbol{B^{*}},\boldsymbol{\pi^{*}})$.

As an example, suppose that $\boldsymbol{X}$ is a VLMC that assumes values in the alphabet $E=\left\{0,1\right\}$ and $k=\max\left\{l(\omega):\omega\in\cal{T}\right\}=2$. Given a sample of the hidden process $\boldsymbol{X}$, $x_{1}^{T}=\left\{0,0,1,0,1,1,0,1,...,0,1\right\}$, then the sample associated to Markov chain $\boldsymbol{X^{*}}$ with order $k=2$, denoted by $X^{*}$, is

[TABLE]

And for an observable sample $z_{1}^{T}=\left\{0,0,1,1,0,0,1,...,1,0\right\}$, of the contaminated process $\boldsymbol{Z}$, we have that the new observable sample of the process $\boldsymbol{Z^{*}}$ is

[TABLE]

Now we are ready to apply the Baum-Welch EM algorithm [1] to the process $(\boldsymbol{X}^{*},\boldsymbol{Z}^{*})$. Given a sample $Z^{*}$, the forward variable is defined as

\alpha_{r}\displaystyle\left(\omega\right)=\mathbb{P}\left(z^{*}_{1},...,z^{*}_{r},X^{*}_{r}=\omega\bigg{|}\boldsymbol{\lambda}^{*}\right).

By induction, we have that

$\alpha_{1}\left(\omega\right)=\pi_{\omega}b_{\omega}\left(z^{*}_{1}\right),\forall\omega\in E^{k},$

$\alpha_{r+1}\left(\omega\right)=\left[\displaystyle\sum_{\nu\in E^{k}}\alpha_{r}\left(\omega\right)p^{*}\left(\omega|\nu\right)\right]b_{\omega}\left(z^{*}_{r+1}\right),\forall\omega\in E^{k},\ \ 2\leq r\leq(T-k)+1$.

Similarly, the backward variable is defined as

$\beta_{r}\left(\omega\right)=\mathbb{P}\left(z^{*}_{r+1},z^{*}_{r+2},...,z_{(T-k)+1}|X^{*}_{r}=\omega,\boldsymbol{\lambda}^{*}\right),$

and by induction follows

$\beta_{(T-k)+1}\left(\omega\right)=1,\forall\omega\in E^{k},$

[TABLE]

Now we define

$\gamma_{r}\left(\omega\right)=\mathbb{P}\left(X^{*}_{r}=\omega|Z^{*},\boldsymbol{\lambda}^{*}\right),$

and

[TABLE]

Given $\alpha_{r}\displaystyle\left(\omega\right)$ and $\beta_{r}\displaystyle\left(\omega\right)$, we can write

[TABLE]

and

[TABLE]

Henceforth, the parameter vector $\boldsymbol{\lambda}^{*}$ can be updated in the follwoing way:

$\hat{\boldsymbol{\pi}}^{*}=\left\{\hat{\pi}_{\omega}^{*}\right\}_{\forall{\omega\in E^{k}}}=\left\{\gamma_{1}(\omega)\right\},$

$\hat{\mathbf{A}}^{*}=\left\{\hat{p}^{*}(\omega|\nu)\right\}_{\forall{\omega,\nu\in E^{k}}},$ where $\hat{p}^{*}(\omega|\nu)=\frac{\displaystyle\sum^{T+k-1}_{r=1}\delta_{r}(\omega,\nu)}{\displaystyle\sum^{T+k-1}_{r=1}\gamma_{r}(\omega)}$,

$\hat{\mathbf{B}}^{*}=\left\{\hat{b}_{\omega}(\nu)\right\}_{\forall{\omega,\nu\in E^{k}}}$, where $\hat{b}_{\omega}(\nu)=\frac{\displaystyle\sum^{T+k}_{r=1}I_{\left\{z^{*}_{r}=\nu\right\}}\gamma_{r}(\omega)}{\displaystyle\sum^{T+k}_{r=1}\gamma_{r}(\omega)}.$

Remark 4.1**.**

If $\boldsymbol{X}$ has infinite order, for a finite sample $z_{1}^{T}$ it is possible to estimate only a truncated tree \cal{T}\big{|_{k}}, where $k\in\mathbb{N}$ is as big as possible. We apply the same methodology proposed for finite trees to obtain \hat{\cal{T}}\big{|_{k}} using a truncated likelihood function, given the sample of the observable process $\boldsymbol{Z}$.

4.1 Estimation Methodology

First of all, we stress that, for estimation purposes, the main difference between a VLHMM and an HMM is that in an HMM the states of the original (hidden) process are known while in a VLHMM they are unknown. Then the estimation procedure has also to learn the contexts of the hidden tree $\cal{T}$ associated to the VLMC $\boldsymbol{X}$. This fact makes the estimation procedure much more complex. In order to overcome this difficulty, we propose a methodology that has two steps described as follows. We describe the procedure considering a finite order VLMC. The infinite order case is analogous.

*First step: * Given a sample $z_{1}^{T}$ of the obervable process, for $k$ fixed, which depends on the sample size $T$, estimate the hidden $k_{full}$ Markov chain associated to the order $k$ VLMC, $\hat{\cal{T}}_{k}$, applying the Baum-Welch algorithm considering that the process is the HMM $\left(\boldsymbol{X^{*}},\boldsymbol{Z^{*}}\right)$. 2. 2.

Second step: Generate a bootstrap sample from the estimated tree $\hat{\cal{T}}_{k}$, with transition matrix $\hat{A}^{*}$, and apply a pruning procedure based on a Bayesian Information Criterion to obtain the estimated tree $\hat{\cal{T}}$ of the $\cal{T}$.

In the first step, since the Baum-Welch (BW) algorithm is an EM algorithm, its convergence to a local maximum of the likelihood function is guaranteed [1]. Our propose to make the BW algorithm to reach a global maximum is to vary the required initial guess $\boldsymbol{\lambda}^{*}_{0}$ over a set

$\Lambda_{0}^{*}=\left\{\boldsymbol{\lambda}^{*1}_{0},\ldots,\boldsymbol{\lambda}^{*N}_{0}\right\},$

for a fixed $N$, as big as possible.

Then, for each initial guess, $\boldsymbol{\lambda}^{*j}_{0}$, the Baum-Welch algorithm returns an estimate $\tilde{\boldsymbol{\lambda}}_{j}^{*}\in\Lambda=\{\tilde{\boldsymbol{\lambda}}_{1}^{*},\ldots,\tilde{\boldsymbol{\lambda}}_{N}^{*}\}$ of the vector of parameters $\boldsymbol{\lambda}^{*}=\left(\boldsymbol{A^{*}},\boldsymbol{B^{*}},\boldsymbol{\pi}^{*}\right)$.

Finally, our propose to estimate $\boldsymbol{\lambda^{*}}$ is to take the $\tilde{\boldsymbol{\lambda}}^{*}_{j}$ that maximizes the likelihood $\mathbb{L}\left(\tilde{\boldsymbol{\lambda}}^{*}|Z^{*}\right)$, given the observed sample $Z^{*}=z^{*}_{1},...,z^{*}_{(T-k)+1}$, that is

[TABLE]

In this way, if the likelihood function $\mathbb{L}\left(\tilde{\boldsymbol{\lambda}}^{*}|Z^{*}\right)$ has a finite number of local maxima, for $N$ large enough, our estimator $\hat{\boldsymbol{\lambda}}^{*}$ is the Maximum Likelihood Estimator (MLE) of ${\boldsymbol{\lambda}}^{*}$.

In the estimation procedure, for each $\boldsymbol{\lambda}^{j*}_{0}$, $j\in\{1,\ldots,N\}$, we consider that the values of the noise parameter $\epsilon_{j},j\in\{1,\ldots,N\}$ in the initial emission distribution $\mathbf{B}_{0}(\epsilon_{j})$ assumes a distinct increasing value in the interval $(0,1)$. Then, for each value of the noise parameter, $\epsilon_{j}$, it corresponds an initial vector $\boldsymbol{\lambda}^{*}_{0}(\epsilon_{j})$. Henceforth, $\hat{\boldsymbol{\lambda}}^{*}$ also provides an estimator $\hat{\epsilon}$ of the noise parameter since we can rewrite (9) in terms of $\epsilon_{j}$. This estimator is defined in the following way

[TABLE]

We keep fixed the initial guess of the transition probabilities of the $k_{full}$ (matrix) $\boldsymbol{A^{*}}$, for each $j\in\{1,\ldots,N\}$, as the empirical transition matrix of the observable values $\boldsymbol{Z^{*}}$, truncated at order $k$. And $\boldsymbol{\pi}^{*j}_{0}$ is considered as an uniform distribution.

Now we proceed with the second step of our estimation procedure. Once we have the estimate $\hat{\boldsymbol{\lambda}^{*}}=\left(\hat{\boldsymbol{A^{*}}},\hat{\boldsymbol{B^{*}}},\hat{\boldsymbol{\pi}^{*}}\right)$ of $\boldsymbol{\lambda^{*}}=\left(\boldsymbol{A^{*}},\boldsymbol{B^{*}},\boldsymbol{\pi}^{*}\right)$, associated to the HMM $(\boldsymbol{X^{*}},\boldsymbol{Z^{*}})$, we propose a pruning in the estimated $k_{full}$ tree, $\hat{\cal{T}}_{k}$, with transition matrix $\hat{\boldsymbol{A}^{*}}$, to obtain an estimator for the parameters of the context tree $\cal{T}$, associated to the hidden process $\boldsymbol{X}$. To this end, we apply an adaptation of the Bayesian Information Criterion (BIC) estimator for VLMC, proposed in [6], which is explained in this section. Under some mild conditions, [6] showed that the BIC pruning provides a consistent estimator for a VLMC when the sample comes from a VLMC. However, this is not the case here, since we do not have a sample directly from the hidden process. Then we propose a pruning procedure based on a sample of the estimated $k_{full}$ tree $\hat{\cal{T}}_{k}$, with transition matrix $\hat{\boldsymbol{A}^{*}}$. Thus, our proposal is a bootstrap version of the BIC algorithm, where we replace a sample of the true VLMC by a bootstrap sample ${\hat{x}_{1}^{m}}:=\hat{x}_{1},\ldots,\hat{x}_{m},\ m=O(T),$ drawn from the estimated transition matrix $\hat{\mathbf{A}}^{*}$.

Following [6], we need to define some auxiliary variables. Let $\hat{N}_{m}(\omega,a)$ be the number of occurrences of the string $\omega\in\cup_{j=1}^{k}E^{j}$ followed by the symbol $a\in E$ in the bootstrap sample $\hat{x}_{1}^{m}$, that is

[TABLE]

and the number of occurrences of $\omega$ in $\hat{x}_{1}^{m}$ is

[TABLE]

A feasible bootstrap context tree is such that, given a bootstrap sample $\hat{x}_{1}^{m}$, $d({\cal{T}})\leq$ $D(m)$, $\hat{N}_{m}(\omega)\geq 1$ for all $\omega\in{\cal{T}}$, and $\omega^{\prime}$ is a suffix of some $\omega\in\cal{{T}}$ with $\hat{N}_{m}(\omega)^{\prime}\geq 1$. The set of boostrap feasible context trees is denoted by \mathcal{F}\mbox{\scriptsize(\hat{x}_{1}^{m},D(m))}.

We define the Bootstrap Bayesian Information Criterion (BIC) for a set of feasible trees as

[TABLE]

where $ML_{\cal{T}}(\hat{x}_{1}^{m})=\displaystyle\prod_{\omega\in\tau:\hat{N}_{m}(\omega)\geq 1}\displaystyle\prod_{a\in E}\left(\frac{\hat{N}_{m}(\omega,a)}{\hat{N}_{m}(\omega)}\right)^{\hat{N}_{m}(\omega,a)}$.

For a finite bootstrap sample, $\hat{x}_{1}^{m}$, the BIC estimator of $\cal{T}$ is defined by

[TABLE]

with $D(m)=o(\log m)$.

Since we have replaced the sample of the VLMC by a bootstrap sample we need to show that the BIC estimator is still consistent. To this end, we present the following proposition which plays a key role to prove this consistency.

Proposition 4.1.1**.**

Let $\hat{\mathbf{A}}^{*}$ be a strongly consistent estimator of the transition probability matrix of the markovian process $\boldsymbol{X}^{*}$, with law $\hat{P}$ and transition probabilities $\hat{p}(a,\omega)$, $\forall a\in E,\omega\in\hat{\cal{T}}_{k}$. And let $\hat{x}_{1}^{m}$ be a bootstrap sample of size $m=O(T)$, drawn from $\hat{P}$ fixed, where $T$ is the size of the hidden sample, $x_{1}^{T}$. Then, conditionally on $\hat{P}$, $\forall a\in E,\omega\in\hat{\cal{T}}_{k}$ and for almost all realizations of the VLHMM $\left(\boldsymbol{X},\boldsymbol{Z}\right)$,

[TABLE]

Proof in Appendix.

Now we are ready to state the main result of this work.

Theorem 4.1.1**.**

Let $\hat{x}_{1}^{m}$ be a bootstrap sample of size $m=O(T)$ drawn from $\hat{P}$ fixed. For $d(\cal{T})<\infty$, for the BIC estimator of $\cal{T}$, given by equation (11) we have that

$\hat{\cal{T}}_{BIC}\left(\hat{x}_{1}^{m}\right)=\cal{T},$

almost surely when $m\rightarrow\infty$.

In the general case, we have

\hat{\cal{T}}_{BIC}\big{|_{k}}\left(\hat{x}_{1}^{m}\right)=\cal{T}\big{|_{k}},

almost surely when $m\rightarrow\infty.$

Proof in Appendix.

4.1.1 Computation of the Bootstrap BIC Estimator

The direct application of the BIC procedure is impracticable due to a large number of possible trees to be checked in the likelihood function. Then the estimation of the tree $\cal{T}$, associated to the hidden process $\boldsymbol{X}$, is made in the same manner proposed in [6] through a recursive procedure, based on the CTM algorithm [6], which assigns a value and a binary indicator to each node. The difference in our procedure is that the sample is not directly generated from a VLMC but from an estimated tree $\hat{\cal{T}}_{k}$. In the following we adapt the definitions given in [6], replacing the original sample by the bootstrap sample $\hat{x}_{1}^{m}$. Let

[TABLE]

Then the estimator $\hat{\cal{T}}_{BIC}\left(\hat{x}_{1}^{m}\right)$, defined in equation (11), can be written as

[TABLE]

where $\tilde{\mathbb{P}}_{\omega}\left(\hat{x}_{1}^{m}\right)=m^{-\frac{\left|E\right|-1}{2}}\tilde{\mathbb{P}}_{L,\omega}(\hat{x}_{1}^{m})$.

Definition 4.1.1.1**.**

Given a sample $\hat{x}_{1}^{m}$, let $S_{d}$ be the set of all contexts of maximum size $d=D(m)=o(logm)$ and such that $\hat{N}_{m}(\omega)\geq 1$. For each string $\omega\in S_{d}$ with $\hat{N}_{m}(\omega)\geq 1$, we assign recursively starting from the leaves of the $d_{full}$ tree $\hat{\cal{T}}_{d}$, the value

[TABLE]

and the indicator function

[TABLE]

Definition 4.1.1.2**.**

*For each $\omega\in S_{d}$ the estimated tree $\hat{\cal{T}}$ is the set of contexts $\nu\succeq\omega$ such that

[TABLE]

Proposition 4.1.1.1**.**

The bootstrap context tree estimator $\hat{\cal{T}}_{BIC}\left(\hat{x}_{1}^{m}\right)$ equals the maximizing tree assigned to the root,

[TABLE]

Proof See Appendix.

Remark 4.1.1.1**.**

Once we have the estimate $\hat{\lambda}$, associated with the tree $\cal{T}$, it is straightforward to obtain a Viterbi Algorithm version, adapted to TSCM and TPCM, to estimate an optimal sequence of the hidden process $\boldsymbol{X}$.

5 Simulations and sensitivity analysis concerning the random noise

In this section, we present some simulations to evaluate the methodology proposed in this work. In these simulations, we are interested in evaluating the impact on the estimation of the parameters of the hidden stochastic process $\boldsymbol{X}$, $\boldsymbol{\lambda}$, as we increase the degree of contamination $\epsilon$, for both models TSCM and TPCM. We consider simulations with sample sizes $T=10.000$ and $30.000$ with $100$ Monte Carlo replications. The contamination parameter $\epsilon$ ranged from $0.01$ to $0.99$ with steps of $0.01$. In order to allow such a refinement in the parameter space of the random noise, we decided to use a binary alphabet to decrease the time of simulations. But we stress that there is no restriction on the methodology concerning the use of larger alphabets. According to the proposed methodology, we used the Baum-Welch algorithm for the estimation of the parameters $\boldsymbol{\lambda}^{*}$ and the BIC bootstrap algorithm.

The section is organized as follows: we present two simulation scenarios with very different trees structures of branches. For each scenario, we apply TSCM and TPCM, in order to evaluate the parameter estimates behavior as we increase the degree of contamination of the sample.

5.1 First Scenario

In this scenario, we chose a VLMC $\boldsymbol{X}$ of order $k=3$ with context tree $\cal{T}$ showed in Figure 5.1 and probability transition matrix given in Table 5.1. We chose very different values ​​for the transition probabilities, ranging from 0.05 to 0.87, in order to observe that the behavior of the estimates does not change depending on the value chosen.

For sample sizes 10000 and 30000, we obtain accurate estimates of the parameters. We also notice that, as the contamination noise increases, the estimate of the noise parameter becomes closer to the true value, even for small samples. On the other hand, the variability of the estimates decreases as the sample size increases, as expected, and the estimates become more precise. For the TPCM model, we observe that the transition probabilities of the observed process are increasingly close to zero as the noise parameter increases, even for large samples, which makes the estimation of the noise parameter difficult.

Table 5.3 shows the estimates of transition probabilities for a very small noise, $\epsilon=0.01$, and sample sizes $T=10000$ and $T=30000$. In both cases, the estimates are close to the true values. This was also expected since there was little change in the symbols of the original VLMC since $\epsilon=0.01$. And also note that as the sample size increases, the estimates becomes closer to the true one and variability decreases.

Figure 5.2 shows evidence of normality in the behavior of the estimates of transition probabilities as the sample size increases. We only present results for estimates of transition probabilities $0.27$ and $0.87$ but this behavior remains the same for all other transition probabilities estimates. More than that, we notice that our methodology was able to recover the true tree.

Table 5.4 shows estimates of the transition probabilities of the hidden process $\boldsymbol{X}$ for TSCM. We observe that the variability decreases as the sample size increases and the estimates become more accurate.

Figure 5.3 shows clearly the impact of the increasing of the random noise on the estimates of transition probabilities. For noise values close to $0.50$, although the noise parameters are well estimated, estimates of transition probabilities tend to be distant from the true ones and closer to $0.50$. As a consequence, when the noise is between $0.40$ to $0.60$, the bootstrap BIC algorithm estimates an independent model, ie, a tree with just a root. The problem occurs in the first step of the estimation procedure and not in the boostrap BIC, the Baum-Welch algorithm fails to recover the true transition probabilities in this range of random noise. It is intuitive that in TSCM with $E=\{0,1\}$ the variability of the estimators should attain higher values for noise perturbation around $0.50$ since the emission distribution is bernoulli. The high variability in this interval can lead the Baum-Welch to fail.

Nevertheless, if the value of the estimated noise belongs to the interval $40\%$ to $60\%$ we can conclude that its value is well estimated, but estimates of transition probabilities are far from the true ones. But outside this range, the proposed methodology is able to provide accurate estimates for transition probabilities and also for the context tree.

Table 5.5 shows simulations considering TPCM regime. We observe that the estimates of the transition probabilities are very close to the true ones and become increasingly accurate by increasing the sample size, as in the TSCM. Figure 5.4 shows that, for the TPCM regime, if the contamination is smaller than $0.25$, the estimates of transition probabilities are accurate.

We notice that, as the contamination increases, the estimates become distant from the true value, even for large samples. This is because of higher the noise in this model, most inflated zeros is the contaminated sample and more difficult is to obtain accurate estimates. Again, we only present estimates of transition probabilities $0.27$ and $0.87$ but results are similar for all other values of transition probabilities.

For a sample with size $T=10.000$ we are able to estimate noise values at most $75\%$ but the estimates are not very accurate. The methodology is able to recover the true tree for noise values smaller than $40\%$. After this value, the bootstrap BIC algorithm estimates an independent model (only a root) because all transition probabilities become closer to $50\%$. However, for noise values less than $75\%$, the methodology was able to estimate the parameters of the model and to recover the true context tree $\cal{T}$.

5.2 Second scenario

In this scenario, we chose a tree with a larger order and a more complex tree structure, $\boldsymbol{X}$ is a renewal process. And as in the first scenario, the values of the transition probabilities are very different, ranging from $0.1$ to $0.83$.

Table 5.6 shows the transition probability matrix associated with the process $\boldsymbol{X}$.

The context tree $\cal{T}$ associated to $\boldsymbol{X}$ is shown in Figure 5.5 (order $k=4$).

Table 5.7 shows that we obtain accurate estimates of the parameters when $T=30000$.

According to the Figure 5.6 we observe that the estimates of the transition probabilities are accurate if the random noise is outside the interval $40\%$ to $60\%$, as in the first scenario. Regarding the variability of the estimates, we note that there is a range where the variability increases for a fixed sample size, but it decreases as the sample sizes increases.

Figure 5.7 shows that the estimates of the random noise and transition probabilities present the same behavior shown in the first scenario for the TPCM.

Although the tree structure considered in scenario 2 is more complex, the results were similar to those in scenario 1.

6 Conclusions

In this paper, we have presented a methodology to estimate the parameters of some stochastically contaminated models. These models can be viewed as a bivariate process $(\boldsymbol{X},\boldsymbol{Z})$ where the original, hidden, process $\boldsymbol{X}$ is as a VLMC and $\boldsymbol{Z}$ is the contaminated observed process. We considered two contamination regimes, one regime in which a random noise is added to the original value and the other contamination regime where the original value of the process is multiplied by a random noise. Our inference methodology for the parameters of these contaminated models has two steps. If the tree associated to the VLMC is finite, in the first step we rewrite $\boldsymbol{X}$ as order $k$ Markov chain ($k_{full}$ tree) and apply the Baum-Welch EM algorithm to estimate the parameters of this transformed model. In the second step, we proposed a bootstrap BIC in order to prune the branches of the estimated $k_{full}$ tree in the first step and, in this way, to obtain an estimate of the transition matrix of the hidden VLMC. If the tree associated to the hidden VLMC is infinite, we apply the same methodology to obtain an estimate of the parameters of the hidden tree truncated at some order $k$.

We have shown that our bootstrap BIC estimator for the context tree associated to $\boldsymbol{X}$ is strongly consistent under some mild conditions. We have presented simulations showing that our methodology is capable of recovering the hidden tree and the noise parameters from a contaminated sample. For samples sizes above $10.000$ the accuracy of the estimate of the noise parameter is quite satisfactory and the estimates of the transition probabilities associated to the hidden VLMC are close to the true values, with low variability, in a reasonable range of random noises, namely out of $40$ to $60\%$, in the additive model, and up to $25\%$, in the multiplicative model. Hence, if the estimate of the noise parameter is outside these ranges we can conclude that the estimates of transition probabilities associated to the hidden VLMC are reliable.

Although the simulations have been made considering an alphabet $E=\{0,1\}$, in order to decrease the time of simulations, the method can be applied to any type of emission distribution with any discrete alphabet.

7 Appendix: **Proofs **

Proof of Proposition 3.1.1

Proof.

Let $\boldsymbol{Z}$ be a contaminated process according to a TSCM. Without loss of generallity, we take $t=0$ and for some $\omega=a_{-l(\omega)+1}^{0}\in\cal{T}$, with $l(\omega)<T$, we have that

[TABLE]

The event $\left\{Z_{0}=z_{0}\right\}$ can be written in terms of $\boldsymbol{X}$ and $\boldsymbol{\xi}$, according to a TSCM, as

$\left\{Z_{0}=z_{0}\right\}=\displaystyle\bigcup_{x_{0},b_{0}=0:\atop z_{0}=x_{0}\oplus b_{0}}^{\left|E\right|-1}\left\{X_{0}=x_{0},\xi_{0}=b_{0}\right\}.$

Henceforth

[TABLE]

Note that $\left\{X_{0}=x_{0},X_{0}=a_{0}\right\}$ are empty sets if $x_{0}\neq a_{0}$, then

[TABLE]

Hence, by independence of $\boldsymbol{X}$ and $\boldsymbol{\xi}$, we have that

[TABLE]

On the other hand, we have that

[TABLE]

Since the events $\left\{X_{0}=x_{0},X_{0}=a_{0}\right\}$ are empty for all $x_{0}\neq a_{0}$, then the only remaining event is $\left\{X_{0}=a_{0}\right\}$. We notice that the events $\left\{X_{0}=a,\xi_{0}=b_{0}\right\}$, for each $a_{0},b_{0}\in E$ fixed, are mutually exclusive and $\boldsymbol{X}$ is independent of $\boldsymbol{\xi}$, then

[TABLE]

This concludes the proof of item (i).

(ii) We want to show that the likelihood function of the observed process $\boldsymbol{Z}$, for a sample $z_{1}^{T}$, is:

[TABLE]

Like in item (i) in Proposition 3.1.1 we can write the events $\left\{Z_{t}=z_{t}\right\}$ in terms of $\boldsymbol{X}$ and $\boldsymbol{\xi}$,

[TABLE]

By the distributive property $A\cap\left\{B\cup C\right\}=\left\{A\cap B\right\}\cup\left\{A\cap C\right\}$, we have that

[TABLE]

Since $\left\{X_{t}=a_{t},\xi_{t}=b_{t}\right\}$ are mutually exclusive,

[TABLE]

Finally, the claim follows by independence of $\boldsymbol{X}$ and $\boldsymbol{\xi}$. ∎

Proof of Proposition 3.2.1

Proof.

Proofs of items (i),(ii) are analogous to the proof of proposition 3.1.1, but changing the indicator function of $a\oplus b$ by $a\cdot b,\ \forall a,b\in E$. ∎

Proof of Proposition 4.1.1

Proof.

Let $\hat{\mathbf{A}}^{*}$ be a strongly consistent estimator of the transition probability matrix of the markovian process $\boldsymbol{X}^{*}$, with law $\hat{P}$. We observe that if each entry of the transition probability matrix $\hat{\mathbf{A}}^{*}$ , $\hat{p}(a|\omega)$, is a MLE of the transition probability of the hidden Markov chain $\boldsymbol{X^{*}}$, ${p}(a|\omega)$, $\forall a\in E,\omega\in\cal{T}$. Then for almost all realizations of the process $(\boldsymbol{X},\boldsymbol{Z})$, we have that $\hat{p}(a|\omega)\longrightarrow p(a|\omega)$ almost surely as $m=O(T){\rightarrow}\infty$, since the regularity conditions $A_{1}$ to $A_{6}$ in [12] are satisfied. Hence, we only have to show that for a bootstrap sample, $\hat{x}_{1}^{m}$, of size $m=O(T)$, drawn from $\hat{P}$ fixed, for almost all realizations of the process $(\boldsymbol{X},\boldsymbol{Z})$, the following holds

[TABLE]

almost surely as $m\longrightarrow\infty$.

But since we can write

[TABLE]

then the random variable $\frac{\hat{N}_{m}(\omega,a)}{m}$, conditionally in $\hat{P}$, converges almost surely to

$E(1_{\left\{\hat{x}_{t}^{t+k}=\omega,\hat{x}_{t+k+1}=a\right\}}|\hat{P})=\hat{P}({\omega a})$, as $m\longrightarrow\infty$

by the Ergodic Theorem, where $\hat{P}(\omega a)$ is the measure of the string $\omega a$ given $\hat{P}$. Analogously, we have that

[TABLE]

almost surely as $m\longrightarrow\infty$.

Then, from 17 and 18 the result in 16 follows.

∎

Proof of the Theorem 4.1.1

Proof.

Lemmas 3.1 and 3.2 presented in [6] guarantee the consistency of the BIC estimator $\hat{\cal{T}}$ in the case where the sample is obtained directly from a VLMC with tree $\cal{T}$. The only difference in our case is that we have replaced the variable $N_{T}(\omega,a),\omega\in\cal{T}$, $a\in E$ in [6], which counts the frequency of the string $\omega$ followed by the symbol $a$ in the sample ${x_{1}^{T}}$, by its bootstrap version $\hat{N}_{m}(\omega,a)$, but applying Proposition 4.1.1 to Lemmas 3.1 and 3.2, with this replacement, Theorem 4.1.1 follows. ∎

Proof of the Theorem 4.1.1.1

Proof.

Analogously, Propositions 4.3 and Lemma 4.4 presented in [6] guarantee the that $\hat{\cal{T}}_{BIC}(x_{1}^{n})={\cal{T}}^{d}_{\emptyset}$ if the sample is obtained directly from a VLMC with tree $\cal{T}$. Proposition 4.1.1 allows us to replace the original sample by a bootstrap sample, implying proposition 4.1.1.1, by replacing $N_{T}(\omega,a)$, $\omega\in\cal{T}$, $a\in E$, in [6] by its bootstrap version $\hat{N}_{m}(\omega,a)$ in Propositions 4.3 and Lemma 4.4 presented in [6]. ∎