1. 概述
語音是人類之間溝通交流的最直接也是最快捷方便的一種手段，而實現人類與計算機之間暢通無阻的語音交流，一直是人類追求的一個夢想。
伴隨着移動智能設備的普及，各家移動設備的廠家也開始在自家的設備上集成了語音識別系統，像Apple Siri、Microsoft Cortana、Google Now等語音助手的出現，使得人們在使用移動設備的同時，也能夠進行語音交流，極大的方便了人們的生活。但是此類助手也存在一些尷尬的瞬間，例如在一些工作場合或者聚會的場合，某人的一句“Hey Siri”就可能喚醒多台蘋果設備，使用者難免尷尬困惑。
而此類予語音助手背后，均是一種被稱作“聞聲識人”的計算機技術，稱為語音識別。語音識別技術屬於生物認證技術，而其中的說話人識別（speaker recognize，SR）是其中的一種，該技術通常也被稱為聲紋識別技術，該技術是一項通過語音波形中反映說話人生理特征和行為特征的一組語音參數，自動識別說話人身份的技術。其核心是通過預先錄入說話人的聲音樣本，提取出說話人獨一無二的語音特征並存入數據庫，應用的時候將待驗證的語音進行特征提取並與數據庫中的特征進行匹配，以確定說話人的身份。
1.1 什么是聲紋？
聲紋（voiceprint）是用電聲學儀器顯示的攜帶者言語信息的聲波頻譜，是由波長、頻率以及強度等百余種特征維度組成的生物特征，具有穩定性、可測量性以及唯一性等特點。
人類語言的產生是由人體語言中樞與發生器官之間進行的一個復雜的生物物理反應過程。發聲器官如舌頭、牙齒、喉嚨、肺、鼻子在尺寸和形態上因人而異，所有任何兩個人的聲波圖譜都有一定的差異性。
每個人的語音聲學特征既有相對穩定性，又有個體差異性。這種差異可能來自生理、病理、心理、模擬、偽裝等，也可能會周圍環境的干擾相關。
由於每個人的發生器官都有其獨特性，因此在一般情況下，人們仍然能區別不同的人的聲音或者判斷是否是同一個人的聲音。
聲紋不像圖像那樣的直觀，在實際的分析中，可以通過波形圖和語譜圖進行繪制展現，例如下圖是一段從1到10的讀數語音文件對應的波形圖和語譜圖（上部分為聲音波形圖，下部分為聲音語譜圖）：
import wave import numpy as np import matplotlib.pyplot as plt fw = wave.open('test.wav','r') soundInfo = fw.readframes(-1) soundInfo = np.fromstring(soundInfo,np.int16) f = fw.getframerate() fw.close() plt.subplot(211) plt.plot(soundInfo) plt.ylabel('Amplitude') plt.title('Wave from and spectrogram of test.wav') plt.subplot(212) plt.specgram(soundInfo,Fs = f, scale_by_freq = True, sides = 'default') plt.ylabel('Frequency') plt.xlabel('time(seconds)') plt.show()
- 語譜圖更簡單的繪制方法，可參考 scipy.signal.spectrogram。
- 語譜圖繪制的原理，可參考 Create audio spectrograms with Python。
與其他的生物認證技術如指紋識別、人臉識別、虹膜識別等相同，聲紋識別具有不會遺忘、無需記憶和使用方便等優點。在生物認證技術領域，說話人識別技術以其獨特的方便性、經濟性和准確性收到人們的廣泛關注，並日益成為人們日常生活和工作中重要且普及的安全認證方式。
但是，說話人識別有着其他生物認證技術所不具有的優勢：
用戶接受度高：以聲音作為識別特征，因其非接觸性和自然醒，用戶易接受。用戶不用刻意的用手指觸摸相應的傳感器上，也不用將眼睛湊向攝像頭，只需要簡單的說一兩句話即可完成識別認證。
設備成本低：對輸入設備如麥克風，攝像頭等沒有特別的要求，特征提取，模型訓練和匹配只需要普通的計算機即可完成。
其他生物認證特征技術各有其劣勢：指紋識別需要特殊的傳感器芯片，虹膜識別精確度較高，但是設備較為昂貴。
在遠程應用和移動互聯網環境下優勢明顯：通過電話、移動設備進行身份認證，聲音是最具優勢的生物特征，語音控制也逐漸成為流行的交互形式，以聲音為特征的身份鑒別技術也越發重要。
1.2 聲紋識別技術的歷史
聲紋識別技術的研究始於20世紀30年代，早期的工作主要集中於人耳聽辨實驗和探討聽音識別的可能性方面。隨着研究手段和計算機技術的發展，研究工作逐漸脫離了單純的人耳聽辨，使得通過機器自動識別人的聲音稱為可能。在這個過程中也出現了很多不同的計算機技術，從早期的模板匹配到最新的深度學習技術，均在不斷的刷新着語音識別技術手段。整體來看，聲紋識別技術的發展經歷了七個技術演進之路，詳見下圖（下圖來自speakin）：
1.3 聲紋識別的種類
聲紋識別根據實際應用的范疇可以分為 1:1識別 和 1:N識別兩種：
- 1:1識別：指確定待識別的一段語音是否來自其所聲明的目標說話人，即確認目標說話人是目標說話人的過程。通常應用於電子支付、智能硬件、銀行證券交易等。1:1識別有兩個系統的性能評價參量，分別為
- 錯誤接受率(False Acceptation Rate, FAR)：將非目標說話人判別為目標說話人造成的錯誤率
- 錯誤拒絕率(False Rejection Rate, FRR)：將目標說話人誤識成非目標說話人造成的錯誤率
- 對安全性要求越高，則設定閾值越高，此時接受目標說話人的條件越嚴格，即FRR越高，FAR越低；對用戶體驗要求越高，則設定閾值越低，此時接受目標說話人的條件越寬松，即FAR越高，FRR越低。在聲紋系統中，可以通過設定不同的閾值來平衡FAR和FRR。
- 1:N識別：指判定待識別語音屬於目標說話人模型集合中的哪一個人，即在N個人中找到目標說話人的過程。通常應用於公安司法、軍隊國防等。
2. 語音的特征提取方法概述
語音是一種數字信號，其數字⾳頻的采樣率為44100Hz（根據乃奎斯特取樣定理得出的結果，在模擬訊號數字化的過程中，如果保證取樣頻率大於模擬訊號最高頻率的2倍，就能100%精確地再還原出原始的模擬訊息。音頻的最高頻率為20kHz，所以取樣率至少應該大於40kHz，為了留一點安全系數，再考慮到工程上的習慣，最終選擇了44.1kHz這個數值）。通常情況下使用傅里葉變換將信號在時域與頻域之間進行轉換，而頻譜圖可以顯示傅里葉變換后的振幅與時間和頻率的對應關系。
2.1 特征提取方法
對於語音識別系統而言，所提取的特征參數需要能夠反映特定發信的信息，在說話人無關的系統中，更要求參數能夠反映不同說話人相同發音的信息，要求說話人的特征參數要能夠代表特定的說話人，能夠區分不同說話人相同語音之間的差異，最好能夠做到與具體的發音內容無關，也稱為文本無關。
在語音特征參數提取技術的發展歷程中，線性預測編碼（Linear Predictive Coding, LPC）被廣泛應用於語音特征參數的提取，其中包括LPC系數、反射LPC系數、面積函數和LPC倒譜系數，能夠很好的反映語音的聲道特征，但是卻對語音的其他特征無能為力。
不同於LPC等通過對人的發聲機理進行研究而得到的聲學特征，Mel倒譜系數MFCC是受人的聽覺系統研究成果推出而導出的聲學特征。根據人耳聽覺機理的研究發現，人耳對不同頻率的聲波有不同的聽覺靈敏度。從200Hz到5000Hz的語音信號對語音的清晰度影響最大。人們從低頻到高頻這一段頻帶內按臨界帶寬的大小由密到疏安排一組帶通濾波器，對輸入信號進行濾波。將每個帶通濾波器輸出的信號能量作為信號的基本特征，對此特征經過進一步處理后就可以作為語音的輸入特征。由於這種特征不依賴於信號的性質，對輸入信號不做任何的假設和限制，又利用了聽覺模型的研究成果。因此，這種參數比基於聲道模型的LPC相比具有更好的魯棒性，更符合人耳的聽覺特性，而且當信噪比降低時仍然具有較好的識別性能。
MFCC（MeI-Freguency CeptraI Coefficients）是需要語音特征參數提取方法之一，因其獨特的基於倒譜的提取方式，更加的符合人類的聽覺原理，因而也是最為普遍、最有效的語音特征提取算法。MFCC是在Mel標度頻率域提取出來的倒譜系數，Mel標度描述了人耳對頻率感知的非線性特性。
2.2 MFCC語音特征提取
MFCC 語音特征的提取過程，如下圖：
需要對語音信號進行預加重、分幀、加窗等等處理，而這些處理的方式均是為了能夠最大化語音信號的某些信息，以達到最好特征參數的提取。
2.2.1 預加重
預加重其實就是將語音信號通過一個高通濾波器，來增強語音信號中的高頻部分，並保持在低頻到高頻的整個頻段中，能夠使用同樣的信噪比求頻譜。在本實驗中，選取的高通濾波器傳遞函數為： 
式中a的值介於0.9-1.0之間，我們通常取0.97。同時，預加重也是為了消除發生過程中聲帶和嘴唇的效應，來補償語音信號受到發音系統所抑制的高頻部分，也為了突出高頻的共振峰。
def pre_emphasis(signal, coefficient=0.97): '''對信號進行預加重''' return numpy.append(signal[0], signal[1:] - coefficient * signal[:-1])
2.2.2 分幀
分幀是指在跟定的音頻樣本文件中，按照某一個固定的時間長度分割，分割后的每一片樣本，稱之為一幀，這里需要區分時域波形中的幀，分割后的一幀是分析提取MFCC的樣本，而時域波形中的幀是時域尺度上對音頻的采樣而取到的樣本。
分幀是先將N個采樣點集合成一個觀測單位，也就是分割后的幀。通常情況下N的取值為512或256，涵蓋的時間約為20-30ms。也可以根據特定的需要進行N值和窗口間隔的調整。為了避免相鄰兩幀的變化過大，會讓兩相鄰幀之間有一段重疊區域，此重疊區域包含了M個取樣點，一般M的值約為N的1/2或1/3。
語音識別中所采用的信號采樣頻率一般為8kHz或16kHz。以8kHz來說，若幀長度為256個采樣點，則對應的時間長度是256/8000×1000=32ms。本次實驗中所使用的采樣率(Frames Per Second)16kHz，窗長25ms（400個采樣點），窗間隔為10ms（160個采樣點）。
def audio2frame(signal, frame_length, frame_step, winfunc=lambda x: numpy.ones((x,))): '''分幀''' signal_length = len(signal) frame_length = int(round(frame_length)) frame_step = int(round(frame_step)) if signal_length <= frame_length: frames_num = 1 else: frames_num = 1 + int(math.ceil((1.0 * signal_length - frame_length) / frame_step)) pad_length = int((frames_num - 1) * frame_step + frame_length) zeros = numpy.zeros((pad_length - signal_length,)) pad_signal = numpy.concatenate((signal, zeros)) indices = numpy.tile(numpy.arange(0, frame_length), (frames_num, 1)) + numpy.tile(numpy.arange(0, frames_num * frame_step, frame_step),(frame_length, 1)).T indices = numpy.array(indices, dtype=numpy.int32) frames = pad_signal[indices] win = numpy.tile(winfunc(frame_length), (frames_num, 1)) return frames * win
2.2.3 加窗
在對音頻進行分幀之后，需要對每一幀進行加窗，以增加幀左端和右端的連續性，減少頻譜泄漏。在提取MFCC的時候，比較常用的窗口函數為Hamming窗。
假設分幀后的信號為 S(n),n=0,1,2…,N-1，其中N為幀的大小，那么進行加窗的處理則為：
W(n)的形式如下：
不同的a值會產生不同的漢明窗，一般情況下a取值0.46。進行值替換后，W(n)則為：
對應的漢明窗時域波形類似下圖：
def deframesignal(frames, signal_length, frame_length, frame_step, winfunc=lambda x: numpy.ones((x,))): '''加窗''' signal_length = round(signal_length) frame_length = round(frame_length) frames_num = numpy.shape(frames)[0] assert numpy.shape(frames)[1] == frame_length, '"frames"矩陣大小不正確，它的列數應該等於一幀長度' indices = numpy.tile(numpy.arange(0, frame_length), (frames_num, 1)) + numpy.tile(numpy.arange(0, frames_num * frame_step, frame_step),(frame_length, 1)).T indices = numpy.array(indices, dtype=numpy.int32) pad_length = (frames_num - 1) * frame_step + frame_length if signal_length <= 0: signal_length = pad_length recalc_signal = numpy.zeros((pad_length,)) window_correction = numpy.zeros((pad_length, 1)) win = winfunc(frame_length) for i in range(0, frames_num): window_correction[indices[i, :]] = window_correction[indices[i, :]] + win + 1e-15 recalc_signal[indices[i, :]] = recalc_signal[indices[i, :]] + frames[i, :] recalc_signal = recalc_signal / window_correction return recalc_signal[0:signal_length]
2.2.4 對信號進行離散傅立葉變換 (DFT)
由於信號在時域上的變換通常很難看出信號的特性，所有通常將它轉換為頻域上的能量分布來觀察，不同的能量分布，代表不同語音的特性。所以在進行了加窗處理后，還需要再經過離散傅里葉變換以得到頻譜上的能量分布。對分幀加窗后的各幀信號進行快速傅里葉變換得到各幀的頻譜。並對語音信號的頻譜取模平方得到語音信號的功率譜。設語音信號的DFT為：
能量的分布為：
在本次實驗中，采用DFT長度 N=512，結果值保留前257個系數。
def deframesignal(frames, signal_length, frame_length, frame_step, winfunc=lambda x: numpy.ones((x,))): '''加窗''' signal_length = round(signal_length) frame_length = round(frame_length) frames_num = numpy.shape(frames)[0] assert numpy.shape(frames)[1] == frame_length, '"frames"矩陣大小不正確，它的列數應該等於一幀長度' indices = numpy.tile(numpy.arange(0, frame_length), (frames_num, 1)) + numpy.tile(numpy.arange(0, frames_num * frame_step, frame_step),(frame_length, 1)).T indices = numpy.array(indices, dtype=numpy.int32) pad_length = (frames_num - 1) * frame_step + frame_length if signal_length <= 0: signal_length = pad_length recalc_signal = numpy.zeros((pad_length,)) window_correction = numpy.zeros((pad_length, 1)) win = winfunc(frame_length) for i in range(0, frames_num): window_correction[indices[i, :]] = window_correction[indices[i, :]] + win + 1e-15 recalc_signal[indices[i, :]] = recalc_signal[indices[i, :]] + frames[i, :] recalc_signal = recalc_signal / window_correction return recalc_signal[0:signal_length]
下圖是有頻譜到功率譜的轉換結果示意圖：
2.2.5 應用梅爾濾波器 (Mel Filterbank)
MFCC考慮到了人類的聽覺特征，先將線性頻譜映射到基於聽覺感知的Mel非線性頻譜中，然后轉換到倒譜上。 在Mel頻域內，人對音調的感知度為線性關系。舉例來說，如果兩段語音的Mel頻率相差兩倍，則人耳聽起來兩者的音調也相差兩倍。Mel濾波器的本質其實是一個尺度規則，通常是將能量通過一組Mel尺度的三角形濾波器組，如定義有M個濾波器的濾波器組，采用的濾波器為三角濾波器，中心頻率為 f(m),m=1,2…M，M通常取22-26。f(m)之間的間隔隨着m值的減小而縮小，隨着m值的增大而增寬，如圖所示：
從頻率到Mel頻率的轉換公式為：
其中 f 為語音信號的頻率，單位赫茲（Hz）。
def hz2mel(hz): '''把頻率hz轉化為梅爾頻率''' return 2595 * numpy.log10(1 + hz / 700.0) def mel2hz(mel): '''把梅爾頻率轉化為hz''' return 700 * (10 ** (mel / 2595.0) - 1)
假如有10個Mel濾波器（在實際應用中通常一組Mel濾波器組有26個濾波器。），首先要選擇一個最高頻率和最低頻率，通常最高頻率為8000Hz，最低頻率為300Hz。使用從頻率轉換為Mel頻率的公式將300Hz轉換為401.25Mels，8000Hz轉換為2834.99Mels，由於有10個濾波器，每個濾波器針對兩個頻率的樣點，樣點之間會進行重疊處理，因此需要12個點，意味着需要在401.25和2834.99之間再線性間隔出10個附加點，如：
$$m(i) = 401.25,622.50,843.75,1065.00,1286.25,1507.50, 1728.74,1949.99,2171.24,2392.49,2613.74,2834.99$$
現在使用從Mel頻率轉換為頻率的公式將它們轉換回赫茲：
$$h(i) = 300,517.33,781.90,1103.97,1496.04,1973.32,2554.33, 3261.62,4122.63,5170.76,6446.70,8000$$
將頻率映射到最接近的DFT頻率：
$$f(i) = 9,16,25,35,47,63,81,104,132,165,206,256$$
於是，我們得到了一個由10個Mel濾波器構成的Mel濾波器組。
def get_filter_banks(filters_num=20, NFFT=512, samplerate=16000, low_freq=0, high_freq=None): '''計算梅爾三角間距濾波器，該濾波器在第一個頻率和第三個頻率處為0，在第二個頻率處為1''' low_mel = hz2mel(low_freq) high_mel = hz2mel(high_freq) mel_points = numpy.linspace(low_mel, high_mel, filters_num + 2) hz_points = mel2hz(mel_points) bin = numpy.floor((NFFT + 1) * hz_points / samplerate) fbank = numpy.zeros([filters_num, NFFT / 2 + 1]) for j in xrange(0, filters_num): for i in xrange(int(bin[j]), int(bin[j + 1])): fbank[j, i] = (i - bin[j]) / (bin[j + 1] - bin[j]) for i in xrange(int(bin[j + 1]), int(bin[j + 2])): fbank[j, i] = (bin[j + 2] - i) / (bin[j + 2] - bin[j + 1]) return fbank
2.2.6 對頻譜進行離散余弦變換 (DCT)
在上一步的基礎上使⽤離散余弦變換，即進⾏了⼀個傅⽴葉變換的逆變換，得到倒譜系數。
由此可以得到26個倒譜系數。只取其[2:13]個系數，第1個用能量的對數替代，這13個值即為所需的13個MFCC倒譜系數。
def lifter(cepstra, L=22): '''升倒譜函數''' if L > 0: nframes, ncoeff = numpy.shape(cepstra) n = numpy.arange(ncoeff) lift = 1 + (L / 2) * numpy.sin(numpy.pi * n / L) return lift * cepstra else: return cepstra
2.2.7 動態差分參數的提取（包括一階微分系數和加速系數）
標准的倒譜參數MFCC只反映了語音參數的靜態特性，語音的動態特性可以用這些靜態特征的差分譜來描述。通常會把動、靜態特征結合起來以有效提高系統的識別性能。差分參數的計算可以采用下面的公式：
上式中，d(t)表示第t個一階微分，c(t)表示第t個倒譜系數，Q表示倒譜系數的階數，K表示一階導數的時間差，可取1或2。將上式的結果再代入就可以得到加速系數。
⾄此，我們計算到了了⾳頻⽂件每⼀幀的39個Mel頻率倒譜系數（13個MFCC+13個一階微分系數+13個加速系數），這些即為一個語音文件的特征數據，這些特征數據可以運用在之后的分類中。
def derivate(feat, big_theta=2, cep_num=13): '''計算一階系數或者加速系數的一般變換公式''' result = numpy.zeros(feat.shape) denominator = 0 for theta in numpy.linspace(1, big_theta, big_theta): denominator = denominator + theta ** 2 denominator = denominator * 2 for row in numpy.linspace(0, feat.shape[0] - 1, feat.shape[0]): tmp = numpy.zeros((cep_num,)) numerator = numpy.zeros((cep_num,)) for t in numpy.linspace(1, cep_num, cep_num): a = 0 b = 0 s = 0 for theta in numpy.linspace(1, big_theta, big_theta): if (t + theta) > cep_num: a = 0 else: a = feat[row][t + theta - 1] if (t - theta) < 1: b = 0 else: b = feat[row][t - theta - 1] s += theta * (a - b) numerator[t - 1] = s tmp = numerator * 1.0 / denominator result[row] = tmp return result
3. 總結
本文針對語音數據的特征提取方法—MFCC進行了簡單的概述和實踐，MFCC是音頻特征處理中比較常用而且很有效的方法。當特征數據提取出來之后，就可以進一步的進行數據的歸一化、標准化，然后應用於機器學習、神經網絡等等模型訓練算法中，以得到能夠識別語音類別的模型。在實際的應用中，可能還需要考慮很多的其他因素，例如源語音數據的采集方法、采集時長、模型的構建方式、模型的部署方式等等因素，因此需要根據業務的具體場景，來進行平衡取舍，以達到識別的時效性、准確性等。
目前關於語音識別相關的研究還在持續中，目標是能夠最小化成本的在移動端部署語音識別相關的功能，提高SDK在人工智能方便的能力等。