在任意一个Automatic speech recognition 系统中，第一步就是提取特征。换句话说，我们需要把音频信号中具有辨识性的成分提取出来，然后把其他的乱七八糟的信息扔掉，例如背景噪声啊，情绪啊等等。

      搞清语音是怎么产生的对于我们理解语音有很大帮助。人通过声道产生声音，声道的shape（形状？）决定了发出怎样的声音。声道的shape包括舌头，牙齿等。如果我们可以准确的知道这个形状，那么我们就可以对产生的音素phoneme进行准确的描述。声道的形状在语音短时功率谱的包络中显示出来。而MFCCs就是一种准确描述这个包络的一种特征。

       MFCCs（Mel Frequency Cepstral Coefficents）是一种在自动语音和说话人识别中广泛使用的特征。它是在1980年由Davis和Mermelstein搞出来的。从那时起。在语音识别领域，MFCCs在人工特征方面可谓是鹤立鸡群，一枝独秀，从未被超越啊（至于说Deep Learning的特征学习那是后话了）。

       好，到这里，我们提到了一个很重要的关键词：声道的形状，然后知道它很重要，还知道它可以在语音短时功率谱的包络中显示出来。哎，那什么是功率谱？什么是包络？什么是MFCCs？它为什么有效？如何得到？下面咱们慢慢道来。

 

一、声谱图（Spectrogram）

         我们处理的是语音信号，那么如何去描述它很重要。因为不同的描述方式放映它不同的信息。那怎样的描述方式才利于我们观测，利于我们理解呢？这里我们先来了解一个叫声谱图的东西。

       这里，这段语音被分为很多帧，每帧语音都对应于一个频谱（通过短时FFT计算），频谱表示频率与能量的关系。在实际使用中，频谱图有三种，即线性振幅谱、对数振幅谱、自功率谱（对数振幅谱中各谱线的振幅都作了对数计算，所以其纵坐标的单位是dB（分贝）。这个变换的目的是使那些振幅较低的成分相对高振幅成分得以拉高，以便观察掩盖在低幅噪声中的周期信号）。

  

       我们先将其中一帧语音的频谱通过坐标表示出来，如上图左。现在我们将左边的频谱旋转90度。得到中间的图。然后把这些幅度映射到一个灰度级表示（也可以理解为将连续的幅度量化为256个量化值？），0表示黑，255表示白色。幅度值越大，相应的区域越黑。这样就得到了最右边的图。那为什么要这样呢？为的是增加时间这个维度，这样就可以显示一段语音而不是一帧语音的频谱，而且可以直观的看到静态和动态的信息。优点稍后呈上。

         这样我们会得到一个随着时间变化的频谱图，这个就是描述语音信号的spectrogram声谱图。

      下图是一段语音的声谱图，很黑的地方就是频谱图中的峰值（共振峰formants）。

      那我们为什么要在声谱图中表示语音呢？

      首先，音素（Phones）的属性可以更好的在这里面观察出来。另外，通过观察共振峰和它们的转变可以更好的识别声音。隐马尔科夫模型（Hidden Markov Models）就是隐含地对声谱图进行建模以达到好的识别性能。还有一个作用就是它可以直观的评估TTS系统（text to speech）的好坏，直接对比合成的语音和自然的语音声谱图的匹配度即可。

通过对语音进行分帧进行时频变换，得到每一帧的FFT频谱再将各帧频谱按照时间顺序排列起来，得到时间-频率-能量分布图。很直观的表现出语音信号随时间的频率中心的变化。

 

 

