遺伝的アルゴリズムで、YouTube利用者が最も好む音声を生成してみた【matlab】

# 0.遺伝的アルゴリズムとは

1. 適当に数十個データを作って、
2. よかったやつ上位2つを選び
3. それらを適当に混ぜ合わせたり
4. 突然変異させたりする

の繰り返しです。手順2~4の結果、数十個のデータが再び出来上がるので、手順1に戻るわけです。

# 1.アイデア
本研究では、データをYouTube上で公開して、手順2(よかったやつ上位2つを選ぶ)をYouTube上のコメント欄における投票とします。
専門用語でいうなら、「評価関数」に人の主観的評価を使うわけです。

# 2.各手順の詳細説明

## 2-1. 手順1(適当に数十個のデータを作る)

音声はそれぞれ5秒とします。

### 2-1-1. フォーマットはwavがよい
mp3など、圧縮のかけられた音声フォーマットは、ファイル構造が複雑であるため乱数による「ランダムな音声」の生成に不向きです。
そこで、波形(をDA変換した)データをそのままの形で保有するフォーマットであるwavを利用します。

wavファイルは波形データ以外の情報(ファイルの大きさや、ファイルがwavデータであることを宣言する文字列など)を書いておく部分(メタデータと呼ぶことにします)と、波形データが置かれている部分(メインデータと呼ぶことにします)に分かれます。

メタデータを自分で作るのは面倒なので、ソフトで正弦波wavを生成して、そのメタデータをそのまま使ってしまいましょう。
正弦波にこだわる理由は、「繰り返し現れるデータ⇔メインデータ」という命題が成り立つため、メタデータとメインデータの分離を安全に行えるからです。

### 2-1-2. 正弦波wavの生成

まず、WavePad(無料版)をインストールし、Tone Generatorを開き、5000ms, 440Hzの正弦波を作ってWaveデータに書き出しました。(拡張子はwavです)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/319477/9d586d76-1c67-df85-ab6c-a676564dec94.png)
図1 Tone Generator上の設定

このデータは[GitHub](https://github.com/17ec084/genetic-algorithm/blob/master/5000ms_440hz_sine.wav)で公開しています。

再生すると、音叉のような音が聞こえます。(ブラウザ上で再生されない場合、ダウンロードしてみてください)

### 2-1-3. バイナリエディタによるメタデータとメインデータの分離

次にバイナリエディタ(stirling)にてこのwavデータを開きました。
[wavデータの構造を解説しているホームページ](http://home.a00.itscom.net/hatada/asp/wav.html)と見比べながら、本当に波形にかかわっているのはどの部分か特定しました。
要は'data'の文字列(16進表示で`64	61	74	61`)の先頭から8バイト先以降が波形データのようです。
つまり、図2で黒く選択した部分(441000 Byte)です。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/319477/8a0bacc0-f177-aaa4-5ef4-fa56e4a54bd4.png)
図2 黒塗りのメインデータと、塗っていないメタデータ

今、波形データの先頭27byteである
`FF 07 FC 0F E9 17 BE 1F 73 27 01 2F 5F 36 86 3D 70 44 14 4B 6D 51 74 57 24 5D 75`
をデータの中から検索していきました。前半でのみ、何件か見つかりました。後半では1件も見つかりませんでした。
このような27byteのデータが、別の箇所で偶然に出てくることは考えられません※。
したがって、この27byteがくり返し出てきたのは、**「同じデータを繰り返している」ことをほぼ確実に示しています**ので、
おそらくちゃんと正弦波の繰り返しを表現してくれているのでしょう。(後半では見つからなかったのは、サンプルレートが44100Hzであるため440Hzでは割り切れず、このことが「細かいずれ」が生じさせたものと考えられます)
<sub><sub>
※1bitを「2通り」、1byteを「256通り」と表現するとき、27byteはおよそ「$10^{65}$通り」となります。
つまり、これと全く同じ27byteのデータが偶然作られる確率は$10^{-65}$未満です。
今回生成したwavファイルは441046byteでした。
441046-27=441019回の試行にて、確率$10^{-65}$のことが1度以上起きる確率は
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/319477/29b15fcc-dacd-60bc-634d-4be0c280b328.png)
≃$4.41019 × 10^{-60}$
となるそうです。
参考:http://www.math.kobe-u.ac.jp/HOME/saji/mathyomi/probability.html
偶然ではまず起きないような確率ですよね。</sub></sub>

したがって、図2で黒く示したメインデータの部分(ちょうど441000 Byte)を書き換えていけば、Waveファイルを動的に生成することができます。

