python 处理音频并提取特征可视化

深度学习和机器学习在图像和视觉方面有很多的研究，而且取得了很好的效果，但是在音频处理方面还是比较弱，还有很大的提升空间。图像上有openCV这个强大的库，而在音频方面有哪些比较合适的库呢？这里说的是python支持的库，因为本人接触的主要还是python。所以今天就介绍一下python处理音频的库。

进入正题：Python有一些很棒的音频处理库，比如Librosa和PyAudio，今天着重学习一下librosa库，这个库的安装非常简单，和其它库安装一样，直接pip install librosa就行了。以下代码在jupyter notebook中完成。

In[1]:

import librosa as lr
 audio_path = './data/js.wav'
 x, sr = lr.load(audio_path) #sr-sample rate采样率
 print(type(x),type(sr))
 print(x.shape,sr)

Out[2]:

<class 'numpy.ndarray'> <class 'int'>
 (127200,) 22050

#这会将音频时间序列作为numpy数组返回，默认采样率sr为22KHZ mono我们可以通过以下方式更改：
 lr.load(audio_path, sr = 44100)
 #以44.1KHZ重新采样，或者禁重新采样。采样率是每秒传输的音频样本数。
 lr.load(audio_path, sr=None)

Out[3]:

(array([-3.0517578e-04, -2.7465820e-04, -3.6621094e-04, ...,
 -1.2207031e-04, -9.1552734e-05, -9.1552734e-05], dtype=float32), 44100)

In [4]:

#使用Ipython.display.Audio播放音频
 import IPython.display as ipd
 ipd.Audio(audio_path)

可视化音频：

可视化音频的目的主要让你更直观的理解音频特征

In [5]:

#波形：我们可以绘制音频数组librosa.display.waveplot:
 %matplotlib inline 
 import matplotlib.pyplot as plt
 import librosa.display
plt.figure(figsize = (14,5))
 librosa.display.waveplot(x, sr = sr)

In[8]

#谱图：谱图是通过视觉表示频谱的频率，声音或其它信号，因为它随时间变化。频谱图有时被称为声纹或语音图。当数据在3D中表示时，它们可以称为waterfalls。
 #在二维阵列中，第一轴是频率，第二轴是时间。我们可以使用显示频谱图：librosa.display.specshow:
 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()

In[9]

#从轴表示频率（从0到10kHZ），横轴表示剪辑的时间，由于我们看到的所有动作都发生在频谱的底部，我们可以将频谱轴转换为对数轴。
 plt.figure(figsize=(14,5))
 librosa.display.specshow(Xdb, sr=sr, x_axis='time', y_axis='log')
 plt.colorbar()

In[10]:

#写音频：librosa.output.write_wav 将Numpy数组保存到wav文件, 下面举例说明它的用法。
 librosa.output.write_wav('example.wav', x, sr)

#创建音频信号：现在让我们创建一个220Hz的音频信号，音频信号是一个numpy数组，所以我们将创建一个并将其传递给音频函数

 import numpy as np
 #sample rate
 sr = 22050
 #seconds
 T = 5.0
 #time variable,linspace()产生一个固定间隔的一个序列的一维数组
 t = np.linspace(0, T, int(T*sr), endpoint=False)
 #pure sine wave at 220Hz
 x = 0.5*np.sin(2*np.pi*220*t)
#播放这个音频
 ipd.Audio(x, rate=sr)

In [11 ]:

#保存这个音频信号
 librosa.output.write_wav('tone_220.wav', x, sr)

音频特征提取：

每个音频信号都包含许多特征。但是，我们必须提取与我们试图解决的问题相关的特征，提取要使用它们进行分析的特征的过程称为特征提取 接下来我们详细了解一下我们感兴趣的特征：

过零率：

该特征在语音识别和音乐信息检索中被经常使用。

过零率（Zero Crossing Rate，ZCR）是指在每帧中，语音信号通过零点（从正变为负或从负变为正）的次数。 这个特征已在语音识别和音乐信息检索领域得到广泛使用，比如敲击声音的分类，过零率是分类的关键特征。

特性：

1. 一般而言，清音（unvoiced sound）和环境噪音的ZCR都大于浊音（voiced sound）；
2. 由于清音和环境噪音的ZCR大小相近，因而不能够通过ZCR来区分它们；
3. 在实际当中，过零率经常与短时能量特性相结合来进行端点检测，尤其是ZCR用来检测清音的起止点；
4. 有时也可以用ZCR来进行粗略的基频估算，但这是非常不可靠的，除非有后续的修正（refine）处理过程。