二、倒谱分析（Cepstrum Analysis）

       下面是一个语音的频谱图。峰值就表示语音的主要频率成分，我们把这些峰值称为共振峰（formants），而共振峰就是携带了声音的辨识属性（就是个人身份证一样）。所以它特别重要。用它就可以识别不同的声音。

        既然它那么重要，那我们就是需要把它提取出来！我们要提取的不仅仅是共振峰的位置，还得提取它们转变的过程。所以我们提取的是频谱的包络（Spectral Envelope）。这包络就是一条连接这些共振峰点的平滑曲线。

      我们可以这么理解，将原始的频谱由两部分组成：包络和频谱的细节。这里用到的是对数频谱，所以单位是dB。那现在我们需要把这两部分分离开，这样我们就可以得到包络了。

      那怎么把他们分离开呢？也就是，怎么在给定log X[k]的基础上，求得log H[k] 和 log E[k]以满足log X[k] = log H[k] + log E[k]呢？

      为了达到这个目标，我们需要Play a Mathematical Trick。这个Trick是什么呢？就是对频谱做FFT。在频谱上做傅里叶变换就相当于逆傅里叶变换Inverse FFT (IFFT)。需要注意的一点是，我们是在频谱的对数域上面处理的，这也属于Trick的一部分。这时候，在对数频谱上面做IFFT就相当于在一个伪频率（pseudo-frequency）坐标轴上面描述信号。

         由上面这个图我们可以看到，包络是主要是低频成分（这时候需要转变思维，这时候的横轴就不要看成是频率了，咱们可以看成时间），我们把它看成是一个每秒4个周期的正弦信号。这样我们在伪坐标轴上面的4Hz的地方给它一个峰值。而频谱的细节部分主要是高频。我们把它看成是一个每秒100个周期的正弦信号。这样我们在伪坐标轴上面的100Hz的地方给它一个峰值。

         把它俩叠加起来就是原来的频谱信号了。

      在实际中咱们已经知道log X[k]，所以我们可以得到了x[k]。那么由图可以知道，h[k]是x[k]的低频部分，那么我们将x[k]通过一个低通滤波器就可以得到h[k]了！没错，到这里咱们就可以将它们分离开了，得到了我们想要的h[k]，也就是频谱的包络。

       x[k]实际上就是倒谱Cepstrum（这个是一个新造出来的词，把频谱的单词spectrum的前面四个字母顺序倒过来就是倒谱的单词了）。而我们所关心的h[k]就是倒谱的低频部分。h[k]描述了频谱的包络，它在语音识别中被广泛用于描述特征。

      那现在总结下倒谱分析，它实际上是这样一个过程：

1）将原语音信号经过傅里叶变换得到频谱：X[k]=H[k]E[k]；

　　只考虑幅度就是：|X[k] |=|H[k]||E[k] |；

2）我们在两边取对数：log||X[k] ||= log ||H[k] ||+ log ||E[k] ||。

3）再在两边取逆傅里叶变换得到：x[k]=h[k]+e[k]。

       这实际上有个专业的名字叫做同态信号处理。它的目的是将非线性问题转化为线性问题的处理方法。对应上面，原来的语音信号实际上是一个卷性信号（声道相当于一个线性时不变系统，声音的产生可以理解为一个激励通过这个系统），第一步通过卷积将其变成了乘性信号（时域的卷积相当于频域的乘积）。第二步通过取对数将乘性信号转化为加性信号，第三步进行逆变换，使其恢复为卷性信号。这时候，虽然前后均是时域序列，但它们所处的离散时域显然不同，所以后者称为倒谱频域。

      总结下，倒谱（cepstrum）就是一种信号的傅里叶变换经对数运算后再进行傅里叶反变换得到的谱。它的计算过程如下：

以下部分还未整理 

---

 

 

三、Mel频率分析（Mel-Frequency Analysis） 

         好了，到这里，我们先看看我们刚才做了什么？给我们一段语音，我们可以得到了它的频谱包络（连接所有共振峰值点的平滑曲线）了。但是，对于人类听觉感知的实验表明，人类听觉的感知只聚焦在某些特定的区域，而不是整个频谱包络。

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

       人的听觉系统是一个特殊的非线性系统，它响应不同频率信号的灵敏度是不同的。在语音特征的提取上，人类听觉系统做得非常好，它不仅能提取出语义信息, 而且能提取出说话人的个人特征，这些都是现有的语音识别系统所望尘莫及的。如果在语音识别系统中能模拟人类听觉感知处理特点，就有可能提高语音的识别率。

        梅尔频率倒谱系数（Mel Frequency Cepstrum Coefficient, MFCC）考虑到了人类的听觉特征，先将线性频谱映射到基于听觉感知的Mel非线性频谱中，然后转换到倒谱上。

        将普通频率转化到Mel频率的公式是：
 
      由下图可以看到，它可以将不统一的频率转化为统一的频率，也就是统一的滤波器组。

      在Mel频域内，人对音调的感知度为线性关系。举例来说，如果两段语音的Mel频率相差两倍，则人耳听起来两者的音调也相差两倍。

 

