随笔- 9 文章- 0 评论- 5 阅读- 5183

HTS语音合成训练（三）-理论基础

HTS语音合成训练（三）- 总结

背景简介：

前两个月由于项目需要，进行了语音拼接合成相关技术的调研和HTS技术实践，从最初的看各种论文和找各种资料来了解相关技术，到实践HTS语音合成，调用hts_engine API实现模型参数提取供同事进行loss计算来选取最优片段，走了很多弯路，在此记录和总结TTS相关经验

1.技术发展路线简介

1.1.概述

语音合成（text to speech），简称TTS，主要作用是将文字转语音，先后经历了共振峰合成、LPC 合成、PSOLA 技术合成、HMM(Hidden Markov Model)语音合成、神经网络合成等方法

1.2.发展过程

相关资料：

论文：

共振峰语音合成算法研究与实现

1.2.1.相关概念

共鸣：物体因共振而发声的现象，例如两个频率相近的音叉靠近，其中一个振动，另一个也会发声，如下图，（好像中学时候做过的实验？）

**共振(resonance)：**物理学上指特定频率和波长下，比其他频率和波长以更大的振幅做振动的情形，这些特定频率和[波长](https://baike.baidu.com/item/波长/829184)称之为共振频率和共振波长，一个系统有多个共振频率，在共振频率上振动会更加容易，在其他频率上振动更加困难

元音：英语音标就是记录音素的符号，是音素的标写符号，它的制定原则是：一个音素只用一个音标表示，一个音标只表示一个音素。英语是拼音文字，字母就表达语音。然而英语有40多个发音，却仅有26个字母。为了能准确标识发音，引入了语音符号，这些语音符号称作音标（Phonetic symbol）。国际音标分为两种：英式音标（DJ音标）和美式音标（K.K.音标）

浊音：当声带处于收紧状态，流经的气流使声带振动，这时产生的音是浊音，浊音是有具体意义的，例如元音，

浊辅音：浊辅音则不是由气流发出的声音/b/

清音：当声带处于放松状态，不伴有声带振动的音，称为清音，清音是无具体意义的音，如爆破音、摩擦音

清辅音：清辅音则是由于气流发出的声音，如/p/是清辅音，/b/是浊辅音，两者发音的口型相同，但是将手放于嘴前，会感受到/p/发出的气流，而/b/则不会

激励源与滤波器：气流与声门，等效为激励源，对应的是激励频率，也就是基音频率F0，声道等效为滤波器，对应是谐振频率

谐振频率：由于电容和电感充放电呈现出周期性特性，这个周期的频率就是”谐振频率“，在串联震荡回路中谐振频率下阻抗会迅速减小，并联振荡回路中，谐振频率下阻抗会迅速增大，由此可以形成各种不同的滤波器

以图9—3(a)所示的滤波电路来说，当有信号从左至右传输时，L对低频信号阻碍小，对高频信号阻碍大C则对低频信号衰减小，对高频信号衰减大。因此该滤波电路容易通过低频信号，称为低通滤波电路。其特点可用图中的幅一频(UF特性f}}I线表示
对于图9—3(b)所示的滤波电路来说，容易通过高频信号，所以称为高通滤波电路
对于图9—3(c)所示的滤波电路，它利用C l和L1串联对谐振信号阻抗小、C2和L7并联对谐振信号阻抗大的特性，能让谐振信号f容易通过，而阻碍其他频率信号通过，所以称为带通滤波电路。该电路的这种特点可用图中的幅一频(U-F特性曲线概括
对于图9—3(d)所示的滤波电路，它利用Cl和Ll并联对谐振信号阻抗大、C，和L，串联对谐振信号阻抗小的特点，容易让谐振频率以外的信号通过，而抑制谐振信号厂F通过，所以称为带阻滤波电路。该电路的特点可用图中的幅一频(U-F性曲线来概括

阻抗：主要由三部分组成，电阻，电容容抗和电感感抗，电阻对直流电起阻碍作用，电容和电感对交流电起阻碍作用，容抗隔直流，通交流，感抗通直流，阻交流

容抗计算公式：C表示电容，单位是法拉（F），f表示频率

X_{C} = 1 / 2 π f C

$X_{\mathrm{C}}=1 / 2 \pi f C$

感抗计算公式：L表示电感，单位是亨（H），f表示频率

X_{L} = 2 π f L

$X_{\mathrm{L}}=2 \pi f L$

1.2.2.声学基础和共鸣原理

1.2.2.1.声道共鸣机制

咽腔、口腔和鼻腔共同构成了声道，他们都是共鸣的器官，喉音（声门气流）自声带产生后，向上进入声道，经过对气体分子的压缩和稀释，声学特性（声道共鸣曲线）得到改变并产生输出声波

1.2.2.2.共振峰

通过构造器官的活动改变声道的大小和形状，声道的共鸣特性发生变化，因而在声音频谱中的一些频率得到了共振加强，而另一些则被削弱减幅，这些被加强的共振频率被称为共振峰，了解共振峰之间的相互关系对理解元音尤其重要，不同的元音对应不同的共振峰，对应不同的声道形状

1.2.2.3.构音三要素

在声道内，构音主要由三个要素构成：

舌的收缩位置
舌的收缩程度
嘴唇的收缩

所有的元音和辅音都是根据声道形状和嘴唇的变化而形成

1.2.3.源-滤波器

参考论文：

《基于HMM的单元挑选语音合成方法研究》

源-滤波器语音合成方法是基于声音是通过激励和滤波器产生的声学理论。整个过程分为：经过激励得到声源以及在声道中传播信号。激励源主要包括周期性激励和非周期性激励，其中周期性激励用来产生浊音信号，非周期性激励也就是类似于噪声激励，用来产生清音信号，当同时使用这两种激励时会产生某些浊辅音信号，合成声音时根据语音的特征不同，选择不同的激励源，通过改变声道参数进而改变合成语音的音色。源-滤波器语音合成主要有两种实现方法：共振峰合成与线性预测分析合成，二者的区别在于声道模型不同

1.2.4.共振峰合成

参考论文：

《基于HMM的单元挑选语音合成方法研究》
《共振峰语音合成算法研究与实现》

20 世纪出现电子器件以来，语音合成技术就得到了快速的发展。1939 年 H.Dudley 制作出了第一个电子合成器，它是利用共振峰原理制作的合成器，通过一些类似于白噪声的激励产生出非浊音信号，而通过周期性的激励产生出浊音信号，其中的模拟声道是一个10阶的带通滤波器，最后通过人为控制可以合成一些简单易懂的语音，共振峰模型是基于对声道的一种比较准确的模拟，该模型中语音生成的过程是在激励信号的激励下，声波经谐振腔（声道），由嘴或鼻辐射声波，声道传输频率响应上的极点称之为共振峰，而语音的共振峰频率（极点频率）的分布决定着声音的音色，音色各异的语音具有不同的共振峰模式，因此，以每个共振峰频率及其带宽作为参数，可以构成共振峰滤波器，再用若干个这种滤波器的组合来模拟声道传输特性（声道模型），对激励源发出的信号进行调制，再经过辐射模型就可以得到合成语音

虽然使用共振峰合成语音可以描述语音中最基本的部分以及描述自然流中的各种现象，但是它却有明显的缺点，因为它将对声道的模拟作为基础，如果对声道模型的模拟不够准确，将会大大的影响合成语音的质量，因为共振峰参数并不能表征语音的自然度，所以合成出来的语音自然度不高，另外对于共振峰合成器的控制也比较复杂，一个好的合成器的控制参数高达几十个，实际操作起来也非常困难

1.2.5.拼接合成

参考论文：

《基于HMM的单元挑选语音合成方法研究》

最开始的拼接合成是采用语音的编码技术，将一些语音编码或者直接存储起来得到语料库，在合成时解码组合输出即可，由于存放的是单个字的语音，语音逐字发音，语音单一，而由于人实际说话的自然流并不是各个孤立语音的简单拼接，所以整体效果的自然度不高，一般主要用于自动报时、报站等应用

而真正的波形拼接合成方法并不是单纯的将语音拼接在一起，它还注重语音合成时的韵律控制与修改，尤其是提出基音同步叠加方法后，使得拼接合成技术的质量提高，该方法在拼接语音之前，对拼接单元按照上下文环境进行韵律的调整，使得语音的韵律符合上下文要求，大大提高了语音合成的自然度和可懂度

随着计算机技术的迅猛发展以及海量数据库的提出，近年来大量使用基于大语料库的波形拼接合成方法，其核心思想是将海量的原始语音信号以及韵律信息存储到机器中，构成一个庞大的语音库，在合成时根据待合成信息选取适当的算法规则，从语音库中挑选合适的波形单元，经过波形编辑拼接以及平滑处理，并为其添加合适的上下文韵律信息，从而使得合成的语音具有较高的自然度。大型数据库的引入使波形拼接语音合成技术有了很大的进步，合成出的语音效果也有了质的突破，使用该技术既保持了原始语音的主要音段特征，又能在拼接时灵活地调整语音的音高、时长等超音段特征，使得合成的语音非常接近人的原始声音，因此该技术被人们广泛研究。

虽然波形拼接合成方法使得语音的自然度有了一定的提高，但是它也存在很大的缺陷，即，要想合成出大量的高质量语音，则需要一个庞大的存储空间来存放各学科的语音信息以及韵律信息，所以不易在嵌入式设备上实现，而且工作量非常大。另外，合成出来的语音非常依赖于音库的制作质量，存在一定的稳定性问题，要想改变语音的特征，则需要建立新的语音库，因此就加大了成本的需要，除此之外合成单元调整策略的不完备还会导致合成的语音的不连续

1.2.6.参数合成

参考论文：

《基于HMM的单元挑选语音合成方法研究》

由于波形拼接合成方法存在的一些缺点，因此基于参数的语音合成方法再次得到了重视，尤其是 20 世纪末以来，基于统计声学建模的语音合成技术迅速发展起来并且得到了广泛的应用，其中隐马尔可夫模型（HMM）是最常用的声学统计模型形式。目前，HMM 是语音领域中使用最广泛的方法，早期它主要被用于语音识别方面，并且技术已经非常成熟。后来提出将 HMM 用于合成领域中也得到了成功的示例，该模型有自动化构建程度高以及合成效果稳定等优点，已经成为语音合成技术研究的热门

基于 HMM 的语音合成系统结构简单，系统的训练基本上不需要语言学知识的指导，构建时间短，构建过程也不需要人为过多干涉。整个系统主要分为模型的训练和语音的合成，要训练 HMM 首先需要对语料库中的语音进行自动切分，然后进行参数（如，基频、时长或频谱参数）提取，建立隐马尔可夫模型，合成语音时首先要对合成文本进行文本分析处理，再根据 HMM 中的参数生成算法得到合成文本的参数，最后通过参数合成器得到语音。由于系统中并没有直接存储语音波形，而是存储语音的特征参数，所以需要的存储空间较小，故该系统也非常适合应用于资源受限的嵌入式设备上。虽然基于 HMM 的参数合成法相比传统的参数合成法简单易用，但是由于最终的语音还是通过参数合成得到的，因此生成的语音与自然语音相比有着明显的差距，合成结果具有很浓的“机器味”，人们难以接受

1.2.7.神经网络合成

参考：

牛津大学TextToSpeech

在HMM参数合成中，特征提取是使用传统决策树方法，声学模型是使用GMM-HMM，有了神经网络之后，对特征提取可以使用DNN网络来替代，HMM由RNN或者LSTM来代替，由此HMM被神经网络所取代

1.2.8.端到端合成

端到端的TTS深度学习模型tacotron(中文语音合成)

在引入神经网络之后，语音合成的基本流程和结构没有被改变，需要前端处理->声学模型->声码器这一套流程，但是，直到google推出tacotron端到端模型，将前端处理部分全部包含在模型中，模型输入为一段文字向量，输出为频谱，然后频谱通过Griffin-Lim算法生成音频

2.原理篇

2.1.特征提取

在HMM语音合成系统HTK中，需要用到两大特征来进行训练：

语言学特征
声学特征

语言学特征主要是文本的上下文关系如前后单词是什么，前后音素是什么等等，还有强制对齐之后音素的边界，声学特征主要是基频、MFCC特征和bap(Band Aperiodicity)特征

2.1.1.基音周期与基音频率

参考：

语音基音周期估计基础

基音决定了语音的音调。汉语是一种有调的语音，同一句话，语调不同，意义差别非常大，准确的检测语音的基音周期，有利于“理解”话语中的意思

按照经典的语音激励模型，声带周期性的张开和闭合所需要的时间，就是指基音的周期（假如声带某次从开启到闭合的时间为4ms，那么，它的基音频率为250Hz）。很显然，要准确的检测基音周期，需要尽可能的先去掉声道模型和辐射模型的干扰

一般来说，男性的基音频率约为50Hz_{250Hz，女性的基音频率约为：100Hz}500Hz，所以我们通常把基音频率的上限限制为500Hz，同时，又由于50Hz对应的是交流电的工频频率，难以去掉它的噪声干扰，所以这个频率检测时也不予考虑。综上，基音频率的检测范围在60Hz_{500Hz或者70Hz}500Hz比较好。

基音频率：即后面所要讲到的傅里叶变换后的一次谐波分量的频率，可在文档中搜索基频查看

2.1.2.MFCC(Mel Frequency Cepstral Coefficents)特征

参考文献：

此部分涉及理论比较多，需要弄懂FTT（快速傅里叶变换）、同态信号处理、DCT(discrete cosine transform，离散余弦变换)等概念，下面是对参考内容的简单总结，如需详细了解可以看参考文献

2.1.2.1.FTT快速傅里叶变换

我们通常看到的语音波形图都是在时域上，也就是x轴是时间t，y是对应的值，跟通常函数类似，我们知道函数f(x)是可以通过泰勒公式进行展开的，如下

f (x) = f (x_{0}) + f^{'} (x_{0}) (x - x_{0}) + o (x - x_{0})

$f(x)=f\left(x_{0}\right)+f^{\prime}\left(x_{0}\right)\left(x-x_{0}\right)+o\left(x-x_{0}\right)$

如果 $x_0$ 取0的话，就变成了关于x的n阶多项式了，相应的x的一次、二次一直到n次方都会有一个系数如 $(x-x_0)$ 的系数为 $f^{\prime}(x_0)$ ，如果我们用 $x^1$ 到 $x^n$ 次来代替时间轴上的1、2、3等时间轴的值，就得到了x轴上是x的n次方的表示形式，y的值也就变成了x的次方前面的系数，这样就将f(x)转变为了非时域形式，但是同样可以表征f(x)，这里的f(x)是以 $x$ 、 $x^1$ ... ${x^n}$ 来表示，称为幂级数

f(x)如果是周期函数的话，还可以用三角函数如sin,cos来进行表示：

f (t) = A_{0} + \sum_{n = 1}^{\infty} [a_{n} \cos (n ω t) + b_{n} \sin (n ω t)]

$f(t)=A_{0}+\sum_{n=1}^{\infty}\left[a_{n} \cos (n \omega t)+b_{n} \sin (n \omega t)\right]$

也可以写作：

f (x) = a_{0} + \sum_{n = 1}^{\infty} (a_{n} \cos (\frac{2 π n}{T} x) + b_{n} \sin (\frac{2 π n}{T} x)), a_{0} \in R

$f(x)=a_{0}+\sum_{n=1}^{\infty}\left(a_{n} \cos \left(\frac{2 \pi n}{T} x\right)+b_{n} \sin \left(\frac{2 \pi n}{T} x\right)\right), a_{0} \in \mathbb{R}$

也叫做傅里叶级数，类比幂级数，这里展开成了关于cos和sin相关的式子，也就是x轴上的值变成了sin和cos，也叫做基，与x的系数类似，sin和cos前面都有对应的系数，如：

f (x) \approx 1 + \frac{4}{π} \sin (x) + 0 \sin (2 x) + \frac{4}{3 π} \sin (3 x) + 0 \sin (4 x) + \frac{4}{5 π} \sin (5 x)

$f(x) \approx 1+\frac{4}{\pi} \sin (x)+0 \sin (2 x)+\frac{4}{3 \pi} \sin (3 x)+0 \sin (4 x)+\frac{4}{5 \pi} \sin (5 x)$

此时对周期函数的拟合效果：

此时基，也就是展开后的与幂级数中 $x$ ... $x^n$ 对应的值为

{1, \sin (x), \sin (2 x), \sin (3 x), \sin (4 x), \sin (5 x)}

$\{1, \sin (x), \sin (2 x), \sin (3 x), \sin (4 x), \sin (5 x)\}$

他们对应的系数为

(1, \frac{4}{π}, 0, \frac{4}{3 π}, 0, \frac{4}{5 π})

$\left(1, \frac{4}{\pi}, 0, \frac{4}{3 \pi}, 0, \frac{4}{5 \pi}\right)$

在频域中表示为

则将周期函数f(x)时域图，通过傅里叶级数变换成了频域图，不管在频域还是时域，都表示同一个函数

那么上面的傅里叶级数和傅里叶变换有什么关系呢？

由于傅里叶级数只能表示周期函数，如下的非周期函数将不能表示

但是可以将此函数的周期T看作是无穷大

T = 2 π \to T = \infty

$T=2 \pi \rightarrow T=\infty$

则在T变化的同时函数图像变化如下

由于函数的基与周期的关系为 $cos(\frac{2n \pi}{T}x)$ 或者 $sin(\frac{2n \pi}{T}x)$ ，当 $T$ 逐渐增加时，基是逐渐减小的，则基的变化则如下图所示

这些频率会变得稠密，直到连续，形成连续曲线

上面这个对非周期函数让 $T \rightarrow \infty$ 的形成连续频域的过程就是傅里叶变换

由以上的傅里叶级数公式，通过欧拉公式变换（欧拉公式可以将f(x)转为复数域，推导过程可以见参考文献），复数形式的傅里叶级数如下：

f (x) = \sum_{n = - \infty}^{\infty} c_{n} \cdot e^{i \frac{2 π n x}{T}}

$f(x)=\sum_{n=-\infty}^{\infty} c_{n} \cdot e^{i \frac{2 \pi n x}{T}}$

其中 $c_n$ 为

c_{n} = \frac{1}{T} \int_{x_{0}}^{x_{0} + T} f (x) \cdot e^{- i \frac{2 π n x}{T}} d x

$c_{n}=\frac{1}{T} \int_{x_{0}}^{x_{0}+T} f(x) \cdot e^{-i \frac{2 \pi n x}{T}} d x$

以上是傅里叶级数的复数形式，推导过程可见参考文献傅里叶级数到三角级数的详细推导

同幂级数时一样分析，现在基变成了 $e^{i \frac{2 \pi n x}{T}}$ ，也就是x轴坐标不再是时间t，也不是 $x^1$ ... $x^n$ ，也不是 $\cos \frac{2 \pi n}{T} x$ 和 $\sin \frac{2 \pi n}{T} x$ ，而是 $e^{i \frac{2 \pi n x}{T}}$ ，而基所对应的系数变为了 $c_n$ ，也就是频谱图中所对应的纵坐标的值（频率HZ），再将 $c_n$

表示为关于 $w=2 \pi/T$ 的函数 $F(w)$

\begin{array}{r} c_{n} = \frac{1}{T} \int_{x_{0}}^{x_{0} + T} f (x) \cdot e^{- i \frac{2 π n x}{T}} d x \\ T = \infty \end{array}} ⟹ F (ω) = \frac{1}{2 π} \int_{- \infty}^{\infty} f (x) e^{- i ω x} d x

$\left.\begin{array}{r} c_{n}=\frac{1}{T} \int_{x_{0}}^{x_{0}+T} f(x) \cdot e^{-i \frac{2 \pi n x}{T}} d x \\ T=\infty \end{array}\right\} \Longrightarrow F(\omega)=\frac{1}{2 \pi} \int_{-\infty}^{\infty} f(x) e^{-i \omega x} d x$

$F(w)$ 就是傅立叶变换(FT)，得到的就是频域曲线

而反傅里叶变换(IFT)公式则为

\begin{array}{r} f (x) = \sum_{n = - \infty}^{\infty} c_{n} \cdot e^{i \frac{2 π n x}{T}} \\ T = \infty \end{array}} ⟹ f (x) = \int_{- \infty}^{\infty} F (ω) e^{i ω x} d ω