GitHub上で、[メタデータ](https://github.com/17ec084/genetic-algorithm/blob/master/5000ms_meta)と[メインデータ](https://github.com/17ec084/genetic-algorithm/blob/master/5000ms_main)に分離したwavファイルを公開しています。

### 2-1-4. メインデータの構造
結論から先に言うと、

- メインデータは1つ2byteの標本の塊である
- 1つの標本は1/44100秒間の音を表現する
- 1つの標本はリトルエンディアン(「128～1の位」→「-32768～256の位」の順)で記録されている

です。

5秒間の音声を441000 Byteで表現しているため、440Hz正弦波1周期、すなわち(1/440)秒間の音声は200 byteほどで表現されているはずです。
一方サンプリングレートは44100Hzなので、1秒間に44100回、(1/440)秒間に100回ほどサンプリングされているはずです。
つまり、100個の標本が200byteで表されていることになりますから、1個の標本は2byteとなっているはずです。
実際、波形データ先頭400byteを2byteずつ読み込んでプロットしてみると、図3のように波数2の正弦波を確かに得ることができます。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/319477/beb7d085-f021-aa52-b554-0b8d62752c5b.png)
図3 音声データ(先頭およそ4.5ms)

図3は、400バイト分だけバイナリエディタから「data.txt」に16進数をコピペしたのち、matlabで次のコードを実行するだけで得られます。

```matlab
bytes = fscanf(fopen('data.txt','r'), '%x'); %符号なし16進数2桁(0～255)として読む

for i=[1:2:400]
    idx = floor((i-1)/2)+1; %4桁の16進数を合成する「two_bytes」用のインデックス 
    two_bytes(idx) = bytes(i+1)*256+bytes(i); %順番が直感と逆なので注意(「128～1の位」→「256～32768の位」)

    if two_bytes(idx) >= (256^2)/2
        two_bytes(idx) = two_bytes(idx) - (256^2); %符号あり形式に変換 
    end
end

scatter([1:200], two_bytes);
```

### 2-1-5. メインデータの書き換え
メインデータをランダムに置き換えていきましょう。
とはいえ、完全なランダムは避けた方がいいでしょう。
完全なランダムですと、1秒に44100回も「音の傾向」が変化してしまうこととなります。つまり1秒間の中に44100個の確率変数が入っていることになります。すると統計学的には「巨視的に」音を聞くことになっているわけです。
巨視的にみる、聞くと、何らかの傾向が強く表れてしまいます。
たとえ話でいうなら、サイコロを44100回も振ったら、どの目が出る回数も大体44100÷6 回となり、**違いが現れなくなってしまいます**よね。これではランダムにした意味が失われてしまいます。

そこで、0.25秒間ずつ、同じ正弦波による出力をさせてみたいと思います。

$220\times 1.054766076^{scale(0以上39以下の一様分布する整数)}　　　　ヘルツ$の正弦波0.25秒,音量amp(0～32767)をベクトルで出力してくれるmatlab関数は次の通りです。

```matlab
function [vector, angle] = generate_sine(scale, amp, angle)
    vector = [];
    sample_rate = 44100;
    hz = round(220*(1.054766076^scale));
    sec = 0.25;
    sample_cnt = round(sample_rate*sec); 
    each_wave_sample_cnt = sample_rate/hz; %1周期の標本数
    while(size(vector,2)<sample_cnt)
        new_vec = amp*sin(angle+linspace(0, 2*pi, each_wave_sample_cnt));%ここでsinを変えればいろんな音色に置き換え可
        new_vec(end) = [];
        vector = [vector, new_vec]; 
    end
    vector = round(vector(1:sample_cnt));
    angle = mod(angle + 2*pi*(mod(sample_cnt, each_wave_sample_cnt-1))/each_wave_sample_cnt, 2*pi);
    %位相ずれ情報も返却
    return;
end  
```

ボリューム常に最大、音程は0.25秒ごとに4段階上がるのを5回繰り返すことで5秒間の音源を作るとしたら、次のようなコードで出来ます。

```matlab
scales =        [4,8,16,32,4,8,16,32,4,8,16,32,4,8,16,32,4,8,16,32];
volume = 32767.*[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1];

vector = [];
angle=0;

for i=[1:size(scales,2)]
    [new_vector, angle] = generate_sine(scales(i), volume(i), angle);
    vector = [vector, new_vector];
end

plot(vector);
```

もしくは

```matlab
function vector = generate_waves(scales, volumes)
    vector = [];
    angle = 0;
    for i=[1:size(scales,2)]
        [new_vector, angle] = generate_sine(scales(i), volumes(i), angle);
        vector = [vector, new_vector];
    end
    return;
end 
```
としてから`plot(generate_waves([4,8,16,32,4,8,16,32,4,8,16,32,4,8,16,32,4,8,16,32], 32767.*[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]));`とします。


プロットし、一部に注目したものが図4です。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/319477/5f381132-f1b3-d1ba-bc54-a73411f37560.png)
図4. テスト波形の一部

最後に、このファイル出力機能を使って、wavファイルを出力する関数も書きます。

```matlab
function wav_write_5000ms(filename, scales, volumes)
    fwrite(fopen(filename,'w'),fread(fopen('5000ms_meta')));
    fwrite(fopen(filename,'a'),generate_waves(scales, volumes),'integer*2');
end
```