四、Mel频率倒谱系数（Mel-Frequency Cepstral Coefficients） 

       我们将频谱通过一组Mel滤波器就得到Mel频谱。公式表述就是：log X[k] = log (Mel-Spectrum)。这时候我们在log X[k]上进行倒谱分析：

1）取对数：log X[k] = log H[k] + log E[k]。

2）进行逆变换：x[k] = h[k] + e[k]。

      在Mel频谱上面获得的倒谱系数h[k]就称为Mel频率倒谱系数，简称MFCC。

         现在咱们来总结下提取MFCC特征的过程：（具体的数学过程网上太多了，这里就不想贴了）

1）先对语音进行预加重、分帧和加窗；（加强语音信号性能（信噪比，处理精度等）的一些预处理）

2）对每一个短时分析窗，通过FFT得到对应的频谱；（获得分布在时间轴上不同时间窗内的频谱）

3）将上面的频谱通过Mel滤波器组得到Mel频谱；（通过Mel频谱，将线形的自然频谱转换为体现人类听觉特性的Mel频谱）

4）在Mel频谱上面进行倒谱分析（取对数，做逆变换，实际逆变换一般是通过DCT离散余弦变换来实现，取DCT后的第2个到第13个系数作为MFCC系数），获得Mel频率倒谱系数MFCC，这个MFCC就是这帧语音的特征；（倒谱分析，获得MFCC作为语音特征）

       这时候，语音就可以通过一系列的倒谱向量来描述了，每个向量就是每帧的MFCC特征向量。

     这样就可以通过这些倒谱向量对语音分类器进行训练和识别了。

 

function varargout = GUI(varargin)
% GUI MATLAB code for GUI.fig
%      GUI, by itself, creates a new GUI or raises the existing
%      singleton*.
%
%      H = GUI returns the handle to a new GUI or the handle to
%      the existing singleton*.
%
%      GUI('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in GUI.M with the given input arguments.
%
%      GUI('Property','Value',...) creates a new GUI or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before GUI_OpeningFcn gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to GUI_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
%      instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help GUI

% Last Modified by GUIDE v2.5 15-Mar-2020 17:37:45

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @GUI_OpeningFcn, ...
                   'gui_OutputFcn',  @GUI_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin && ischar(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT


% --- Executes just before GUI is made visible.
function GUI_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to GUI (see VARARGIN)

% Choose default command line output for GUI
handles.output = hObject;

% Update handles structure
guidata(hObject, handles);

% UIWAIT makes GUI wait for user response (see UIRESUME)
% uiwait(handles.figure1);


% --- Outputs from this function are returned to the command line.
function varargout = GUI_OutputFcn(hObject, eventdata, handles) 
% varargout  cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure
varargout{1} = handles.output;


% --- Executes on button press in pushbutton1.
function pushbutton1_Callback(hObject, eventdata, handles)
fprintf('\n识别中...\n\n');
%加载训练好的GMM模型
load speakerData;
load speakerGmm;

waveDir='trainning\';  %导入测试集
Test_speakerData = dir(waveDir);   %获取测试集中的结构体数据，这是一个char类型的结构体
Test_speakerData(1:2) = [];
Test_speakerNum=length(Test_speakerData);
Test_speakerNum
count=0;
%%%%%%%%%%%%%%%%for i=1:Test_speakerNum
%%%读取语音
[filename,filepath]=uigetfile('*.wav','选择音频文件');
set(handles.text1,'string',filepath)
   filep=strcat(filepath,filename); 
[testing_data, fs]=wavread(filep);
sound(testing_data, fs);
save testing_data
load testing_data
y=testing_data
axes(handles.axes1)
plot(y);
xlabel('t');ylabel('幅值');
title('时域图');
%频域
%幅频图
N=length(y);   
fs1=100;  %采样频率
n=0:N-1;
t=n/fs;   %时间序列
yfft =fft(y,N);
mag=abs(yfft);     %取振幅的绝对值
f=n*fs/N;    %频率序列
axes(handles.axes2)
plot(f(1:N/2),mag(1:N/2)); %绘出Nyquist频率之前随频率变化的振幅
xlabel('频率/Hz');
ylabel('振幅');title('频域图');

%相谱
A=abs(yfft);
ph=2*angle(yfft(1:N/2));
ph=ph*180/pi;
axes(handles.axes3);
plot(f(1:N/2),ph(1:N/2));
xlabel('频率/hz'),ylabel('相角'),title('数字0-9的相位谱');
% 绘制功率谱
Fs=1000;
n=0:1/Fs:1;
xn=y;
nfft=1024;
window=boxcar(length(n)); %矩形窗
noverlap=0; %数据无重叠
p=0.9; %置信概率
[Pxx,Pxxc]=psd(xn,nfft,Fs,window,noverlap,p);
index=0:round(nfft/2-1);
k=index*Fs/nfft;
plot_Pxx=10*log10(Pxx(index+1));
plot_Pxxc=10*log10(Pxxc(index+1));
axes(handles.axes4)
plot(k,plot_Pxx);
title('数字0-9的功率谱');

axes(handles.axes5)
surf( speakerData(1).mfcc);  %绘制MFCC的三维图
title('第一个人语音的三维MFCC');   %第一个人说话的mfcc的特征   mfcc是指梅尔倒谱系数

%绘制第一个人的MFCC的全部二维图
axes(handles.axes6)

plot(speakerData(1).mfcc);
title('第一个人语音的全部二维MFCC');

%为了观看效果，绘制一个特征吧
axes(handles.axes7)
plot(speakerData(1).mfcc(:,2));
title('第一个人语音中第一个数字0的二维MFCC');


% hObject    handle to pushbutton1 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)


% --- Executes on button press in pushbutton2.
function pushbutton2_Callback(hObject, eventdata, handles)
load testing_data
load speakerGmm
match= MFCC_feature_compare(testing_data,speakerGmm);
axes(handles.axes8)
plot(match);
hold on
title('匹配结果展示');
[max_1 index]=max(match);
fprintf('\n %s',match);
fprintf('\n\n\n说话人是%s。',speakerData(index).name(1:end-4))
set(handles.text3,'string',speakerData(index).name(1:end-4))
if index==i
  count=count+1;
end


% hObject    handle to pushbutton2 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)