$\left.\begin{array}{r} f(x)=\sum_{n=-\infty}^{\infty} c_{n} \cdot e^{i \frac{2 \pi n x}{T}} \\ T=\infty \end{array}\right\} \Longrightarrow f(x)=\int_{-\infty}^{\infty} F(\omega) e^{i \omega x} d \omega$

如果需要知道详细推导过程，请看参考文献傅里叶级数（FS）、傅里叶变换（FT）快速傅里叶变换（FFT）及量子傅里叶变换（QFT）之间推导关系

由以上的傅里叶变换式子可以看出，如果想要得到 $F(\omega)$ ，需要对 $f(x)$ 的值与 $e^{-i \omega x}$ 的值做积，由于上面的假设是 $d x$ 取 ${-\infty}$ 到 ${+\infty}$ ，而真正计算机在计算的时候并不能取到无穷，所以上述式子会变成离散形式

F (w) = \sum_{t = 0}^{T - 1} f (t) e^{- i w t}

$F(w)=\sum_{t=0}^{T-1} f(t) e^{-i w t}$

也就是离散傅里叶变换(DFT)的函数式

如果计算上述式子，例如 $T$ 取8，则需要计算8个 $f(t)$ 与8个 $e^{-i w t}$ 相乘才能得出 $F(w)$ 的值，即频谱图，此时计算效率为 $O(n^2)$ ，需要消耗大量的计算机资源，改进点就是利用 $e^{-i w t}$ 在单位圆上的性质，利用分治法进行改进，具体改进方式可查看参考文献中十分简明易懂的FFT（快速傅里叶变换，改进后的时间复杂度为 $O\left(n \log _{2} n\right)$ ，也就是快速傅里叶变换

