FFTNetを読んで，メルケプストラム係数（Mel Cepstral Coefficients）が何かわからなかったので調べてみると，メル周波数ケプストラム係数（MFCC）ばかりが検索上位に来てしまい困った．なので，今後のためにもメルケプストラム係数のメモと，Pythonで出力する方法をまとめる．
また，メルケプストラムと表現している文献が多いので，以下ではメルケプストラムと省略する．
参考
#メルケプストラムとは
統計的声質変換を行うための知識と手法では以下のように説明されている．
スペクトル包絡は音色をよく表現できる音響特徴ですが、次元が大きく扱いづらいという問題があります。 スペクトル包絡は周波数方向に周期性があるので、周波数方向に逆フーリエ変換することで、表現能力をあまり損なわずに低次元な特徴量を得ることができます。 これがケプストラムです。 （ケプストラム【cepstrum】はスペクトル【spectrum】のアナグラムです。）
人間の聴覚を考慮することで、ケプストラムよりも表現能力を保持させた、メルケプストラムという音響特徴量があります。 統計的声質変換にはこの特徴量がよく用いられます。 サンプリング周波数によりますが、スペクトル包絡の次元数は512や1024と大きい一方、メルケプストラムの次元数は60ほどに小さくしても元と同程度の品質で復元できます。 （別の概念でメル周波数ケプストラム係数MFCCというものがありますが、これはメルケプストラムとは異なるので注意。）
##つまり
- スペクトル包絡だとデータが大きいので，表現力をある程度保ちつつ情報量を圧縮したのがケプストラム
- 人間の聴覚特性を考慮したケプストラムがメルケプストラム
#Pythonでの実装
メルケプストラムはSPTK（Signal Processing Toolkit）を使うと簡単に求められる．しかし，SPTKはLinux環境でないと動作しないため，pysptkを利用する．
##pysptkでメルケプストラムを出力
pysptk.mcepで計算できるようだが，イマイチ使い方がわからなかったので，今後調べる．他に，pysptk.sp2mc(powerspec, order, alpha)という関数で出力が可能．pysptkのリポジトリを参考にそれぞれの入力を示す．
###powerspec
スペクトル包絡（spectral envelope）を入力すれば良い．スペクトル包絡はWORLDを使えば算出可能．
こちらもPythonラッパーがあるのでこれを利用する．
###order
メルケプストラムの次数．機械学習等のモデルパラメータに合わせて設定を行う
###alpha
all-pass constantと呼ばれるもの．SPTK Reference manualを参照すると，サンプリング周波数に対して
- 16kHz：0.42
- 10kHz：0.35
- 8kHz：0.31
に設定すると良いとのこと
#出力の確認
試しにLJ Speechを使って，メルケプストラムを出力してみる．
##alphaの設定
SPTKのマニュアルには16kHzまでしか書いていなかったが，LJ Speechのサンプリング周波数は22.5kHzである．そこで，既知の数値から適切な値を予測した．サンプリング周波数に対して，alphaの値が対数っぽく変化してるなと考え，少し乱暴に計算してみた所$$alpha = \frac{1}{2.9}\log fs$$でそれっぽい値が算出できた．ただし，$fs[kHz]$とする．この式でalphaを求めると，0.46となったので，とりあえずこれを使ってメルケプストラムを求める
##サンプルコード
サンプルコードを以下に示す．メルケプストラムの次数は25とした．
import numpy as np
import librosa
import pyworld as pw
import pysptk as sptk
import matplotlib.pyplot as plt
wav, fs = librosa.load('LJ001-0001.wav')
wav = wav.astype(np.float)
#スペクトル包絡算出にWORLDを使用
_f0, t = pw.dio(wav, fs)
f0 = pw.stonemask(wav, _f0, t, fs)
sp = pw.cheaptrick(wav, f0, t, fs)
ap = pw.d4c(wav, f0, t, fs)
#メルケプストラムに変換前のスペクトル包絡を可視化
plt.figure()
plt.imshow(np.log10(sp).T, origin='lower')#対数を取ることで見やすくなる
plt.colorbar()
plt.title('spectral envelope')
#メルケプストラムの算出及び可視化
mcep = sptk.sp2mc(sp, order=25, alpha=0.46)
plt.figure()
#スペクトル包絡に合わせてアスペクト比の調節
plt.imshow(mcep.T, aspect=mcep.shape[0]/mcep.shape[1]/3.7, origin='lower')
plt.colorbar()
plt.title('Mel-cepstrum coefficients')
#メルケプストラムからスペクトル包絡に変換及び可視化
sp_from_mcep = sptk.mc2sp(mcep, alpha = 0.46, fftlen = 1024)
plt.figure()
plt.imshow(np.log10(sp_from_mcep).T, origin='lower')
plt.colorbar()
plt.title('spectral envelope from mcep')
#WORLDでそれぞれ音声波形を合成し，保存
synthesized = pw.synthesize(f0, sp, ap, fs)
librosa.output.write_wav('./world.wav', synthesized, fs)
syn_mcep = pw.synthesize(f0, sp_from_mcep, ap, fs)
librosa.output.write_wav('./world_mcep.wav', syn_mcep, fs)
##出力
###WORLDから算出したスペクトル包絡
###メルケプストラム
###メルケプストラムから変換したスペクトル包絡
###出力を見て
- メルケプストラムは，スペクトル包絡と比べてかなり情報圧縮されているのがわかる．
- メルケプストラムは，この情報量で音声波形の合成ができると考えると凄い
- メルケプストラムから変換したスペクトル包絡は，少し平滑化されたようになっており，当たり前だが完全な復元は厳しい
- 音声を聞いたところ，元のスペクトル包絡で合成したものより，心なしか全体的に歪んで聞こえる
#まとめ
メルケプストラムは低次の音声の特徴量であり，ニューラルボコーダや声質変換で使われる．Pythonでは音声処理のラッパーも充実しているため，簡単なコードで求められるのが良い．しかし，音声の特徴量は種類が多いため，原理まで理解するのは大変である．なので，今後は原理的な部分の理解も行っていきたい．