那我们看看如何计算音频片段的过零率：

In [12]:

#过零率的计算
#导入音频文件
 x, sr = librosa.load('./data/js.wav')
 plt.figure(figsize=(14,5))
 librosa.display.waveplot(x, sr=sr)
#放大局部查看
 n0 = 9000
 n1 = 9100
 plt.figure(figsize=(14,5))
 plt.plot(x[n0:n1])
 plt.grid()

#这里举一个不需要librosa库来实现过零率的计算方法
import math
import numpy as np
def ZeroCR(waveData,frameSize,overLap):
 wlen = len(waveData)
 step = frameSize - overLap
 frameNum = math.ceil(wlen/step)
 zcr = np.zeros((frameNum,1))
 for i in range(frameNum):
 curFrame = waveData[np.arange(i*step,min(i*step+frameSize,wlen))]
 curFrame = curFrame - np.mean(curFrame) # 过零点判断
 zcr[i] = sum(curFrame[0:-1]*curFrame[1::]<=0)
 return zcr
#通过上面定义的函数来计算过零率
import math
import wave
import numpy as np
import pylab as pl
# ============ 测试一下算法 =============
# 读取音频
fw = wave.open('./data/js.wav','rb')
params = fw.getparams()
print(params)
nchannels, sampwidth, framerate, nframes = params[:4]
str_data = fw.readframes(nframes)
wave_data = np.frombuffer(str_data, dtype=np.short)
wave_data.shape = -1, 1
#wave_data = wave_data.T
fw.close()
# 计算过零率
frameSize = 256
overLap = 0
zcr = ZeroCR(wave_data,frameSize,overLap)
# 可视化结果
time = np.arange(0, len(wave_data)) * (1.0 / framerate)
time2 = np.arange(0, len(zcr)) * (len(wave_data)/len(zcr) / framerate)
pl.figure(figsize=(14,10))
pl.subplot(211)
pl.plot(time, wave_data)
pl.ylabel("Amplitude")
pl.subplot(212)
pl.plot(time2, zcr)
pl.ylabel("ZCR")
pl.xlabel("time (seconds)")
pl.show()

In[13]:

#验证一下有几个过零点, 本例中的片段含有3个
 zero_crossings = librosa.zero_crossings(x[n0:n1], pad=False)
 print(sum(zero_crossings))
3

光谱质心：

谱质心(Spectral Centroid)是描述音色属性的重要物理参数之一，是频率成分的重心，是在一定频率范围内通过能量加权平均的频率，其单位是Hz。它是声音信号的频率分布和能量分布的重要信息。在主观感知领域，谱质心描述了声音的明亮度，具有阴暗、低沉品质的声音倾向有较多低频内容，谱质心相对较低，具有明亮、欢快品质的多数集中在高频，谱质心相对较高。该参数常用于对乐器声色的分析研究。

如果有两首曲子，一首来自布鲁斯类型，另一首属于金属。与长度相同的布鲁斯流派歌曲相比，金属歌曲在最后有更多的频率。因此，布鲁斯歌曲的光谱质心将位于其光谱中间附近，而金属歌曲的光谱质心将朝向它的末端。 librosa.feature.spectral_centroid计算信号中每帧的光谱质心：

In [14]:

#计算光谱质心
 spectral_centroids = librosa.feature.spectral_centroid(x, sr=sr)[0]
 spectral_centroids.shape

(249,)

In [15]:

#可视化每一帧的质心

 import sklearn
frames = range(len(spectral_centroids))
 t = librosa.frames_to_time(frames)
#归一化光谱质心
 def normalize(x, axis=0):
 return sklearn.preprocessing.minmax_scale(x, axis=axis)
#沿着音频波形绘制光谱质心
 plt.figure(figsize=(14,5))
 librosa.display.waveplot(x, sr=sr, alpha=0.4)
 plt.plot(t, normalize(spectral_centroids), color='r')

光谱衰减：

它是信号形状的度量。librosa.feature.spectral_rolloff计算信号中每帧的滚降系数：

In [16]:

#计算光谱衰减：
 spectral_rolloff = librosa.feature.spectral_rolloff(x+0.01, sr=sr)[0]
 plt.figure(figsize=(14,5))
 librosa.display.waveplot(x, sr=sr, alpha=0.4)
 plt.plot(t, normalize(spectral_rolloff), color='r')

梅尔频率倒谱系数：