2.1.2.2.同态信号处理

同态信号处理作用之一是从混响中提取真实的声音，那什么是混响呢？

混响主要用于唱卡拉OK,增加话筒声音的延时,产生适量的回声,使唱歌的声音更圆润更优美,歌声不那么“干”

我们的目的在ktv唱歌的时候提取出真实的声音，而不是混响

如果在时域中，可以将真实声音信号设为 $s(n)$ ，混响信号设为 $x(n)$ ，由于在房间内会产生回波，所以回波信号设为 $h(n)$ ，假设回波信号的方程已知，混响信号由回波信号和真实信号组成，方程如下

x (n) = s (n) * h (n)

$x(n)=s(n) * h(n)$

其中，*代表卷积，卷积公式如下

(f ⋆ g) (t) = \int_{τ} f (τ) g (t + τ) d τ

$(f \star g)(t)=\int_{\tau} f(\tau) g(t+\tau) d \tau$

时域中如果已知 $x(n)$ 和 $h(n)$ 求 $s(n)$ 不好求得，所以将其通过傅里叶变换转为频域中，因为时域中卷积对应这频域中乘积，如果不懂可以见参考文档中为什么频域中相乘对应这时域中卷积？，会有详细推导，于是公式变成了如下

X (z) = S (z) \times H (z)

$X(z)=S(z) \times H(z)$

然后对频域中函数取对数得到

