ガウス過程と非線形変換は，それぞれ様々な分野で活用されています．

また，その組み合わせもたくさんの応用先が考えられます．

この記事では，Reed-Solomon符号などの業績で有名なIrving S. Reedにより IEEE Transactions on Information Theory へ投稿された論文 On the use of Laguerre polynomials in treating the envelope and phase components of narrow-band Gaussian noise を紹介します．

この記事はD言語 Advent Calendarの18日目の記事です．

この記事はTUT Advent Calendar 2018の17日目の記事です．

この記事はD言語くんAdvent Calendarの16日目の記事です．

概要

D言語くんの生態に学術・産業界から注目が集まっている．

本報告では，D言語くんが支配する都市の存在について，またその都市における人々やD言語くんの生態について報告する．

世の中には専門家が書いた精度保証付きの高機能・高性能な数値積分ライブラリが多数あります．

しかし，実際はそこまで機能・性能を求めていなかったり，関数の形がキレイで積分が簡単だったりで，そのようなライブラリがオーバースペックな場合が結構あります．

しかも，そのようなライブラリは高度な技法を使っているため，専門的な知識がないとコードの理解ができないという問題があります．

でも，自分で台形公式やシンプソン則を書いてもあまり精度が出ないような場合もあります．

そのようなときに役に立つ二重指数関数型積分公式（Double Exponential: DE公式）を用いた数値積分ライブラリdeintを少し前に公開しましたので，利用方法を紹介します．

deintのGitHubリポジトリ

私が現在勉強しながら研究している非線形性に関する理論を簡単にまとめておきます．

もちろん，まだ一回も発表もしていない研究の内容をそのまま書くのはさすがにいけないので，ここではすでに世間ではよく知られているような内容のみをまとめておきます．

表記について

• \(j\)：虚数単位
• \(\mathcal{F}[\cdot]\)：フーリエ変換

\[S(f) = \mathcal{F}[s(t)] = \int_{-\infty}^{\infty} s(t) e^{-j2\pi f t} dt\]

ここで，\(S(f)\)を\(s(t)\)のスペクトル，\(\left\lvert S(f) \right\rvert^2\)を電力スペクトル密度とそれぞれ呼びます．

• \(\delta(x)\)：ディラックのデルタ関数

• \(\text{sinc}(x)=\frac{\sin \pi x}{\pi x}\)：sinc関数

正弦波のフーリエ級数展開

まず，フーリエ変換について復習しておきます． 正弦波は周期関数なので，正弦波のフーリエ変換（スペクトル）はフーリエ級数に展開したときの係数で表せます．

ここで，周期\(T_0 = \frac{1}{f_0}\)の任意の信号\(s(t)\)の複素フーリエ係数が\(c_n\)であれば，\(s(t)\)のスペクトルは次のように表せます．

\[S(f) = \sum_{n=-\infty}^{\infty} c_n \delta(f + n f_0)\]

そのため，まずは複素フーリエ係数\(c_n\)を求めてみます． ここで，\(c_n\)とは複素正弦波\(e^{j 2 \pi n f_0 t }\)と信号\(s(t)\)の内積で計算されます．

\[c_n = f_0 \int_{-\frac{1}{2 f_0}}^{\frac{1}{2 f_0}} s(t) e^{-j 2\pi nf_0 t} dt\]

\(s(t)=e^{j2\pi f_0 t}\)のとき，\(c_1=1\)で\(c_{n\not=1}=0\)になりますね？

では，話を戻して正弦波，つまり\(s(t)=\cos(2\pi f_0 t)\)のときにはどうなるでしょうか？ ここでは，オイラーの公式より

\[\cos(2\pi f_0 t) = \frac{1}{2} \left(e^{j2\pi f_0 t} + e^{-j2\pi f_0 t} \right)\]

になることから直ちに

\[c_n = \begin{cases} \frac{1}{2}, & n = 1 \\ \frac{1}{2},& n = -1 \\ 0, & \text{otherwise} \end{cases}\]

がわかります． つまり，正弦波のフーリエ変換はは次のようになります．

