ある講義で紹介された「フーリエ級数」
どんな周期関数も三角関数の無限の足し算で表せる！
そして、別の講義で紹介された「最小二乗法」
サンプルを、できるだけ誤差が小さくなるように多項式で表せる！

この2個のツールを紹介された俺は、あることを思いつきました。
…そうだ！フーリエ級数の式にwavファイルの波形データを最小二乗法で当てはめれば、
名状しがたい不可逆圧縮のようなものができるのじゃないか？

実際にやってみた(トリビアの泉風)

紹介されたフーリエ級数はA1*sin(2πx/p)+B1*cos(2πx/p)+A2*sin(4πx/p)+B2*cos(4πx/p)+…の形(pは周期)だったが、
実験してみるとf(x)=1が表せないっぽいことに気づき、定数項を加えることにした。
そして、これだけだと掛ける値が2,4,6,8,...なのか2,4,8,16,...なのかわからなかったのでネットで調べたところ、
前者らしいことがわかった。

wavファイル全体をフーリエ級数で表そうとすると多量の項が必要であることがわかったので、
wavファイルを適当なサンプル数ごとに区切り、フーリエ級数の式にあてはめることにした。
(もともとフーリエ級数は無限個の項が必要らしいし…)

とりあえず、結果をご覧頂こう。
[youtube][/youtube]
項数が64のあたりから聞き取れるようになり、項数が128くらいでだいぶましになってくる。
しかし、今のサンプル数は512であり、512サンプルを項数64で表現しているのである。
項数64となっているが、ここに表示している項数はAn*sin(2*n*πx/p)+Bn*cos(2*n*πx/p)という組の個数のことであり、
定数項もあるので、実際の最小二乗法のパラメータの数は「項数」*2+1、すなわち129個である。
このパラメータ1個を4バイトのfloat型で表すと、516バイトになる。
一方、16ビットモノラルのwavファイルの512サンプルは1024バイトである。
よって、それなりの音質で約半分のデータ容量にすることができると言える。
ただし、プププププ…というノイズが入ってしまうのが難点である。これを消すのは今後の課題とする。

また、最小二乗法をかけるサンプル数(ブロックサイズ)を変えてみると、このようになる。
[youtube][/youtube]
音質は、だいたい項数÷ブロックサイズと相関があることがわかる。(比例しているかはわからない)
ブロックサイズと項数をともに2倍にすると、ノイズ以外の音声はほとんど変わらずに、ノイズの間隔を広くすることができるようだ。

ここまで音質のみを気にしてきたが、処理時間も気にしたほうがいい。
最小二乗法の時間計算量は、サンプル数をN、推定するパラメータ数をMとすると、
だいたいO((N+M)*M^2)であることがわかる。
例えば1回の処理でN=1024、M=256とすると、(N+M)*M^2==83,886,080となり、だいたい8秒程度かかると見積もれる。
実測してみると、自分の環境では1回の測定で約3.7秒であった。
これをwavファイルの長さだけ繰り返すので、かなり時間がかかることがわかる。

この処理は実用的な不可逆圧縮には向いていなそうであるが、簡単なものを自作して遊ぶぶんには
ある程度の成果が得られそうだ、ということがわかった。

おまけ
この実験で扱ったのは「フーリエ級数」である。
そういえば、波形を処理する方法として「離散フーリエ変換」というものがあった。
同じ「フーリエ」という言葉が入っているし、何か関係あるかも？
そういえば、離散フーリエ変換のコードでもsin(2*n*πx/p)とかcos(2*n*πx/p)に相当するものがあったな…
ということで、とりあえず同じデータを使い、それぞれの変換結果のグラフを描いてみた。 …なるほど、わからん。
よく考えたら、同じ文字列が入っていても「ガウスの消去法」と「ガウス記号」のように、関係ないものもありますね。

最後に、実験に使用したプログラムとデータを用意しておきます。
renhenkan.plがあるディレクトリをカレントディレクトリにし、 というコマンドで最小二乗法→逆変換ができ、
逆変換の結果が「出力フォルダ」の中の"riv_result.wav"に格納されます。

みけCAT さんが書きました: 紹介されたフーリエ級数はA1*sin(2πx/p)+B1*cos(2πx/p)+A2*sin(4πx/p)+B2*cos(4πx/p)+…の形(pは周期)だったが、
実験してみるとf(x)=1が表せないっぽいことに気づき、定数項を加えることにした。

一般のFourier級数は
f(x) = a_0/2 + Σ(k=1..∞){a_k cos(kx) + b_k sin(kx)}
と表したもののことを指します。上の形のものはf(0)=0や周期の条件等をおいたものでしょう。
(cf. フーリエ級数 - wikipedia)

みけCAT さんが書きました: この実験で扱ったのは「フーリエ級数」である。
そういえば、波形を処理する方法として「離散フーリエ変換」というものがあった。
同じ「フーリエ」という言葉が入っているし、何か関係あるかも？

離散Fourier変換とは、離散的な(連続でない)函数をFourier変換するアルゴリズムのことで、これは高速で計算できることが知られています。
Fourier級数は、函数をFourier変換したときに得られる級数のことを指しますから、関係はもちろんありますね。
(cf. 離散フーリエ級数 - wikipedia)

日記の内容だけではどうやって変換をかけたのかがあまりよく分からなかったのですが、最小二乗法で誤差を小さくするという考え方は面白いと思いました。私は圧縮や音声加工の技術には詳しくないので分かりませんが、ケースバイケースで役に立つこともあるんでしょうか。
単純なFourier級数展開をして項を切り捨てて圧縮をするときと比べてどれくらい差が出るのかも調べてみると面白いかもしれないですね〜。