X^{'} (z) = \log X (z) = \ln S (z) + \log H (z)

$X^{\prime}(z)=\log X(z)=\ln S(z)+\log H(z)$

这样就将频域中的乘法变换为加法，然后再进行反离散傅里叶变换

IDFT (\log (X (k))) = IDFT (\log (S (k))) + I D F T (\log (H (k)))

$\operatorname{IDFT}(\log (X(k)))=\operatorname{IDFT}(\log (S(k)))+I D F T(\log (H(k)))$

假设得到的时域信号如下

x^{'} (n) = s^{'} (n) + h^{'} (n)

$x^{\prime}(n)=s^{\prime}(n)+h^{\prime}(n)$

此时的时域信号 $x^{\prime}(n)$ 则为倒谱，由于 $x^{\prime}(n)$ 和 $h^{\prime}(n)$ 是已知的，所以可以得到 $s^{\prime}(n)$ 的值，也就分离出了我们想要的真实声音信号的倒谱

2.1.2.3.DCT离散余弦变换

本节可以看参考文献DCT离散余弦变换

这篇文章主要让我懂了些关于高频和低频的知识，当我们在音频信号处理中，需要对频谱进行低频滤波处理获得音频信号的包络，因为音频信号的低频部分代表了包络而高频部分代表了频谱的细节，为什么呢？

个人理解，例如我们在进行频谱转换的时候，如果是采用泰勒公式展开，低频部分也就相当于x的一次方、二次方、三次方这样的值，高数里面一般对某个值做极限近似处理的时候，一般都是取前面的值，如 $sinx$

\sin x = x - \frac{x^{3}}{3!} + \frac{x^{5}}{5!} - \dots + (- 1)^{- 1} \frac{x^{2 m - 1}}{(2 m - 1)!} + R_{2 m}

$\sin x=x-\frac{x^{3}}{3 !}+\frac{x^{5}}{5 !}-\cdots+(-1)^{-1} \frac{x^{2 m-1}}{(2 m-1) !}+R_{2 m}$

当极限x->0时，我们可以将x约等于x，也就是说x已经代表了 $sinx$ ，或者说由于 $sinx$ 满足收敛条件，所以展开式中后面展开项是越来越小的，只用前面的几项就能代表 $sinx$ 函数的大致位置了，同理，由于傅里叶级数需要满足狄利克雷条件，即

给定连个实数级数{ $a_n$ }和{ ${b_n}$ }，如果级数满足

$a_{n} \geq a_{n+1}>0$
$\lim _{n \rightarrow \infty} a_{n}=0$
$\left|\sum_{n=1}^{N} b_{n}\right| \leq M$ ，对于所有正整数 $N$

其中 $M$ 是某个常数，那么级数 $\sum_{n=1}^{\infty} a_{n} b_{n}$ 收敛

也就是说傅里叶级数是递减的，取级数的前几个值，就可以代表周期函数 $f(x)$ ，而后面的几项由于值较小，更多的是更细节上的微调，也就理解了为什么说低频部分代表包络，高频部分代表细节了

还有一个启发就是对基频和谐波的理解，刚开始学习并不知道什么是基频什么是谐波，通过查阅文献傅里叶变换（一）——认识傅里叶变换和结合自己理解有了新的发现

再看傅里叶级数公式

f (t) = A_{0} + \sum_{n = 1}^{\infty} [a_{n} \cos (n ω t) + b_{n} \sin (n ω t)]

$f(t)=A_{0}+\sum_{n=1}^{\infty}\left[a_{n} \cos (n \omega t)+b_{n} \sin (n \omega t)\right]$

而这里的 $A_{0}$ 也就是常说的直流分量，当n=1时，式子变成了

f (t) = A_{0} + a_{1} \cos (ω t) + b_{1} \sin (ω t)

$f(t)=A_{0}+a_{1} \cos (\omega t)+b_{1} \sin (\omega t)$

如果我们使用三角函数的知识，设三角形的两直角边为 $a_1$ 和 $b_1$ ，则斜边为 $\sqrt{a_1^{2}+b_1^{2}}$ ，设为 $c_1$ ，则此时 $sin\alpha=a_1/c_1$ ， $cos\alpha=b_1/c_1$ ，则分别将 $a_1$ 和 $b_1$ 替换则 $f(t)$ 变成

f (t) = A_{0} + c_{1} \sin (φ) \cos (ω t) + c_{1} \cos (φ) \sin (ω t)

$f(t)=A_{0}+c_1 \sin (\varphi) \cos (\omega t)+c_{1} \cos (\varphi) \sin (\omega t)$

化简后可以得到

f (t) = A_{0} + c_{1} \sin (ω t + φ)

$f(t)=A_{0}+c_1 \sin (\omega t + \varphi)$

这里的 $c_1sin(wt + \varphi)$ 称为一次谐波分量，也成为基波，它的频率也就是基频或基准频率，一次代表傅里叶级数中的n取1，谐波也指正弦或余弦波，其中 $\varphi$ 为相位，也说明了傅里叶变换成频谱图只包含了 $c_1$ 到 $c_n$ 的值，并没有包含相位 $\varphi$ ，于是完整的转换还需要一个相位谱，则完整的图如下

2.1.2.4.MFCC原理

上面讲完了需要了解的知识点，就可以讲MFCC原理了，结合梅尔频率倒谱系数（MFCC）学习笔记，很全面很详细的一篇入门文章，我们可以有更深入的理解

首先语音信号是时域函数，与我们通常看到的波形图一致

通过FTT快速傅里叶变换可以将时域转为频域

频域中的主要频率成分也就是共振峰(formants)，而共振峰决定着音色，所以可以用来分辨不同人的声音

而这些共振峰是在频谱图的包络中，所以我们需要将频谱图中的包络提取出来

我们在离散余弦变换中讲过，包络主要是频谱的低频部分，频谱细节是在高频部分，我们只需要将包络提取出来后经过低频滤波器进行滤波就可以得到低频部分，而怎么进行提取呢？用到的就是前面讲的同态信号处理的方法，假设频谱由两部分组成：

包络
频谱细节

假设两者是通过频域中乘积的形式（时域中卷积）

至于为什么要这样假设本人也不知道，有兴趣可以去找找资料

X [k] = H [k] E [k]

$X[k]=H[k] E[k]$

然后取对数，就变成了两者相加，再进行滤波就得到了低频的倒谱(cepstrum)(如有疑惑请参考同态信号处理)

Mel频率分析就是基于人类听觉感知实验的。实验观测发现人耳就像一个滤波器组一样，它只关注某些特定的频率分量（人的听觉对频率是有选择性的）。也就说，它只让某些频率的信号通过，而压根就直接无视它不想感知的某些频率信号。但是这些滤波器在频率坐标轴上却不是统一分布的，在低频区域有很多的滤波器，他们分布比较密集，但在高频区域，滤波器的数目就变得比较少，分布很稀疏

总结就是我们要将普通频谱转换为Mel频谱，于是有了如下公式：

mel (f) = 2595 * \log_{10} (1 + f / 700)

$\operatorname{mel}(f)=2595 * \log _{10}(1+f / 700)$

得到Mel频谱之后，由于我们的目标依然是分离出包络，所以还是继续进行取对数，然后进行反傅里叶变换变成时域信号，这样就得到了Mel倒谱(Mel Frequency Cepstral)

那么什么是Mel倒谱系数(Mel Frequency Cepstral Coefficents)呢？

上面的步骤是对频谱取对数，然后进行反傅里叶变换，如果将反傅里叶变换改成DCT离散余弦变换

c_{i} = \sqrt{\frac{2}{N}} \sum_{j = 1}^{N} m_{j} \cos (\frac{π i}{N} (j - 0.5))

$c_{i}=\sqrt{\frac{2}{N}} \sum_{j=1}^{N} m_{j} \cos \left(\frac{\pi i}{N}(j-0.5)\right)$

如果用反傅里叶变换将频域变为时域，那么为什么不是反DCT变换而是DCT变换呢？

一般取12-16阶系数（HTK中为13），也就是前第1阶到第13阶系数（系数为13）或其他如2到13阶系数（系数为12），即cos前面的系数作为MFCC系数，也就是倒谱系数（Mel Frequency Cepstral Coefficents）

前面还有预加重、分帧、加窗等操作这里就不过多介绍了，请参考声学特征（二） MFCC特征原理

2.1.2.4.bap特征

bap全称为Band Aperiodicity，可以参考语音合成vocoder（四） aperiodicity参数，由于做的的拼接合成，所以此参数没有详细研究

2.1.3.HMM原理

参考文献：

HMM公式推导

前面介绍了我们训练所需要准备的数据的原理，如基频F0，MFCC特征，忽略了前端处理中语言学特征的原理讲解，由于在HTS训练中，使用的是festival进行前端处理，如果大家感兴趣，可以参考festival文档

现在我们的语言学特征、MFCC、F0、Bap特征都已经准备好了，接下来进行训练则需要知道HMM的工作原理

2.1.3.1.什么是HMM

这里推荐一个讲解非常好的博客一文搞懂HMM（隐马尔可夫模型），里面对HMM模型进行了很生动的讲解，所以这里就不再赘述

总结一下就是有三个骰子和一串从这三个骰子中摇出的数字，我们需要知道这串数字中每个数字是由第几个骰子摇出来的，这个时候我们就可以对此问题用HMM进行建模，将3个骰子看作是3个隐状态(state)，将这串数字看作是观测状态(observation)，通过前后向算法、EM算法或Baum-Welch算法、维特比算法来求解最大概率骰子序列

2.1.3.2.HMM推导

标记说明：

I为隐藏状态序列
O为观测状态序列
$i_t$ 为t时刻的隐藏状态
$o_t$ 为t时刻的观测状态
$Q=\{q_1,q_2,...q_N\}$ 为状态集合，即有规定的 $N$ 个隐藏状态，例如一般在英文语音合成中用3个隐藏状态，中文语音合成中用五个隐藏状态
$V=\{v_1,v_2,...v_m\}$ 为观测集合
$a_{ij}$ 为转移概率，即 $P(i_{t+1}=q_j \mid i_t=q_i)$ ，t时刻隐藏状态为 $q_i$ 的概率转移到t+1时刻隐藏状态为 $q_j$ 的概率
$b_j(k)$ 为发射概率，即 $P(o_{t}=v_k \mid i_t=q_j)$ ，t时刻隐藏状态为 $q_j$ 到观测状态为 $v_k$ 的概率

如果初次学习HMM，建议可以先看看b站的HMM讲解和公式推导

推导过程中主要涉及的理论知识有：

概率论中简单的联合概率公式

P (a, b ∣ c) = P (a ∣ b, c) P (b ∣ c)

$P(a, b \mid c)=P(a \mid b, c) P(b \mid c)$

推导：

\begin{aligned} P (a, b ∣ c) & = \frac{P (a, b, c)}{P (b, c)} \frac{P (b, c)}{P (c)} \\ = P (a ∣ b, c) P (b ∣ c) \end{aligned}

$\begin{aligned} P(a, b \mid c) &= \frac{P(a, b, c)}{P(b, c)} \frac{P(b, c)}{P(c)} \\ &= P(a \mid b, c) P(b \mid c) \end{aligned}$

齐次马尔可夫假设：任意时刻t的隐藏状态只依赖于前一时刻的隐藏状态，与其他时刻隐藏状态及观测状态无关（初始时刻的隐藏状态由参数 $\pi$ 决定)

至于为什么初始状态概率为 $\pi$ ，暂时未知，印象中概率论中对马尔可夫链推导中有用到，这里没有进行详细研究

P (i_{t} ∣ i_{t - 1}, i_{t - 2}, \dots, i_{1}, o_{t}, \dots, o_{1}) = P (i_{t} ∣ i_{t - 1}), t = 2, 3, \dots, T

$P\left(i_{t} \mid i_{t-1}, i_{t-2}, \ldots, i_{1}, o_{t}, \ldots, o_{1}\right)=P\left(i_{t} \mid i_{t-1}\right), \quad t=2,3, \ldots, T$

观测独立假设：任意时刻的观测只依赖于该时刻的马尔科夫链状态，与其他观测及状态无关
$P\left(o_{t} \mid o_{t-1}, o_{t-2}, \ldots, o_{1}, i_{t}, i_{t-1}, \ldots, i_{1}\right)=P\left(o_{t} \mid i_{t}\right), \quad t=1,2, \ldots, T$

HMM中三个问题：

概率计算问题：给定模型参数 $\lambda=(\pi, A, B)$ 和观测序列 $O=\left(o_{1}, o_{2}, \ldots, o_{T}\right)$ ，其中 $\pi$ 为初始转移概率， $A$ 为转移矩阵， $B$ 为观测矩阵，计算在模型在 $\lambda$ 下观测序列 $O$ 出现的概率 $P (O∣\lambda)$ （前后向算法）
参数估计问题：仅已知观测序列 $O$ ，求 $λ$ , 即 $P(λ∣O)$ （EM或Baum-Welch算法）
预测问题/推理问题/decoding：已知 $O$ 和 $λ$ ，求 $I=\left(i_{1}, i_{2}, \ldots, i_{T}\right)$ ，即 $P(I∣O,λ)$ （维特比算法）

一、概率计算问题-求解 $P (O∣\lambda)$ 

P (O ∣ λ) = P (o_{t}, o_{t - 1}, o_{t - 2}, \dots, o_{1} ∣ λ)

$P (O∣\lambda)=P(o_{t},o_{t-1}, o_{t-2}, \ldots, o_{1} \mid \lambda)$

假设我们有 $N$ 个隐状态，则根据边缘概率分布和联合概率分布的关系即联和概率之和为边缘概率，可以得出，对隐状态 $I$ （当作联合概率）求和加入到概率中不影响最终概率（边缘概率）

\begin{aligned} P (O ∣ λ) & = \sum_{I} P (O, I ∣ λ) \\ = \sum_{I} P (O ∣ I, λ) P (I ∣ λ) \end{aligned}

$\begin{aligned} P(O \mid \lambda)&=\sum_{I} P(O, I \mid \lambda) \\ &=\sum_{I} P(O \mid I,\lambda)P(I \mid \lambda) \\ \end{aligned}$

这里的 $P(O, I \mid \lambda)$ 可以使用观测独立假设进行变形，即：

\begin{aligned} P (O, I ∣ λ) & = P (o_{t} ∣ o_{t - 1}, o_{t - 2}, \dots, o_{1}, i_{t}, i_{t - 1}, \dots, i_{1} λ) P (o_{t - 1}, o_{t - 2}, \dots, o_{1}, i_{t}, i_{t - 1}, \dots, i_{1} ∣ λ) \\ = P (o_{t} ∣ i_{t}) P (o_{t - 1}, o_{t - 2}, \dots, o_{1}, i_{t}, i_{t - 1}, \dots, i_{1} ∣ λ) \\ = P (o_{t} ∣ i_{t}) P (o_{t - 1} ∣ i_{t - 1}) P (o_{t - 2}, o_{t - 3}, \dots, o_{1}, i_{t}, i_{t - 1}, \dots, i_{1} ∣ λ) \\ = \dots \\ = \prod_{t = 1}^{T} b_{i_{t}} (o_{t}) \end{aligned}

$\begin{aligned} P(O, I \mid \lambda) &= P(o_{t}\mid o_{t-1}, o_{t-2}, \ldots, o_{1},i_{t}, i_{t-1}, \ldots, i_{1} \lambda)P(o_{t-1}, o_{t-2}, \ldots, o_{1},i_{t}, i_{t-1}, \ldots, i_{1} \mid \lambda) \\ &=P\left(o_{t} \mid i_{t}\right)P(o_{t-1}, o_{t-2}, \ldots, o_{1},i_{t}, i_{t-1}, \ldots, i_{1} \mid \lambda) \\ &=P\left(o_{t} \mid i_{t}\right)P\left(o_{t-1} \mid i_{t-1}\right)P(o_{t-2}, o_{t-3}, \ldots, o_{1},i_{t}, i_{t-1}, \ldots, i_{1} \mid \lambda) \\ &= \ldots \\ &= \prod_{t=1}^{T} b_{i_t}(o_t) \end{aligned}$

则 $P(I \mid \lambda)$ 可使用齐次马尔可夫假设进行变形（这里的初始转移概率为 $\pi_{i}$ ，其中 $\sum_{i=1}^{N} \pi_{i} = 1$ ），即：

\begin{aligned} P (I ∣ λ) & = P (i_{t}, i_{t - 1}, \dots, i_{1} ∣ λ) \\ = P (i_{t} ∣ i_{t - 1}, \dots, i_{1}, λ) P (i_{t - 1}, \dots, i_{1} ∣ λ) \\ = P (i_{t} ∣ i_{t - 1}) P (i_{t - 1}, \dots, i_{1} ∣ λ) \\ = . . . \\ = π_{i} \prod_{t = 2}^{T} a_{i_{t - 1} i_{t}} \end{aligned}

$\begin{aligned} P(I \mid \lambda) &= P(i_{t}, i_{t-1}, \ldots, i_{1} \mid \lambda) \\ &=P(i_{t} \mid i_{t-1}, \ldots, i_{1}, \lambda) P(i_{t-1}, \ldots, i_{1}\mid \lambda) \\ &=P\left(i_{t} \mid i_{t-1}\right)P(i_{t-1}, \ldots, i_{1}\mid \lambda) \\ &=... \\ &=\pi_i \prod_{t=2}^{T} a_{i_{t-1}i_t} \end{aligned}$

带入原式子中得：

\begin{aligned} P (O ∣ λ) & = \sum_{I} P (O ∣ I, λ) P (I ∣ λ) \\ = \sum_{I} \prod_{t = 1}^{T} b_{i_{t}} (o_{t}) π_{i} \prod_{t = 2}^{T} a_{i_{t - 1} i_{t}} \end{aligned}

$\begin{aligned} P(O \mid \lambda)&=\sum_{I} P(O \mid I,\lambda)P(I \mid \lambda) \\ &=\sum_{I}\prod_{t=1}^{T} b_{i_t}(o_t) \pi_i \prod_{t=2}^{T} a_{i_{t-1}i_t} \end{aligned}$

由于 $I=\left(i_{1}, i_{2}, \ldots, i_{T}\right)$ ，则将 $\sum_{I}$ 展开得：

P (O ∣ λ) = \sum_{i_{1}} \sum_{i_{2}} . . . \sum_{i_{T}} \prod_{t = 1}^{T} b_{i_{t}} (o_{t}) π_{i} \prod_{t = 2}^{T} a_{i_{t - 1} i_{t}}

$P(O \mid \lambda) = \sum_{i_1} \sum_{i_2}...\sum_{i_T} \prod_{t=1}^{T} b_{i_t}(o_t) \pi_i \prod_{t=2}^{T} a_{i_{t-1}i_t}$

由于 $i$ 的隐藏集合中有 $N$ 种取法，所以此算法的时间复杂度为 $O(N^T)$ ，会随时间出现指数级增长，所以就需要更加高效的算法，即前向算法和后向算法

前向算法的推导

现在定义一个 $\alpha_{t}(k)$ 表示为：

α_{t} (k) = P (o_{1}, o_{2}, \dots, o_{t}, i_{t} = q_{k} ∣ λ)

$\alpha_{t}(k)=P\left(o_{1}, o_{2}, \ldots, o_{t}, i_{t}=q_{k} \mid \lambda\right)$

则所有的 $i_t$ 的取值集合为边缘概率，即：

P (O ∣ λ) = \sum_{k = 1}^{N} P (o_{1}, o_{2}, \dots, o_{t}, i_{t} = q_{k} ∣ λ) = \sum_{k = 1}^{N} α_{t} (k)

$P(O \mid \lambda) = \sum_{k=1}^{N} P\left(o_{1}, o_{2}, \ldots, o_{t}, i_{t}=q_{k} \mid \lambda\right) = \sum_{k=1}^{N} \alpha_{t}(k)$

现在已知 $P(O \mid \lambda)$ 是由 $\alpha_{t}(k)$ 表示，即与时序t相关，我们推导得出递推公式为：

\begin{aligned} α_{t + 1} (j) & = \sum_{k = 1}^{N} P (o_{1}, \dots, o_{t + 1}, i_{t} = q_{k}, i_{t + 1} = q_{j} ∣ λ) \\ == \sum_{k = 1}^{N} P (o_{t + 1} ∣ o_{1}, \dots, o_{t}, i_{t} = q_{k}, i_{t + 1} = q_{j}, λ) P (o_{1}, \dots, o_{t}, i_{t} = q_{k}, i_{t + 1} = q_{j} ∣ λ) \end{aligned}

$\begin{aligned} \alpha_{t+1}(j) &= \sum_{k=1}^{N} P\left(o_{1}, \ldots, o_{t+1}, i_{t}=q_{k}, i_{t+1}=q_{j} \mid \lambda\right) \\ &= =\sum_{k=1}^{N} P\left(o_{t+1} \mid o_{1}, \ldots, o_{t}, i_{t}=q_{k}, i_{t+1}=q_{j}, \lambda\right) P\left(o_{1}, \ldots, o_{t}, i_{t}=q_{k}, i_{t+1}=q_{j} \mid \lambda\right) \end{aligned}$

则根据观测独立假设

P (o_{t + 1} ∣ o_{1}, \dots, o_{t}, i_{t} = q_{k}, i_{t + 1} = q_{j}, λ) = P (o_{t + 1} ∣ i_{t + 1} = q_{j}) = b_{j} (o_{t + 1})

$P\left(o_{t+1} \mid o_{1}, \ldots, o_{t}, i_{t}=q_{k}, i_{t+1}=q_{j}, \lambda\right) = P(o_{t+1} \mid i_{t+1}=q_{j}) = b_j(o_{t+1})$

则根据齐次马尔可夫假设

\begin{aligned} P (o_{1}, \dots, o_{t}, i_{t} = q_{k}, i_{t + 1} = q_{j} ∣ λ) & = P (i_{t + 1} = q_{j} ∣ o_{1}, \dots, o_{t}, i_{t} = q_{k}, λ) \\ P (o_{1}, \dots, o_{t}, i_{t} = q_{k} ∣ λ) \\ = P (i_{t + 1} = q_{j} ∣ i_{t} = q_{k}) P (o_{1}, \dots, o_{t}, i_{t} = q_{k} ∣ λ) \\ = a_{k j} α_{t} (k) \end{aligned}

$\begin{aligned} P\left(o_{1}, \ldots, o_{t}, i_{t} = q_{k}, i_{t+1}=q_{j} \mid \lambda\right) &= P\left( i_{t+1}=q_{j} \mid o_{1}, \ldots, o_{t}, i_{t}=q_{k},\lambda\right) \\ P\left(o_{1}, \ldots, o_{t}, i_{t}=q_{k} \mid \lambda\right) \\ &= P\left( i_{t+1}=q_{j} \mid i_{t}=q_{k}\right)P\left(o_{1}, \ldots, o_{t}, i_{t}=q_{k} \mid \lambda\right) \\ &= a_{kj}\alpha_{t}(k) \end{aligned}$

带入 $\alpha_{t+1}(j)$ 式子得：

α_{t + 1} (j) = \sum_{k = 1}^{N} b_{j} (o_{t + 1}) a_{k j} α_{t} (k)

$\alpha_{t+1}(j)=\sum_{k=1}^{N}b_j(o_{t+1})a_{kj}\alpha_{t}(k)$

分析时间复杂度，由于时间序列为T，首先需要循环时间T，然后有两个嵌套得对 $N$ 求和，则总共得时间复杂度为 $O(TN^2)$

后向算法的推导

前向算法是从前向后推导，则后向算法是从后向前推导

首先定义后向概率（即上图中的紫色部分）

\begin{aligned} β_{t} (k) & = P (o_{t + 1}, \dots, o_{T} ∣ i_{t} = q_{k}, λ) \\ . . . \\ β_{1} (k) & = P (o_{2}, \dots, o_{T} ∣ i_{1} = q_{k}, λ) \end{aligned}

$\begin{aligned} \beta_{t}(k) &= P\left(o_{t+1},\ldots, o_{T} \mid i_{t}=q_{k}, \lambda\right) \\ ... \\ \beta_{1}(k) &= P\left(o_{2},\ldots, o_{T} \mid i_{1}=q_{k}, \lambda\right) \end{aligned}$

则在模型已知情况下观测概率 $P(O \mid \lambda)$ 与 $\beta_{t}(k)$ 之间的关系

\begin{aligned} P (O ∣ λ) & = P (o_{1}, \dots, o_{T} ∣ λ) \\ = \sum_{k = 1}^{N} P (o_{1}, \dots, o_{T}, i_{1} = q_{k} ∣ λ) \\ = \sum_{k = 1}^{N} P (o_{1}, \dots, o_{T}, ∣ i_{1} = q_{k}, λ) P (i_{1} = q_{k} ∣ λ) \\ = \sum_{k = 1}^{N} P (o_{1}, ∣ o_{2}, \dots, o_{T}, i_{1} = q_{k}, λ) P (o_{2}, \dots, o_{T} ∣ i_{1} = q_{k}, λ) \\ = \sum_{k = 1}^{N} P (o_{1}, ∣ i_{1} = q_{k}, λ) P (i_{1} = q_{k} ∣ λ) P (o_{2}, \dots, o_{T} ∣ i_{1} = q_{k}, λ) \\ = \sum_{k = 1}^{N} b_{k} (o_{1}) π_{k} β_{1} (k) \end{aligned}

$\begin{aligned} P(O \mid \lambda) &= P\left(o_{1},\ldots, o_{T} \mid \lambda\right) \\ &= \sum_{k=1}^{N} P\left(o_{1},\ldots, o_{T},i_{1}=q_{k} \mid \lambda\right) \\ &= \sum_{k=1}^{N} P\left(o_{1},\ldots, o_{T}, \mid i_{1}=q_{k}, \lambda\right) P\left(i_{1}=q_{k} \mid\lambda\right) \\ &= \sum_{k=1}^{N} P\left(o_{1}, \mid o_{2}, \ldots,o_{T}, i_{1}=q_{k}, \lambda\right) P\left(o_{2}, \ldots,o_{T} \mid i_{1}=q_{k}, \lambda\right) \\ &= \sum_{k=1}^{N} P\left(o_{1}, \mid i_{1}=q_{k}, \lambda\right) P\left(i_{1}=q_{k} \mid\lambda\right) P\left(o_{2}, \ldots,o_{T} \mid i_{1}=q_{k}, \lambda\right) \\ &= \sum_{k=1}^{N} b_k(o_1) \pi_k \beta_{1}(k) \end{aligned}$

由于

\begin{aligned} α_{1} (k) & = P (o_{1}, i_{1} = q_{k} ∣ λ) \\ = P (i_{1} = q_{k} ∣ λ) P (o_{1} ∣ i_{1} = q_{k}, λ) \\ = π_{k} b_{k} (o_{1}) \end{aligned}

$\begin{aligned} \alpha_{1}(k) &= P\left(o_{1}, i_{1}=q_{k} \mid \lambda\right) \\ &= P\left(i_{1}=q_{k} \mid \lambda\right) P\left(o_{1} \mid i_{1}=q_{k}, \lambda\right) \\ &= \pi_{k}b_{k}(o_1) \end{aligned}$

即对于t=1的时序

P (O ∣ λ) = \sum_{k = 1}^{N} α_{1} (k) β_{1} (k)

$P(O \mid \lambda) = \sum_{k=1}^{N} \alpha_{1}(k)\beta_{1}(k)$

则时序为t时

P (O ∣ λ) = \sum_{k = 1}^{N} α_{t} (k) β_{t} (k)

$P(O \mid \lambda) = \sum_{k=1}^{N} \alpha_{t}(k)\beta_{t}(k)$

现在假设已知所有的 $\beta_{t+1}(k)$ 来推导 $\beta_{t}(k)$

\begin{aligned} β_{t} (k) & = P (o_{t + 1}, \dots, o_{T} ∣ i_{t} = q_{k}, λ) \\ = \sum_{j = 1}^{N} P (o_{t + 1}, \dots, o_{T}, i_{t + 1} = q_{j} ∣ i_{t} = q_{k}, λ) \\ = \sum_{j = 1}^{N} P (o_{t + 1}, \dots, o_{T} ∣ i_{t + 1} = q_{j}, i_{t} = q_{k}, λ) P (i_{t + 1} = q_{j} ∣ i_{t} = q_{k}, λ) \\ = \sum_{j = 1}^{N} P (o_{t + 1} ∣ o_{t + 2}, \dots, o_{T}, i_{t + 1} = q_{j}, i_{t} = q_{k}, λ) P (o_{t + 2}, \dots, o_{T} ∣ i_{t + 1} = q_{j}, i_{t} = q_{k}, λ) a_{k j} \\ = \sum_{j = 1}^{N} P (o_{t + 1} ∣ i_{t + 1} = q_{j}, λ) P (o_{t + 2}, \dots, o_{T} ∣ i_{t + 1} = q_{j}, i_{t} = q_{k}, λ) a_{k j} \\ = \sum_{j = 1}^{N} b_{j} (o_{t + 1}) β_{t + 1} (j) a_{k j} \end{aligned}

$\begin{aligned} \beta_{t}(k) &= P\left(o_{t+1},\ldots, o_{T} \mid i_{t}=q_{k}, \lambda\right) \\ &= \sum_{j=1}^{N} P\left(o_{t+1},\ldots, o_{T},i_{t+1}=q_{j} \mid i_{t}=q_{k}, \lambda\right) \\ &= \sum_{j=1}^{N} P\left(o_{t+1},\ldots, o_{T} \mid i_{t+1}=q_{j},i_{t}=q_{k}, \lambda\right) P\left(i_{t+1}=q_{j} \mid i_{t}=q_{k}, \lambda\right) \\ &= \sum_{j=1}^{N} P\left(o_{t+1} \mid o_{t+2} ,\ldots, o_{T}, i_{t+1}=q_{j},i_{t}=q_{k},\lambda\right) P\left(o_{t+2} ,\ldots, o_{T} \mid i_{t+1}=q_{j},i_{t}=q_{k},\lambda\right) a_{kj} \\ &= \sum_{j=1}^{N} P\left(o_{t+1} \mid i_{t+1}=q_{j},\lambda\right) P\left(o_{t+2} ,\ldots, o_{T} \mid i_{t+1}=q_{j},i_{t}=q_{k},\lambda\right) a_{kj} \\ &= \sum_{j=1}^{N} b_{j}(o_{t+1}) \beta_{t+1}(j) a_{kj} \end{aligned}$

则推导出 $\beta_{t+1}(k)$ 和 $\beta_{t}(k)$ 之间的关系，如果要求 $P(O \mid \lambda)$ 需要求出 $\beta_{1}(k)$ ，而 $\beta_{1}(k)$ 需要求 $\beta_{2}(k)$ ，依次类推，直到第T个，也就是T*N次循环可以求出 $\beta_{1}(k)$ ，再进行N次循环则得出 $P(O \mid \lambda)$ ，所以时间复杂度为 $O(TN^2)$

二、参数估计问题-求解 $P(λ∣O)$

本节首先需要用到EM相关知识，请结合EM算法原理及推导和（EM算法）The EM Algorithm两篇文章来进行理解和推导，第二篇文章还讲到了在混合高斯模型(GMM)中的应用，GMM也是在后面需要与HMM一起介绍的知识，这里就不重点介绍了

总结上面两篇文章，EM算法主要是用来解决存在隐变量，通过极大似然估计无法直接求参数的问题，通过E步骤先确定一个变量值，然后通过M求出最大化参数，再通过E步骤确定新的变量值，通过M步骤最大化参数，如此循环，直到参数变换小于某个设定的阈值为止，此外，也有使用Baum-Welch算法的，原理EM相同，在此就不做过多介绍

EM中E步骤中的Q函数表达式如下

Q (λ, λ^{(t)}) = E_{I ∣ O, λ^{(t)}} [\log P (O, I ∣ λ)] = \sum_{I} \log P (O, I ∣ λ) P (I ∣ O, λ^{(t)})

$Q\left(\lambda, \lambda^{(t)}\right)=\mathbb{E}_{I \mid O, \lambda^{(t)}}[\log P(O, I \mid \lambda)] =\sum_{I} \log P(O, I \mid \lambda) P\left(I \mid O, \lambda^{(t)}\right)$

通过前面的推导我们得到了如下公式

\begin{aligned} P (O ∣ λ) & = \sum_{I} \prod_{t = 1}^{T} b_{i_{t}} (o_{t}) π_{i} \prod_{t = 2}^{T} a_{i_{t - 1} i_{t}} \end{aligned}

$\begin{aligned} P(O \mid \lambda)&= \sum_{I}\prod_{t=1}^{T} b_{i_t}(o_t) \pi_i \prod_{t=2}^{T} a_{i_{t-1}i_t} \\ \end{aligned}$

可推出

P (O, I ∣ λ) = \prod_{t = 1}^{T} b_{i_{t}} (o_{t}) π_{i} \prod_{t = 2}^{T} a_{i_{t - 1} i_{t}}

$P(O,I\mid \lambda)= \prod_{t=1}^{T} b_{i_t}(o_t) \pi_i \prod_{t=2}^{T} a_{i_{t-1}i_t}$

将式子带入Q得

Q (λ, λ^{(t)}) = \sum_{I} [\log (π_{i}) + l o g (\prod_{t = 2}^{T} a_{i_{t - 1} i_{t}}) + l o g (\prod_{t = 1}^{T} b_{i_{t}} (o_{t}))] P (I ∣ O, λ^{(t)})

$Q\left(\lambda, \lambda^{(t)}\right) =\sum_{I} {[}{\log (\pi_i)} + {log(\prod_{t=2}^{T} a_{i_{t-1}i_t})} + {log(\prod_{t=1}^{T} b_{i_t}(o_t))} {]} P\left(I \mid O, \lambda^{(t)}\right)$

由模型得参数为 $\lambda=(\pi, A, B)$ ， $\pi$ 为初始转移概率，与公式中蓝色部分对应（左）， $A$ 为转移矩阵，与公式中红色部分对应（中）， $B$ 为观测矩阵，与公式中绿色部分对应（右），这里以 $\pi$ 的更新为例进行讲解， $A$ 与 $B$ 的更新雷同，可以参考文章HMM公式推导

求解 $\pi^{t+1}$

\begin{aligned} π^{(t + 1)} & = \arg max_{π} \sum_{I} \log π_{i_{1}} P (O, I ∣ λ^{(t)}) \\ = \arg max_{π} \sum_{i_{1}} \sum_{i_{2}} \dots \sum_{i_{T}} \log π_{i_{1}} P (I, O ∣ λ^{(t)}) \\ = \arg max_{π} \sum_{i_{1}} \log π_{i_{1}} \sum_{i_{2}} \dots \sum_{i_{T}} P (i_{T}, \dots, i_{2}, i_{1}, O ∣ λ^{(t)}) \\ = \arg max_{π} \sum_{i_{1}} \log π_{i_{1}} P (i_{1}, O ∣ λ^{(t)}) \end{aligned}

$\begin{aligned} \pi^{(t+1)}&=\arg \max _{\pi} \sum_{I} \log \pi_{i_1} P\left(O, I \mid \lambda^{(t)}\right) \\ &= \arg \max _{\pi}\sum_{i_{1}} \sum_{i_{2}} \ldots \sum_{i_{T}} \log \pi_{i_1} P\left(I, O \mid \lambda^{(t)}\right) \\ &= \arg \max _{\pi}\sum_{i_{1}} \log \pi_{i_1} \sum_{i_{2}} \cdots \sum_{i_{T}} P\left(i_{T}, \ldots, i_{2}, i_{1}, O \mid \lambda^{(t)}\right) \\ &= \arg \max _{\pi}\sum_{i_{1}} \log \pi_{i_1} P\left(i_{1}, O \mid \lambda^{(t)}\right) \end{aligned}$

上面式子中 $\pi_i$ 由于是初始概率，所以只跟第一个时刻有关，标记改为了 $\pi_{i_1}$ ，从第三步到第四部的变化，是根据边缘概率和联合概率的关系，由于有约束条件 $\sum_{i=1}^{N} \pi_{i} = 1$ ），这里成了条件极值的求解

假设 $f(x,y)$ 满条件 $g(x,y)$ ，求解 $f(x,y)$ 极值的一般步骤为：

构建 $F(x,y)=f(x.y)+ \lambda g(x,y)$
分别求出 $F(x,y)$ 对 $x,y,\lambda$ 的偏微分
令偏微分等于0，由于有三个未知数 $x,y,\lambda$ ，三个偏微分方程，可求出满足条件的 $x, y ,\lambda$

此方法成为拉格朗日乘子法

所以我们要求解的是 $\pi^{(t+1)}$ 在满足 $\sum_{i=1}^{N} \pi_{i} = 1$ 条件下的极值，设构造函数为

L (π, η) = \sum_{i = 1}^{N} \log π_{i_{1}} P (i_{1} = q_{i}, O ∣ λ^{(t)}) + η (\sum_{i = 1}^{N} π_{i} - 1)

$L ( \pi , \eta ) = \sum_{i={1}}^{N} \log \pi_{i_1} P\left(i_{1} = q_{i}, O \mid \lambda^{(t)}\right) + \eta ( \sum _ { i = 1 } ^ { N } \pi _ { i } - 1 )$

对 $\pi_{i}$ 求导，得：

\frac{\partial L (π, η)}{π_{i}} = \frac{1}{π_{i}} P (i_{1} = q_{i}, O | λ^{(t)}) + η

$\frac { \partial L ( \pi , \eta ) } { \pi _ { i } } = \frac { 1 } { \pi _ { i } } P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } ) + \eta$