\[\mathcal{F}[\cos(2 \pi f_0 t)] = \frac{1}{2} \delta(f-f_0) + \frac{1}{2} \delta(f+f_0)\]

グラフで表すと以下の図のようになります． ここでは，ディラックのデルタ関数を表現するために矢印を用いて，その高さを係数\(c_n\)にしています．

矩形波のスペクトル

次に，よくある矩形波のスペクトルを考えてみます． ここでは周期\(T_0=\frac{1}{f_0}\)の次のような矩形波\(r(t)\)を考えます．

\[r(-\frac{1}{2f_0}<t<\frac{1}{2f_0}) = \begin{cases} 1, & \lvert 2 f_0 t \rvert < 0.5 \\ -1, & |2 f_0 t| > 0.5 \\ 0, & 2 f_0 t = 0, 1, -1 \\ \end{cases}\]

ただし，上の式の定義外については\(r(t)\)が周期\(T_0=\frac{1}{f_0}\)の周期関数となるように定義します．

この矩形波は\(0\)を中心に半分の区間では\(1\)に，残り半分の区間では\(-1\)になっています．

この矩形波のスペクトルを調べるために，まずは複素フーリエ係数を調べます． \(r(t)\)が偶関数であることに注目すると，

\[\begin{align} c_n &= f_0 \int_{-\frac{1}{2f_0}}^{\frac{1}{2f_0}} r(t) e^{-j2\pi nf_0 t} dt \\ &= 2 f_0 \int_{0}^{\frac{1}{2f_0}} r(t) \cos(2\pi n f_0 t) dt \\ &= 2 f_0 \left( \int_{0}^{\frac{1}{4f_0}} \cos(2\pi n f_0 t) dt - \int_{\frac{1}{4f_0}}^{\frac{1}{2f_0}} \cos(2\pi n f_0 t) dt \right) \end{align}\]

\(n\not=0\)においては

\[\begin{align} c_n &= 2f_0 \left(\left[\frac{\sin(2\pi n f_0 t)}{2\pi n f_0}\right]_{0}^{\frac{1}{4f_0}} - \left[\frac{\sin(2\pi n f_0 t)}{2\pi n f_0}\right]_{\frac{1}{4f_0}}^{\frac{1}{2f_0}} \right) \\ &= \frac{1}{\pi n} (2\sin(0.5 \pi n) - \sin(\pi n)) \\ &= \frac{\sin(\frac{\pi}{2} n)}{\frac{\pi}{2} n} = \text{sinc}(n/2) \end{align}\]

を得ます． また，\(n=0\)では\(c_0=0\)です． そのため，矩形波のスペクトルは

\[S(f) = \mathcal{F}[r(t)] = \sum_{n=-\infty,n\not=0}^{\infty} \text{sinc}(n/2) \delta(f - n f_0)\]

となります．

非線形変換を受けた正弦波

ここで，次のような非線形関数を考えます．

\[\begin{align} g(x) = \begin{cases} 1, & x > 0\\ -1, & x < 0\\ 0, & x= 0 \end{cases} \end{align}\]

この非線形関数\(g(x)\)は\(x\)が少しでも正に値を取れば\(1\)に，少しでも負に値を取れば\(-1\)に，\(0\)であれば\(0\)になるような関数です．

この関数は矩形波の定義に似ているかもしれませんが，周期関数ではありません．

この非線形関数を用いて矩形関数は次のように定義できます．

\[r(t) = g(\cos(2\pi f_0 t))\]

これを周波数領域で考えてみると，正弦波の2本の線スペクトルは\(g(x)\)な非線形特性を有するシステムを通過することで無限個の線スペクトルに化ける，と考えられます．

では，任意の非線形関数が与えられたとき，その非線形関数のどのような特徴が出力のスペクトルを決定づけるのでしょうか？

つまり，\(g(x)\)の特性だけがわかっているときに，正弦波を実際に適用することなくその出力のスペクトルを予想することはできるのでしょうか？

チェビシェフ多項式による展開

この問に関する答えは，意外かもしれませんが，チェビシェフ多項式という直交多項式にあります．

実は，任意の関数\(g(x)\)をチェビシェフ多項式\(T_k(x)\)で展開したときの係数と，\(g(\cos(2\pi f_0 t))\)のスペクトルにはある関係があります．