% --- Executes on button press in pushbutton3.
function pushbutton3_Callback(hObject, eventdata, handles)
waveDir = uigetdir(strcat(matlabroot,'\work'), '选择训练音频库' );
set(handles.text8,'string','样本库选择完毕！')
speakerData = dir(waveDir);  % 获取测试集中的结构体数据，这是一个char类型的结构体
speakerData(1:2) = [];
speakerNum=length(speakerData);%测试人数；
set(handles.text6,'string',speakerNum)
save speakerNum

% hObject    handle to pushbutton3 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)


% --- Executes on button press in pushbutton4.
function pushbutton4_Callback(hObject, eventdata, handles)
load speakerNum
set(handles.text8,'string','音频文件特征提取进行中。。。')
pause(1) %停留一秒
for i=1:speakerNum
	fprintf('\n正在提取第%d个人%s的特征', i, speakerData(i,1).name(1:end-4));
       
    [y, fs]=wavread(['trainning\' speakerData(i,1).name]);  %读取样本的音频数据
    y=double(y);    %类型转换
    y=y/max(y);
    epInSampleIndex = epdByVol(y, fs);		% 端点检测
    y=y(epInSampleIndex(1):epInSampleIndex(2));	% 消除噪声
    speakerData(i).mfcc=melcepst(y,8000);
end
 set(handles.text8,'string','音频特征提取完成！')
 save speakerData speakerData		% 保存数据
% hObject    handle to pushbutton4 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)


% --- Executes on button press in pushbutton5.
function pushbutton5_Callback(hObject, eventdata, handles)
%GMM训练
fprintf('\n训练每个语者的高斯混合模型...\n\n');
load speakerData.mat
gaussianNum=12;					
speakerNum=length(speakerData);  %获取样本个数


for i=1:speakerNum
	fprintf('\n为第%d个语者%s训练GMM……', i,speakerData(i).name(1:end-4));
	[speakerGmm(i).mu, speakerGmm(i).sigm,speakerGmm(i).c] = gmm_estimate(speakerData(i).mfcc(:,5:12)',gaussianNum,20);  %转置正确
 
end

fprintf('\n');
save speakerGmm speakerGmm;   %保存样本GMM
% hObject    handle to pushbutton5 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)


% --- Executes on button press in pushbutton6.
function pushbutton6_Callback(hObject, eventdata, handles)
clc
close all
% hObject    handle to pushbutton6 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)