【Pythonで音声解析!?】PyWorld, pyreaper入門

音声解析が必要な時があったのでPythonでやってみました.

ちょこちょこ閲覧数も増え, いいねも頂いて有難いのですが, **私は音声解析専門ではありませんのでその辺ご理解ください**.

# wavデータを作成する
どちらのライブラリを使うにせよ, 基本的にwavデータを使うことになります.
その際, ステレオではなく, モノラルでないといけません. 一応, ターミナルでmp4からwavを作成する方法を載せておきます.
もしffmpegがinstallされていなければ次のコマンドを打ってinstall出来ます(brewが必要です).

```shell:Installの方法
$ brew install ffmpeg
```

変換は次のようにして行います.

```shell:実行方法
$ ffmpeg -y -i file.mp4 -ac 1 -ar 44100 file.wav
```

`-i` は入力ファイルの指定です.
`-y` は常に上書きするという指定です.
`-ar 44100` はサンプリング周波数の指定です.
wavは**モノラル形式**でなければなりません. それを `-ac 1` で指定しています.

# dir構造
次のようにdirectoryを構成しています.

```shell
voice_analysis
    ├── data
    │   └── sample.wav
    └── src
        └── voice_analysis.ipynb
```


# PyWorld

## PyWorldとは... ?
C++, matlab用に作られた音声解析ライブラリ[World](https://github.com/mmorise/World)のPython用のwrapperみたいです.

## PyWorldのインストール
[githubのREADME](https://github.com/JeremyCCHsu/Python-Wrapper-for-World-Vocoder/blob/master/README.md)にしたがってinstallすればいいだけです.
と思いきや, installできませんでした.
どうやらPyWorldはanacondaと相性が悪いみたいです...
[Does not work in Anaconda #12](https://github.com/JeremyCCHsu/Python-Wrapper-for-World-Vocoder/issues/12)

自分はpyenvで環境構築しました. 本当はanacondaのpythonで構築しようと思いましたが, なぜか上手く行かなかったので...おそらくはpyenvではなく, systemのpythonでも上手く行くと思います. Anacondaでなければ大丈夫のはずです.

## pyenvでのpython3の導入
brewはinstall済みを想定しています.

参考 : [Mac環境へのPython3系インストール](https://qiita.com/spyc/items/73d1295f8b3dde3b49ca)

ただ, これ通りにやると最新のMojaveではzlibエラーが出ます. どうやらxcodeが追いついていないようです.
参考 : [MacOS Mojave pyenvでpythonのインストールがzlibエラーで失敗した時の対応](https://qiita.com/zreactor/items/c3fd04417e0d61af0afe)

あとは.bashrcなり.zshrcなりにPATHを通せばokです.

```zsh:.zshrc
# anacondaのPATHは通さないようにコメントアウト
# export PATH="/anaconda/bin:$PATH"
# export PATH=/Users/username/Library/Python/2.7/bin:$PATH

export PYENV_ROOT="$HOME/.pyenv"
export PATH="$PYENV_ROOT/bin:$PATH"
eval "$(pyenv init -)"
```

ちゃんとpythonはpyenvの方を示すのか確認してみましょう.

```shell:pythonの環境の確認
$ which python3
/Users/username/.pyenv/shims/python3

$ which pip
/Users/username/.pyenv/shims/pip
```

大丈夫そうです!

## PyWorldのinstall

```shell:scipyとpyworldのinstall
$ pip install scipy
$ pip install pyworld
```
scipyを先にinstallしないとinstall出来ない不具合もあるようなので先にinstallしておきます.

## PyWorld入門
jupyter notebookで試してみましょう.

参考 : [音声合成システム WORLD に触れてみる](https://qiita.com/ohtaman/items/84426cee09c2ba4abc22)

今回はsample2として段階的に高さをあげていった音声を自作してみました.

```python:hello_pyworld.ipynb
from scipy.io import wavfile
import pyworld as pw
import numpy as np
import matplotlib.pyplot as plt

WAV_FILE = "../data/sample2.wav"

# fs : sampling frequency, 音楽業界では44,100Hz
# data : arrayの音声データが入る 
fs, data = wavfile.read(WAV_FILE)

# floatでないとworldは扱えない
data = data.astype(np.float)

_f0, _time = pw.dio(data, fs)    # 基本周波数の抽出
f0 = pw.stonemask(data, _f0, _time, fs)  # 基本周波数の修正

# --- これは音声の合成に用いる(今回は使わない)
# sp = pw.cheaptrick(data, f0, _time, fs)  # スペクトル包絡の抽出
# ap = pw.d4c(data, f0, _time, fs)         # 非周期性指標の抽出
# y = pw.synthesize(f0, sp, ap, fs)    # 合成


# 可視化
plt.plot(data, label="Raw Data")
plt.legend(fontsize=10)
plt.show();

plt.plot(f0, linewidth=3, color="green", label="F0 contour")
plt.legend(fontsize=10)
plt.show();
```

これを入力すると... ちゃんと綺麗にピッチが取れました！
___
- Raw Data

読み込んだデータそのものを可視化したものです.
「あー, あー, あー, あー」と発声しているのがわかります.

![pyworld_rpw.png](https://qiita-image-store.s3.amazonaws.com/0/324488/4fe6f5eb-160d-af44-b342-2ceeafba5493.png)


___
- F0 contour

これはピッチです. 綺麗に取得出来ていそうです.

![pyworld_f0.png](https://qiita-image-store.s3.amazonaws.com/0/324488/1eeb297b-d898-9e78-74bf-242610f7d037.png)




## PyWorldの計算時間
果たしてPyWorldの計算時間はどうなっているのでしょうか.
ということでWAVファイルの時間と解析時間について簡単に可視化してみました.

```python:calc_pyworld.ipynb
from scipy.io import wavfile
import pyworld as pw
import numpy as np
import matplotlib.pyplot as plt
import time

def calc_time(file_name):
    start = time.time()
    
    fs, data = wavfile.read(WAV_FILE)

    data = data.astype(np.float)

    _f0, _time = pw.dio(data, fs)   
    f0 = pw.stonemask(data, _f0, _time, fs)

    elapsed_time = time.time() - start
    print ("elapsed_time:{0}".format(elapsed_time) + "[sec]")
    
    return elapsed_time

def plot_calc_time():
    file_names = os.listdir("../data")
    
    plt.figure()
    t_list = []
    elapsed_time_list = []
    for file_name in file_names:
        t = int(file_name[:-4])
        elapsed_time = calc_time(file_name)
        
        t_list.append(t)
        elapsed_time_list.append(elapsed_time)
    
    t_list.sort()
    elapsed_time_list.sort()
    plt.plot(t_list, elapsed_time_list, marker=".")
    plt.title("Time for calculating f0")
    plt.xlabel("wav_time")
    plt.ylabel("calc_time")
    plt.grid()
    plt.show();

plot_calc_time()
```

file_nameに関しては各々のディレクトリ構造に依存するので適宜修正が必要です.
結果はこんな感じでした.

![pyworld_time.png](https://qiita-image-store.s3.amazonaws.com/0/324488/f2bfcf95-7752-28dc-8b0c-036a916d8344.png)

どちらも秒単位. よくわからない挙動を示しているが一応線形でしょうか?
あまり秒数に依らずに高速に解析できています.

# pyreaper

## pyreaperとは... ?
Googleがgithubにあげている[REAPER](https://github.com/google/REAPER)のPython用のラッパーです.
[A python wrapper for REAPER](https://github.com/r9y9/pyreaper)

## pyreaperのinstall
今回もREADMEにしたがってinstallするだけです.

```shell:pyreaperのinstall
$ pip install pyreaper
```

これでinstall...と思いましたがinstallできません！

>Command "/Users/username/.pyenv/versions/3.7.1/bin/python3.7 -u -c "import setuptools, tokenize;__file__='/private/var/folders/gw/xyk2qtxd61lg89vs0ckk0wkc0000gn/T/pip-install-ury00gf5/pyreaper/setup.py';f=getattr(tokenize, 'open', open)(__file__);code=f.read().replace('\r\n', '\n');f.close();exec(compile(code, __file__, 'exec'))" install --record /private/var/folders/gw/xyk2qtxd61lg89vs0ckk0wkc0000gn/T/pip-record-8wj1vx31/install-record.txt --single-version-externally-managed --compile" failed with error code 1 in /private/var/folders/gw/xyk2qtxd61lg89vs0ckk0wkc0000gn/T/pip-install-ury00gf5/pyreaper/


よくわからないエラーが出てしまいました.
~~諸悪の根源~~anacondaは抹殺したはずなんですが...

どうやらpyreaperの依存関連のエラーのようです.
次を実行して無事installしましょう.

```shell:pyreaperのinstall
$ git clone https://github.com/r9y9/pyreaper
$ cd pyreaper
$ git submodule update --init --recursive
$ python setup.py develop
```

次にpysptkのinstallです.

```shell:pysptkのinstall
$ pip install pysptk
```

## pyreaper入門
jupyter notebookで試してみましょう.

参考 : [How pyreaper works](http://nbviewer.jupyter.org/gist/r9y9/e7d6604757403ec74557ccd93e1f0c8c)


```python:hello_pyreaper.ipynb
%pylab inline
rcParams["figure.figsize"] = (16, 4)

%matplotlib inline
import pysptk
import pyreaper
import matplotlib.pyplot as plt
from scipy.io import wavfile


file = "../data/sample2.wav"
fs, data = wavfile.read(file)

plt.plot(data, label="Raw Data")
plt.legend(fontsize=10)
plt.show();

pm_times, pm, f0_times, f0, corr = pyreaper.reaper(x, fs)

plt.plot(pm_times, pm, linewidth=3, color="red", label="Pitch mark")
plt.legend(fontsize=10)
plt.show();

plt.plot(f0_times, f0, linewidth=3, color="green", label="F0 contour")
plt.legend(fontsize=10)
plt.show();

plt.plot(f0_times, corr, linewidth=3, color="blue", label="Correlations")
plt.legend(fontsize=10)
plt.show();
```
先のwavファイルを入力すると... ちゃんと綺麗にピッチが取れていました！
___
- Raw Data

読み込んだデータそのものを可視化したものです.
「あー, あー, あー, あー」と発声しているのがわかります.

![reaper_1.png](https://qiita-image-store.s3.amazonaws.com/0/279783/13f70ed1-2c32-e868-87bc-bf6d523f8cc0.png)
___
- Pitch Mark

よく理解出来ていません. 音があるときに値を持っていますし, ピッチが検出できているって印なのでしょうか...?

![reaper_2.png](https://qiita-image-store.s3.amazonaws.com/0/279783/2a2f1f68-2d11-ba40-ec13-6831253f3dd5.png)
___
- F0 contour

これはピッチでしょう. 綺麗に取得出来ていそうです.

![reaper_3.png](https://qiita-image-store.s3.amazonaws.com/0/279783/e4487091-6f34-9074-52c6-a86ab4541d05.png)
___
- Correlations

これもまだあまり理解出来ていません. 自己相関の値でしょうか.

おそらくピッチの計算として, フーリエ変換した後にピークを検出する方法を用いていて, そのピーク検出に自己相関関数が用いているということでしょうか.

![reaper_4.png](https://qiita-image-store.s3.amazonaws.com/0/279783/dbd01465-ebf0-c777-a429-1cd09aba8326.png)

## pyreaperの計算時間
果たしてpyreaperの計算時間はどうなっているのでしょう.
ということでWAVファイルの時間と解析時間について簡単に可視化してみました.

```python:calc_pyreaper.ipynb
import pysptk
import pyreaper
import sys,os
import numpy as np
import matplotlib.pyplot as plt
from scipy.io import wavfile


def calc_time(file_name):
    start = time.time()
    
    fs, data = wavfile.read(file_name)

    pm_times, pm, f0_times, f0, corr = pyreaper.reaper(data, fs)

    elapsed_time = time.time() - start
    print ("elapsed_time:{0}".format(elapsed_time) + "[sec]")
    
    return elapsed_time

def plot_calc_time():
    file_names = os.listdir("../data")
    
    plt.figure()
    t_list = []
    elapsed_time_list = []
    for file_name in file_names:
        t = int(file_name[:-4])
        elapsed_time = calc_time(file_name)
        
        t_list.append(t)
        elapsed_time_list.append(elapsed_time)
    
    t_list.sort()
    elapsed_time_list.sort()
    plt.plot(t_list, elapsed_time_list, marker=".")
    plt.title("Time for calculating f0")
    plt.xlabel("wav_time")
    plt.ylabel("calc_time")
    plt.grid()
    plt.show();

plot_calc_time()
```

file_nameに関しては各々のディレクトリ構造に依存するので適宜修正が必要です.
結果はこんな感じでした.

![calc_pyreaper.png](https://qiita-image-store.s3.amazonaws.com/0/279783/b3c34ada-e7f6-a2f2-dc29-e8a934fccc79.png)

どちらも秒単位. これは...計算時間は二次関数なのか指数関数なのか分かりかねます.
調べてみると, (高速)フーリエ変換は標本数 $N$に対して$O(N^2)$または,$O(N logN)$であるらしく, ピッチ解析の基本はFFTとピーク検出のはずなのでだいたいこれ $+\alpha$と言ったところでしょうか...?

参考 : [高速フーリエ変換 by Wikipedia](https://ja.wikipedia.org/wiki/%E9%AB%98%E9%80%9F%E3%83%95%E3%83%BC%E3%83%AA%E3%82%A8%E5%A4%89%E6%8F%9B)