## 2 Materials and Methods
### Glottal Flow Model for the Kalman Filter  
Rosenberg Model for the Glottal Pulse  
A parametric model of the glottal pulse can be obtained from the Rosenberg model [46], which can be formulated as [51]:  
\[g[n]=\left\{ \begin{array}{ll}0.5[1-\cos(\pi(n+1))/N_{1}],&0\leq n\leq N_{1}-1, \\ \cos(0.5\pi(n+1-N_{1})/N_{2}),&N_{1}\leq n\leq N_{1}+N_{2}-1,\\ 0,&\text{otherwise},\end{array} \right.\]  
where \(N_{1}\) is the number of samples of the opening phase and \(N_{2}\) is the number of samples of the closing phase. For a sequence of 96 samples (equivalent to approx. 210 Hz fundamental frequency, pitch period of 4.8 msec., and sampling frequency \(f_{s}=20\) kHz), with \(N_{1}=25\) and \(N_{2}=10\), the \(z\)-transform \(G(z)\) has the form:  
\[G(z)=z^{-33}\prod_{k=1}^{33}(-b_{k}^{-1})\prod_{k=1}^{33}(1-b_{k}z), \tag{12}\]  
where \(b_{k}\) corresponds to the zeros of \(G(z)\), which can also be written in the following form:  
\[ \begin{split} G(z)&=g[0]+g[1]z^{-1}+g[2]z^{-2}+\cdots +g[N-1]z^{-(N-1)},\\ &=\beta_{0}+\beta_{1}z^{-1}+\beta_{2}z^{-2}+\cdots+\beta_{N-1}z^{- (N-1)},\\ &=\sum_{k=0}^{N-1}\beta_{k}z^{-k}.\end{split}  \tag{13}\]  
The glottal pulse time-domain waveform \(g[n]\) and its spectrum are plotted in Figure 1  
The periodic comb excitation \(p[n]\) is modeled as one-sided quasi-periodic impulse train:  
\[p[n]=\sum_{k=0}^{\infty}\gamma^{k}\delta[n-kN_{p}], \tag{14}\]  
which has \(z\)-transform:  
\[P(z)=\sum_{k=0}^{\infty}\gamma^{k}z^{-kN_{p}}=\frac{1}{1-\gamma z^{-N_{p}}}, \tag{15}\]  
where \(N_{p}=f_{s}/f_{0}\) (fundamental period in samples) and \(\gamma\) is a number close to 1 (e.g., 0.999) to make the filter stable. The spectrum of the periodic input \(P(z)\) has a fundamental frequency of \(f_{0}=210\) Hz (\(N_{p}=96\)). Figure 1: Rosenberg model in time domain (only first 50 samples shown, **top** panel) and the magnitude of its spectrum (**bottom** panel). Therefore, \(P(z)G(z)\) is the \(z\)-transform of the glottal flow model (spectrum shown in Figure 2). In the time-domain, the GVV signal can be represented by an ARMA model that can be constructed as a shaping filter (\(sf\)) driving the canonical MA model (see Equation (16)) [43,44]:  
\[x_{sf}(n)=-\sum_{k=1}^{p}\alpha_{k}x_{sf}(n-k)+\sum_{k=0}^{q}\beta_{k}w_{2}(n-k), \tag{16}\]  
where \(x_{sf}(n)\) is the state of the shaping filter, \(\alpha_{k}=-\gamma_{k}\) and \(\beta_{k}\) are the \(k\)th coefficient of the AR and MA model, respectively, and \(w_{2}(n)\) is Gaussian noise with mean 0 and variance \(\sigma_{w_{2}}^{2}\). The state-space equation for this model is:  
\[\mathbf{x_{SF}}(n+1) =\mathbf{A_{SF}}\mathbf{x_{SF}}(n)+\mathbf{B_{SF}}w_{2}(n), \tag{17}\] \[w_{1}(n) =\mathbf{C_{SF}}\mathbf{x_{SF}}(n), \tag{18}\]  
where \(\mathbf{x_{sf}}(n)=(x_{SF}(n-p+1)\quad x_{SF}(n-p+2)\quad\dots\quad x_{SF}(n) )^{T}\) is the state vector and \(p\) is the order of the AR model. Since the periodic input has \(N_{p}\) poles, the order of the AR model is \(p=N_{p}\). \(\mathbf{A_{SF}}\) is the transition matrix \(p\times p\):  
\[\mathbf{A_{SF}}= \begin{bmatrix}0&1&0&\dots&0&0\\ 0&0&1&\dots&0&0\\ \vdots&\vdots&\ddots&\vdots&\vdots&\vdots\\ 0&0&0&\dots&0&1\\ -\alpha_{p}&-\alpha_{p-1}&-\alpha_{p-2}&\dots&-\alpha_{2}&-\alpha_{1}\end{bmatrix}  \in\mathbb{R}^{p\times p},\]  
\(\mathbf{B_{SF}}=[0\quad 0\dots 1]^{T}\in\mathbb{R}^{p\times 1}\) and \(w_{2}(n)\) is a stochastic driving noise with zero mean and variance \(\sigma_{w_{2}}^{2}\). The MA equation Equation (18) contains \(\mathbf{C_{SF}}=[\beta_{q}\quad\beta_{q-2}\dots\beta_{1}\quad\beta_{0}]\in \mathbb{R}^{1\times(q+1)}\) and the colored noise \(w_{1}(n)\in\mathbb{R}\) is the dot product of \(\mathbf{C_{SF}}\) and \(\mathbf{x_{SF}}(n)\). Considering the source-filter theory, the colored noise model can be considered as modeling the GVV for the [PERSON] implementation of IBIF in Equation (19).