読者です 読者をやめる 読者になる 読者になる

どせいたんさき。

ナスダヨー

『Change of Variable』 - フーリエ変換編

目的

Practical Statistics for Astronomers という書籍に以下のような演習問題がある.

Suppose that φ is uniformly distributed between zero and 2π. Find the distribution of sinφ. How could you find the distribution of a sum of sines of independent random angles?

Welcome to the Solutions site for Practical Statistics for Astronomers

これを解く.このエントリではフーリエ変換を用いて分布の足し合わせをおこなう.

解答

足し合わせる前の分布は (1/\pi)/\sqrt{1-x^2} という関数で表現することができる.乱数で生成した分布とあわせてみた.



この関数をフーリエ変換するとハンケル関数で表現できるらしい.とはいえハンケル関数の累乗を逆フーリエ変換するとかちょっと厳しいので数値的に変換する.計算は Octave でおこなった.

w  = linspace(-11,11,301)';  # horizontal axis
dw = diff(w)(1);             # step width

f = @(x) real((1/pi)./sqrt(1-x.*x)); # function

fx = shift(f(w),151)/trapz(w,f(w));  # shift &normalize
gx = fft(fx)*dw;                     # Fourier Transform

for k=1:6
  data{k} = \
  [w/k real(shift(ifft(gx.^k),150))*k/dw]; # Convolution
endfor

思ったような変換結果を得るためには配列をシフトしたりする必要があった.これは分布の中心を配列の先頭に持ってくるため.また,数値的な計算の誤差が無視できなかったので配列 fx を定義するときに規格化している.結果は以下のようになる.分布のピーク付近が微妙に合っていないのは数値的な誤差のせいではないかと思っているが確信はない.


参考資料