ここで，チェビシェフ多項式\(T_k(x)\)とは\(m\not=n\)において次の関係を満たすものです．

\[I_{m,n} = \int_{-1}^{1} T_m(x) T_n(x) \frac{1}{\sqrt{1-x^2}} dx = 0\]

具体的には

\[\begin{align} T_0(x) &= 1 \\ T_1(x) &= x\\ T_2(x) &= 2x^2 - 1\\ T_3(x) &= 4x^3 - 3x\\ T_4(x) &= 8x^4 - 8x^2 + 1 \end{align}\]

となっています． また，面白い特性として

\[T_n(\cos(x)) = \cos(nx)\]

となります．

このチェビシェフ多項式を用いて\(g(\cos(2\pi f_0 t))\)のスペクトルは次のように記述できます．

\[\mathcal{F}[g(\cos(2\pi f_0 t))] = \frac{a_0}{2} \delta(f) + \sum_{n=0}^{\infty} \frac{a_n}{2} \left(\delta(f-nf_0) + \delta(f+nf_0) \right)\] \[a_n = \frac{1}{\pi} \int_{-1}^{1} g(x) T_n(x) \frac{1}{\sqrt{1-x^2}} dx\]

なぜこんなことをするのでしょうか？

それは，\(a_n\)の積分で\(x=\cos(2\pi f_0 t)\)と置換積分することでわかります．

\(\frac{dx}{dt} = -2\pi f_0 \sin(2\pi f_0 t)\) で積分区間は\(\frac{1}{2f_0} \rightarrow 0\)です． また，\(dx/\sqrt{1-x^2} = -2\pi f_0 dt\)になるので，

\[\begin{align} a_n &= \frac{1}{\pi} \int_{\frac{1}{2f_0}}^{0} g(\cos(2\pi f_0 t)) T_n(\cos(2\pi f_0 t)) (-2 \pi f_0) dt \\ &= 2f_0 \int_{0}^{\frac{1}{2f_0}} g(\cos(2\pi f_0 t)) \cos(2\pi n f_0 t) dt \end{align}\]

となります． よく観察してみれば，これは複素フーリエ係数\(c_n\)そのものだということがわかります．

ん？じゃあこれは単に複素フーリエ係数を求める積分を変数変換しただけでは？と思うかもしれません．

そのとおりです！

ですが，もう少し考えていきましょう．

チェビシェフ多項式\(T_n(x)\)と係数\(a_n\)を用いることで，任意の非線形関数\(g(x)\)を表現することができます． つまり，

\[g(x) = \sum_{n=0}^{\infty} a_n T_n(x)\]

となる\(a_n\)があります． そしてこの\(a_n\)はフーリエ係数なのです．

何が嬉しいかというと， ある非線形な入出力特性をもつ素子や電気回路の特性をチェビシェフ多項式で近似したいとします．

この素子や回路の入力に正弦波を入れ，素子の出力をフーリエ変換すると，その係数が素子をチェビシェフ多項式で近似したときの係数になる，というのです！

さらに，逆も成り立ちます． チェビシェフ多項式を用いて特性が近似されたある素子に正弦波を入れたときの出力信号のフーリエ係数はチェビシェフ多項式の係数なのです．

最後に

以上までは，入力信号を正弦波にしたときの結果を述べました． ここではもう少し拡張してみて，入力信号が逆正弦分布という確率分布に従う場合を考えてみます． この逆正弦分布の確率密度関数は

\[p(x) = \frac{1}{\pi} \frac{1}{\sqrt{1-x^2}}\]

となります． \(\theta \sim U(0,2\pi)\)のとき\(\sin \theta\)，や\(\cos \theta\)は逆正弦分布に従います．

\(a_n\)を求めるときに用いた積分は，実は逆正弦分布に従う入力信号\(x\)について\(\mathbb{E}[g(x)T_n(x)]\)を求めていただなのです． そして，逆正弦分布に従う入力信号についてはチェビシェフ多項式が直交多項式となり，その係数はスペクトルと面白い関係を持ちます．