声音信号是连续变化的，为了将连续变化信号简化，我们假设在一个短时间尺度内，音频信号不发生改变。因此将信号以多个采样点集合成一个单位，称为'''讯框'''。一个讯框多为20-40毫秒，如果讯框长度更短，那每个讯框内的采样点将不足以做出可靠的频谱计算，但若长度太长，则每个讯框信号会变化太大。

预强化的目的就是为了消除发声过程中，声带和嘴唇造成的效应，来补偿语音信号受到发音系统所压抑的高频部分。并且能突显高频的共振峰。

由于信号在时域上的变化通常很难看出信号的特性，所以通常透过傅里叶变换将它变换成频域上的能量分布来观察，不同的能量分布，就能代表不同语音的特性。

由于能量频谱中还存在大量的无用讯息，尤其人耳无法分辨高频的频率变化，因此让频谱通过梅尔滤波器。梅尔滤波器，也就是一组20个非线性分布的三角带通滤波器（Triangular Bandpass Filters），能求得每一个滤波器输出的对数能量。必须注意的是：这 20 个三角带通滤波器在'''梅尔刻度'''的频率上是平均分布的。 梅尔频率代表一般人耳对于频率的感受度，由此也可以看出人耳对于频率 f 的感受是呈对数变化的。

MFCC特征在加性噪声的情况下并不稳定，因此在语音识别系统中通常要对其进行归一化处理(normalise)以降低噪声的影响。一些研究人员对MFCC算法进行修改以提升其鲁棒性，如在进行DCT之前将log-mel-amplitudes提升到一个合适的能量(2到3之间)，以此来降低低能量成分的影响。

信号的Mel频率倒谱系数（MFCC)是一小组特征（通常约10-20），其简明地描述了频谱包络的整体形状，它模拟了人声的特征。

In [17]:

#MFCC计算
 x, sr = librosa.load('./data/js.wav')
#通过音频信号计算MFCC
 mfccs = librosa.feature.mfcc(x, sr=sr)
 print(mfccs.shape)
(20, 249)

In [18]:

#显示MFCC图谱
 plt.figure(figsize=(14,5))
 librosa.display.specshow(mfccs, sr=sr, x_axis='time')

In[19]:

#以上计算了249帧的20个MFCC。我们还可以执行特征缩放，使得每个系数维度具有零均值和单位方差：
 #import sklearn
mfccs = sklearn.preprocessing.scale(mfccs, axis=1)
 print(mfccs.mean(axis=1))
 print(mfccs.var(axis=1))
 plt.figure(figsize=(14,5))
 librosa.display.specshow(mfccs, sr=sr, x_axis='time')

[-4.3087693e-09 1.9150086e-09 -4.7875215e-09 -1.1968804e-10

-4.2489257e-09 5.7929011e-08 -9.5750430e-09 4.5481454e-09

2.3937607e-09 1.6905936e-09 4.7875215e-09 -6.8222183e-09

4.2190034e-09 7.1812827e-09 4.0693933e-09 -8.1387865e-09

-9.9640296e-09 5.9844023e-09 1.7833518e-08 5.0268976e-09]

[1.0000001 1.0000001 1.0000001 1.0000001 0.9999998 1.0000001

1.0000001 1.0000001 1.0000001 1.0000001 1.0000001 1.0000001

1.0000001 1.0000001 1. 0.99999994 1.0000002 1.0000001

0.99999976 0.9999999 ]

色度频率

色度频率是音乐音频有趣且强大的表示，其中整个频谱被投影到12个区间，代表音乐八度音的12个不同的半音（或色度）， librosa.feature.chroma_stft用于计算色度频率：

In [20]:

#色度频率
 x, sr = librosa.load('./data/cut.wav')
 hop_length = 512
 chromagram = librosa.feature.chroma_stft(x, sr=sr, hop_length=hop_length)
 plt.figure(figsize=(15,5))
 librosa.display.specshow(chromagram, x_axis='time', y_axis='chroma', hop_length=hop_length, cmap='coolwarm')

In[21]:

#导入pydub库用于音频文件的截取，上面用到的片段就是从piano中截取的5秒的片段。
 from pydub import AudioSegment
 main_wav_path = './data/piano.wav' # 取出.wav后缀的文件名
 sound = AudioSegment.from_wav(main_wav_path)
 start_time = 0
 end_time = 5000
 word = sound[start_time:end_time]
 word.export("./data/cut1.wav", format="wav")

<_io.BufferedRandom name='./data/cut1.wav'>

关于音频特征处理的一般性方法介绍完了，基于这些知识，接下来可以考虑进一步研究声纹匹配或者音乐分类项目。