いま、
`wav_write_5000ms('test.wav', [4,8,16,32,4,8,16,32,4,8,16,32,4,8,16,32,4,8,16,32], 32767.*[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]);`
を実行すると、[test.wav](https://github.com/17ec084/genetic-algorithm/blob/master/test.wav)が生成されました。
(電話の着信音みたい笑と母親が申しておりました。)

実際には、「音の重ね合わせ」にも対応したいと考え、3つの正弦波を重ね合わせる`wav_write_3_superposition_5000ms`関数も作りました。
[GitHubのレポジトリ](https://github.com/17ec084/genetic-algorithm)で公開しています。

### 2-1-6. ランダム音声の生成
3つの正弦波を重ね合わせてみましょう。
`volumes1~3`は0～32767の整数一様乱数、`scale1~3`は0～39の整数一様乱数とします。
100個のデータを作る場合でも、次のようにするだけです！

```matlab
for i=1:100
    wav_write_3_superposition_5000ms...
    (...
        strcat('data', num2str(i), '.wav'),...
        floor(40*rand(1,20)),...
        floor(32768*rand(1,20)),...
        floor(40*rand(1,20)),...
        floor(32768*rand(1,20)),...
        floor(40*rand(1,20)),...
        floor(32768*rand(1,20))...
    );
end
```

## 2-2. 手順2(よかったやつ上位2位)
2-1-6節のコードで生成した100個の第1世代音声データを[YouTubeで公開しました](https://www.youtube.com/watch?v=krZR-9FvvrA)！
Qiitaのこの記事( https://qiita.com/17ec084/items/47743da791811820eb22 )かYouTubeのコメント欄に、**1番いいと思った音声1つ**をコメントしてください。
**ご協力お願いいたします！**

## 2-3. 手順3(上位2つを混ぜ合わせる)

Qiitaの記事およびYouTubeに投票が寄せられ、上位2つが決定可能になったら、その時点から24時間待ちます。
**24時間待った時点においても上位2つが決定可能かつ、上位2つは24時間前と同じものであったなら、そこで投票を締め切ります**。
例えば4月1日17:00の時点で「**データ5に2票、データ61に2票、データ20に1票、それ以外は1票も入っていなかった**」とします。
4月2日17:00の時点で

- 「**データ5に2票、データ61に2票、データ20に2票**」となっていた場合、上位2つを決定できないため、さらに待ちます。
- 「**データ5に2票、データ61に3票、データ20に3票**」となっていた場合、この時点で上位2つが決定可能ですが、4月1日における上位2つと異なるものですので、この場合もさらに待ちます。
- 「**データ5に3票、データ61に2票、データ20に1票**」となっていた場合、投票を締め切り、データ5と61を「上位2つ」とみなします。

その他、明らかなスパム行為による票で締め切りまたは投票そのものを妨害するものなどが発生した場合は、臨機応変に対応します。

このようにして選ばれた上位2つについて、「交叉」させていきます。

2-1-6節に示したコードからわかる通り、本研究における音声データは20次元ベクトル6つで表現できます。
これを120次元ベクトルに変換します。
上位2つのデータについて、2点交叉させたものを24個、一様交差させたものを24個、計48個の120次元ベクトルを新たに生成します。

2点交叉、一様交差については、[北海学園大学のホームページ](http://www.eli.hokkai-s-u.ac.jp/~kikuchi/ma2/chap15.html)が参考になりました。(URLが移転したようなので、もしかしたらそのうちサーバごと止まるのかもしれません。でも[web魚拓を発見したので一安心です](http://web.archive.org/web/20170922113434/http://www.eli.hokkai-s-u.ac.jp:80/~kikuchi/ma2/chap15.html)。)

2点交叉とは、適当に次元を2つ選んで、その間を、上位2つどうしで交換することです。
例えば上位2つのデータが`[1,5,2,9,4]`および`[6,3,5,8,9]`だったとします。(本当は120次元ですがここでは簡単に5次元で書いています。)
2点として2次元目と4次元目を選んだ場合、`[1,3,5,8,4]`および`[6,5,2,9,9]`が作られます。

一様交差とは各次元を半分の確率で入れ替えることです。
例えば、2次元目および5次元目が入れ替わる次元として選ばれたとき、
`[1,5,2,9,4]`と`[6,3,5,8,9]`からは
`[1,3,2,9,9]`および`[6,5,5,8,4]`が得られます。

## 2-4. 手順4(突然変異を起こさせる)
上位2つそのまま、そして2-2節で作った48個、合わせて50個の第2世代データが作られました。残り50個は突然変異で作っていきます。
突然変異とは、任意の次元が乱数に置き換わることです。(2つの次元を入れ替える場合もあるそうですが今回は無視します)
21～36つの次元を突然変異させるのが10個(上位2つから5個ずつ×2)
37～52つの次元を突然変異させるのが10個
53～68つの次元を突然変異させるのが10個
69～84つの次元を突然変異させるのが10個
85～100つの次元を突然変異させるのが10個
としたいと思います。

こうして、第2世代の音声データ100個が完成します。

以降はこれを繰り返していくわけです。