令导数为0得：

P (i_{1} = q_{i}, O | λ^{(t)}) + η π_{i} = 0

$P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } ) + \eta \pi_{i} = 0$

由于 $\pi_{i}$ 在初始状态下的 $N$ 个值都是一样的，所以对 $N$ 个分量求导的结果都是一样的，将这些结果相加得到：

\sum_{i = 1}^{N} P (i_{1} = q_{i}, O | λ^{(t)}) + \sum_{i = 1}^{n} η π_{i} = 0

$\sum _ { i = 1 } ^ { N } P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } ) + \sum _ { i = 1 } ^ { n }\eta \pi _ { i } = 0$

再利用 $\pi$ 的约束条件，就得到：

η \sum_{i = 1}^{N} π_{i} = - \sum_{i = 1}^{N} P (i_{1} = q_{i}, O | λ^{(t)})

$\eta \sum _ { i = 1 } ^ { N } \pi _ { i } = - \sum _ { i = 1 } ^ { N } P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } )$

η = - \sum_{i = 1}^{N} P (i_{1} = q_{i}, O | λ^{(t)})

$\eta = - \sum _ { i = 1 } ^ { N } P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } )$

由于是对 $i$ 的积分，所以由边缘概率和联合概率关系则可以去掉 $i_{1}$

