使用Python對音頻進行特徵提取

寫在前面

由於喜歡玩兒音樂遊戲，因此打算研究一下如何用深度學習的模型生成音遊的譜面。這篇文章主要目的是介紹或者總結一些音頻的知識和代碼。

恩。若是沒玩兒過的話，音樂遊戲大概是下面這個樣子。

下面進入正題。

我Google了一下，找到了這篇文章：Music Feature Extraction in Python。而後這篇文章介紹完了還有一個歌曲分類的實踐：Classification of Music into different Genres using Keras。

下面的內容會主要參考一下這兩篇文章，並加入一些個人理解。內容以下：

• 聲音信號介紹
• 使用Python對音頻進行特徵提取
• 使用Keras對歌曲的題材進行分類

主要涉及的背景知識有：

• 傅里葉變換
• 採樣定理
• Python
• 機器學習

聲音基礎知識

音頻信號

首先先百度一下音頻信號，

音頻信號是（Audio）帶有語音、音樂和音效的有規律的聲波的頻率、幅度變化信息載體。 根據聲波的特徵，可把音頻信息分類爲規則音頻和不規則聲音。其中規則音頻又能夠分爲語音、音樂和音效。規則音頻是一種連續變化的模擬信號，可用一條連續的曲線來表示，稱爲聲波。聲音的三個要素是音調、音強和音色。聲波或正弦波有三個重要參數：頻率 、幅度和相位 ，這也就決定了音頻信號的特徵。

整體來講：音頻信號就是不一樣頻率和相位的正弦波的一個疊加。

通常的聲音大概就是這個樣子。

橫軸是時間，縱軸是聲音的幅度。由於本質上就是正弦波的一個疊加，因此看到實際上是有正有負的。

人耳聽力頻率範圍

生活中存在各類正弦波，但並非全部的波都能被人耳聽到。好比說咱們手機通訊的信號，wifi信號，以及陽光都是一種波，但並不能被人聽見。

正常人耳聽見聲音的頻率範圍是 20Hz 到 2 萬 Hz 。相同強度的聲音如頻率不一樣的話，聽起來的響度是不同的。至敏感的頻率是 3000 和 4000Hz 。

因此聲波的信號基本上只要關注2wHz之內就行了。

奈奎斯特採樣定理

聲音本質上是一種模擬信號，但在計算機或者在其餘數字設備上傳輸時，咱們要把模擬信號轉換爲數字信號，須要進行採樣。

奈奎斯特採樣定理以下：

在進行模擬/數字信號的轉換過程當中，當採樣頻率fs.max大於信號中最高頻率fmax的2倍時(fs.max>2fmax)，採樣以後的數字信號完整地保留了原始信號中的信息。

這個定理描述的很簡單，證實其實也不難，對於聲音信號，只要採樣的頻率大於2*2wHz=4WHz的話，咱們就能夠聽到無損的音質了。

上面說人耳聽力敏感的範圍主要是在4000Hz，因此咱們通常聽到的音樂實際上是使用8000Hz頻率進行採樣的。這裏能夠看下最近比較火的芒種這首歌。

這首歌的時間是3分36秒也就是216秒，它的標準品質的大小是3.3M。這裏能夠計算下使用8000Hz頻率，16bit進行採樣的話，那麼這個文件的大小是：

大概也就是3.3兆了。

使用Python對音頻進行特徵提取

背景知識大概就上面那些，下面開始實踐環節。

準備工做

這裏使用到了librosa，numpy, sklearn與keras。固然爲了方便開發，我這裏還用了IPython NoteBook。

pip install librosa numpy sklearn tensorflow keras
複製代碼

加載音樂

接下來開始在python中加載音樂，這裏我選的是一首我超級喜歡的一首音遊的歌曲。連接在這裏：visions.mp3

import librosa

x , sr = librosa.load("visions.mp3", sr=8000)
print(x.shape, sr)
複製代碼

這裏x是音頻信號的數字信息，能夠看到是一維的，sr是採樣頻率，用8000就行了。

而後能夠看下這首歌在時域上的波形。

%matplotlib inline
import matplotlib.pyplot as plt
import librosa.display
plt.figure(figsize=(14, 5))
librosa.display.waveplot(x, sr=sr)
複製代碼

運行後，效果以下：

橫軸爲時間，縱軸爲幅度（amplitude）。

Spectogram

上面是音樂時域上的信息。接下來須要對音樂的頻率信息進行分析。

Spectogram這裏沒找到好的翻譯， 大概表示的意思就是：時變的頻譜圖。

短時傅里葉變換-STFT

咱們知道對於週期信號，可使用傅里葉變換能夠將時間域的信息變換到頻率域上。

好比說看上面這張圖的這個波實際上是由三個頻率的正弦波構成的，能夠看到從頻率域的原點往上延伸，頻率愈來愈高，能夠看到對應的波形愈來愈密集。

在文章中提到了一個短時傅里葉變換的概念，做爲一名本科學通訊的專業學生，我也不記得老師有講過這個變換，只記得有過離散傅里葉變換。原做者給了個視頻：The Short Time Fourier Transform | Digital Signal Processing

不過在看視頻以前，憑藉個人想象，感受這個變換既然叫短時傅里葉變換，那應該是把音頻信號拆分開，而後進行傅里葉變換。由於若是對整個音樂進行傅里葉變換的話，咱們只能看到全部的頻率信息，但不能肯定每一個頻率的發生的大概時間。只能看到有某些頻率的聲音而已。

後來在網上找了一個我感受比較好的解釋：

短時傅立葉變換基本思想是將信號加滑動時間窗，並對窗內信號作傅立葉變換，獲得信號的時變頻譜。

另外，還有一個東西能夠感覺STFT，這裏能夠打開網易雲音樂，找到一首歌，打開動效。

能夠看到不斷變化的各波形，拆分開來看，其實表明着就當前這一短期內包含的聲音的頻率信息。上面的這張圖其實也就是對當前時刻的音頻信號，進行DFT獲得的。

畫出時變的頻譜圖

接下來咱們調用librosa的stft方法能夠直接獲得短時傅里葉變換的結果。

X = librosa.stft(x)
Xdb = librosa.amplitude_to_db(abs(X))   # 把幅度轉成分貝格式
plt.figure(figsize=(14, 5))
librosa.display.specshow(Xdb, sr=sr, x_axis='time', y_axis='hz')
plt.colorbar()
複製代碼

效果以下：

這裏橫軸是時間，縱軸是頻率，顏色則表明分貝（聲音的響度），能夠看到越紅的地方信號音量越大。

特徵提取

接下來，介紹了一些音樂中的主要特徵，這裏主要翻譯原做者的內容。

過零率（Zero Crossing Rate）

過零率(Zero Crossing Rate,ZCR)是指在每幀中,語音信號經過零點(從正變爲負或從負變爲正)的次數。這個特徵已在語音識別和音樂信息檢索領域獲得普遍使用，是金屬聲音和搖滾樂的關鍵特徵。

回到上面的圖，咱們把一個時間上放大。

n0 = 9000
n1 = 9100
plt.figure(figsize=(14, 5))
plt.plot(x[n0:n1])
plt.grid()
複製代碼

效果以下：

這裏有7個過零點，能夠經過下面的方法進行計算。

zero_crossings = librosa.zero_crossings(x[n0:n1], pad=False)
print(sum(zero_crossings))
複製代碼

頻譜中心（Spectral Centroid）

頻譜中心表明聲音的「質心」，又稱爲頻譜一階距。頻譜中心的值越小，代表越多的頻譜能量集中在低頻範圍內。

#spectral centroid -- centre of mass -- weighted mean of the frequencies present in the sound
import sklearn
spectral_centroids = librosa.feature.spectral_centroid(x[:80000], sr=sr)[0]
# Computing the time variable for visualization
frames = range(len(spectral_centroids))
t = librosa.frames_to_time(frames, sr=8000)
# Normalising the spectral centroid for visualisation
def normalize(x, axis=0):
    return sklearn.preprocessing.minmax_scale(x, axis=axis)
#Plotting the Spectral Centroid along the waveform
librosa.display.waveplot(x[:80000], sr=sr, alpha=0.4)
plt.plot(t, normalize(spectral_centroids), color='r')
複製代碼

這裏用x[:80000]表示音樂的前10秒。

.spectral_centroid 來計算每一幀的頻譜中心。

.frames_to_time將幀轉換爲時間，time[i] == frame[i]。【由於stft是一個窗口(幀)一個窗口取的，和採樣頻率並不對應，因此有個轉換】

爲了更好的可視化，對頻譜中心進行了歸一化。

頻譜滾降點（Spectral Rolloff）

頻譜滾降點的意思，我翻譯過來大概是：比該頻率低的頻率的全部能量大於必定比例的整個頻譜的能量，一般這個比例爲0.85。

感受仍是有點兒彆扭，用公式比較好理解一些。

下面看下代碼，和上面頻譜中心的代碼其實相似，也是一幀一幀的進行特徵抽取。

spectral_rolloff = librosa.feature.spectral_rolloff(x, sr=sr)[0]
librosa.display.waveplot(x, sr=sr, alpha=0.4)
plt.plot(t, normalize(spectral_rolloff), color='r')
複製代碼

.spectral_rolloff 來計算每一幀的頻譜滾降點。

MFCC (梅爾頻率倒譜系數)

MFCC是音頻信號特徵中最重要的一個，基本上處理音頻信號就會用到。（做爲一名通訊學士，我也是才知道的）。

信號的MFCC參數是一個小集合的特徵（通常10-20個），它可以簡潔的表示頻譜的包絡。