η = - P (O | λ^{(t)})

$\eta = - P ( O | \lambda^{ ( t ) } )$

则此时， $\pi$ 为极大值，带入对 $\pi_{i}$ 求导后令其为0的表达式得：

\begin{aligned} π^{(t + 1)} & = \arg max_{π} \sum_{i_{1}} \log π_{i_{1}} P (i_{1}, O ∣ λ^{(t)}) \\ = - \frac{P (i_{1} = q_{i}, O | λ^{(t)})}{η} \\ = \frac{P (i_{1} = q_{i}, O | λ^{(t)})}{P (O | λ^{(t)})} \end{aligned}

$\begin{aligned} \pi^{(t+1)} &= \arg \max _{\pi}\sum_{i_{1}} \log \pi_{i_1} P\left(i_{1}, O \mid \lambda^{(t)}\right) \\ &= - \frac { P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } )}{\eta} \\ &= \frac { P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } )}{ P ( O | \lambda^{ ( t ) } )} \\ \end{aligned}$

由前向概率 $\alpha_{t}(k)$ 定义和 $\beta_{t}(k)$ 定义可得到

\begin{aligned} α_{t} (k) β_{t} (k) & = P (o_{1}, o_{2}, \dots, o_{t}, i_{t} = q_{k} ∣ λ) P (o_{t + 1}, \dots, o_{T} ∣ i_{t} = q_{k}, λ) \\ = P (o_{1}, o_{2}, \dots, o_{t} ∣ i_{t} = q_{k}, λ) P (i_{t} = q_{k} ∣ λ) P (o_{t + 1}, \dots, o_{T} ∣ i_{t} = q_{k}, λ) \\ = P (o_{1}, o_{2}, \dots, o_{t}, o_{t + 1}, \dots, o_{T} ∣ i_{t} = q_{k}, λ) P (i_{t} = q_{k} ∣ λ) \\ = P (o_{1}, o_{2}, \dots, o_{t}, o_{t + 1}, \dots, o_{T}, i_{t} = q_{k} ∣ λ) \\ = P (i_{t} = q_{i}, O | λ^{(t)}) \end{aligned}