mfccs = librosa.feature.mfcc(x, sr=sr)
print(mfccs.shape)
#Displaying the MFCCs:
librosa.display.specshow(mfccs, sr=sr, x_axis='time')
複製代碼

.mfcc 用來計算信號的MFCC參數。

經過打印mfccs.shape，能夠看看每一幀裏面有多少維的MFCC特徵。第一個參數是mfcc參數的維度，第二個參數是幀數。

這裏一共3177幀，每一幀有20維特徵。

使用Keras對歌曲的題材進行分類

接下來，完成了以上的對音頻信號的特徵提取，就能夠進行一波機器學習了。（這裏不會涉及到機器學習的背景知識，默認讀者不是機器學習的小白。）

數據準備

數據集下載

這裏將要使用GTZAN genre collection。這個數據集包含了10個題材，每一個題材包含了100個30s長的音樂。

這裏將使用Tensorflow後端的Keras進行編碼。

加載數據到內存

這裏首先加載這些數據。

import librosa
import numpy as np
import os

genres = 'blues classical country disco hiphop jazz metal pop reggae rock'.split()

data_set = []
label_set = []

label2id = {genre:i for i,genre in enumerate(genres)}
id2label = {i:genre for i,genre in enumerate(genres)}

print(label2id)

for g in genres:
    print(g)
    for filename in os.listdir(f'./genres/{g}/'):
        songname = f'./genres/{g}/{filename}'
        y, sr = librosa.load(songname, mono=True, duration=30)
        chroma_stft = librosa.feature.chroma_stft(y=y, sr=sr)
        rmse = librosa.feature.rms(y=y)
        spec_cent = librosa.feature.spectral_centroid(y=y, sr=sr)
        spec_bw = librosa.feature.spectral_bandwidth(y=y, sr=sr)
        rolloff = librosa.feature.spectral_rolloff(y=y, sr=sr)
        zcr = librosa.feature.zero_crossing_rate(y)
        mfcc = librosa.feature.mfcc(y=y, sr=sr)

        to_append = f'{np.mean(chroma_stft)} {np.mean(rmse)} {np.mean(spec_cent)} {np.mean(spec_bw)} {np.mean(rolloff)} {np.mean(zcr)}'    

        for e in mfcc:
            to_append += f' {np.mean(e)}'

        data_set.append([float(i) for i in to_append.split(" ")])
        label_set.append(label2id[g])
複製代碼

能夠看到除了上面講到的一些特徵以外，原做者這裏還說了一些其餘的特徵，好比說rmse，chroma_stft這種。這裏就再也不一一介紹了。

整體來講，這裏使用chroma_stft，rmse，spec_cent，spec_bw，rolloff，zcr的六個特徵的均值，再加上mfcc的20個特徵的均值，總共26個特徵，來表示一個音樂的。

建立數據集

接下來咱們對數據進行一下歸一化，對標籤進行one-hot編碼。

from sklearn.preprocessing import StandardScaler
from keras.utils import np_utils

scaler = StandardScaler()
X = scaler.fit_transform(np.array(data_set, dtype = float))
y = np_utils.to_categorical(np.array(label_set))
複製代碼

效果以下：

而後把測試集和訓練集進行分割。

from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2)
複製代碼

建立模型

這裏的咱們的特徵不多，因此用幾層神經網絡就能夠了。

from keras import models
from keras.layers import Dense, Dropout

def create_model():
    model = models.Sequential()
    model.add(Dense(256, activation='relu', input_shape=(X_train.shape[1],)))
    model.add(Dense(128, activation='relu'))
    model.add(Dense(64, activation='relu'))
    model.add(Dropout(0.5))
    model.add(Dense(10, activation='softmax'))

    return model

model = create_model()
複製代碼

這裏建立了一個包含三個隱藏層的神經網絡，最後一層輸出的是分類層，由於是10類，因此最後一層是10個單元。（這裏相比原做者的代碼多了一層Dropout減小數據過擬合）

編譯模型

而後是編譯模型，咱們這裏是一個分類問題，因此使用類別交叉熵函數categorical_crossentropy，而後優化器選擇Adam，評價指標選擇正確率。

model.compile(optimizer='adam',
              loss='categorical_crossentropy',
              metrics=['accuracy'])
複製代碼

訓練與評估

接下來使用fit方法進行訓練，訓練50輪。

model.fit(X_train, y_train, epochs=50, batch_size=128)
複製代碼

訓練很快就完了，而後使用evaluate方法進行評估。

test_loss, test_acc = model.evaluate(X_test,y_test)
print('test_acc: ',test_acc)
複製代碼

這裏的測試準確率也就百分之70不到，原做者後面還有使用模型預測的代碼，感受都沒有必要了。

總結

這篇博客主要是介紹了一些聲音方面的特徵提取的知識和代碼，關於最後的案例看着玩玩兒就行，基本不具備實用性。