$\begin{aligned} \alpha_{t}(k)\beta_{t}(k) &= P\left(o_{1}, o_{2}, \ldots, o_{t}, i_{t}=q_{k} \mid \lambda\right) P\left(o_{t+1},\ldots, o_{T} \mid i_{t}=q_{k}, \lambda\right) \\ &= P\left(o_{1}, o_{2}, \ldots, o_{t} \mid i_{t}=q_{k}, \lambda\right) P\left(i_{t}=q_{k} \mid \lambda\right) P\left(o_{t+1},\ldots, o_{T} \mid i_{t}=q_{k}, \lambda\right) \\ &= P\left(o_{1}, o_{2}, \ldots, o_{t},o_{t+1},\ldots, o_{T} \mid i_{t}=q_{k}, \lambda\right) P\left(i_{t}=q_{k} \mid \lambda\right) \\ &= P\left(o_{1}, o_{2}, \ldots, o_{t},o_{t+1},\ldots, o_{T},i_{t}=q_{k} \mid \lambda\right) \\ &= P ( i _ { t } = q _ { i } , O | \lambda^{ ( t ) } ) \end{aligned}$

即得到 $\pi^{(t+1)}$ 的表达式

\begin{aligned} π^{(t + 1)} & = \frac{P (i_{1} = q_{i}, O | λ^{(t)})}{\sum_{i = 1}^{N} P (O | λ^{(t)})} \\ = \frac{α_{t} (k) β_{t} (k)}{P (O | λ^{(t)})} \end{aligned}

$\begin{aligned} \pi^{(t+1)} &= \frac { P ( i _ { 1 } = q _ { i } , O | \lambda^{ ( t ) } )}{\sum _ { i = 1 } ^ { N } P ( O | \lambda^{ ( t ) } )} \\ &= \frac{\alpha_{t}(k)\beta_{t}(k)}{P ( O | \lambda^{ ( t ) } )} \end{aligned}$

自此，就得到了 $\pi^{t+1}$ 的表达式，由于转移矩阵 $A$ 和观测矩阵 $B$ 的表达式与此类似，就不再进行推导，详细的推导可以见HMM公式推导