2 ポイント 投稿者 GN⁺ 2024-04-08 | 1件のコメント | WhatsAppで共有
  • 離散 Fourier 変換(DFT) は通信と信号処理の中核的なツールだが、周波数領域が現実を解釈する唯一の方法というわけではない
  • DCT のように入力サンプルに 基底関数 の値を掛けて周波数 bin を求める構造では、基底を変えるだけで、別の規則に基づく周波数領域を作れる
  • Walsh 行列+1-1 だけを使う方形波基底を提供し、sequency と直交性をそろえれば、時間領域と周波数表現を相互に行き来できる
  • Hadamard 行列は Walsh 行列を並べ替えた形で、Kronecker product やビット演算で作成し、行を sequency 基準で再ソートして WHT に活用する
  • 同じ入力でも、DCT では複数の高調波成分に分散し、Walsh-Hadamard 変換では方形波成分に分かれるため、DFT が「真実」を独占しているわけではないことが分かる

Fourier 周波数領域を見直す

  • 周波数領域は、複雑な信号を正弦波の振幅と位相に変換して表現する数学的空間である
  • この表現のおかげで、時間領域や空間領域で直接扱うのが難しい信号処理タスクを、より簡単に実行できる
  • DFT は通信と信号処理で中心的な役割を果たすが、方形波を奇数次の正弦波高調波の和に変換する解釈が、現実の唯一の解釈なのかは別問題である
  • 正弦波は自然界に広く存在するため Fourier 系のツールは多くのタスクに適しているが、別の規則で動作する きちんと定義された周波数領域 も作ることができる

DCT を基底関数として理解する

  • 離散コサイン変換(DCT)は、DFT を単純化した実数専用版と見なせる
  • DCT-II は、入力値 s_n に特定のコサイン式の値を掛けて合計し、特定の周波数 bin F_k の大きさを求める
  • 核心は、現在の DCT bin 番号に対応する周波数のコサイン波を作る 基底関数 である
  • 一般化すると、B(k, n)kn に応じて multiplier を返し、それを入力サンプルと掛け合わせて合計する構造になる
  • ソフトウェアの観点では B(k, n) をルックアップ配列として、数学的には 行列 として見ることができる
  • N=16 の DCT 基底行列では、最初の行 k=0 は 0Hz に相当する DC 成分であり、すべての値が +1.00 のコサインである
  • 以降の行は、半周期、1 周期、1.5 周期のように、次第に速く変化するコサインの形を持つ

方形波基底と Walsh 行列

  • 正弦波周波数ではなく 方形波 で信号を分解する基底関数は、Walsh 行列で作れる
  • Walsh 行列は、異なる速度で動く方形波で構成され、すべての multiplier は +1 または -1 である
  • 計算は、入力データの一部の符号を反転して合計する方式に単純化される
  • 単純に見える行列でも、2 つの条件を満たす必要がある
    • 各行は、前の行より符号反転が 1 つ多い sequency 順でなければならない
    • 時間領域データと周波数表現の間を滑らかに往復するには 直交性 を維持する必要がある
  • Walsh 行列を直接作るには N×N 配列から始め、N は 2 のべき乗でなければならない
    • 左端の最初の列に、すべての行の +1 を入れる
    • 新しい列を既存値のミラーコピーとして作り、新たに追加した領域を複数の水平区間に分けて、一部区間の符号を反転する
    • 反復ごとに列をコピーし、行区間の数を増やして交互に符号を反転する

Hadamard 行列から Walsh 配列を作る

  • 文献やオープンソースコードでは、Walsh 配列を直接作るより Hadamard 行列 から導出することが多い
  • Hadamard 行列は、Walsh 配列の行順を変えた形である
    • たとえば N=16 では、Walsh の row #15 は Hadamard では #1 に移動し、Walsh の row #1 は #8 に置かれる
  • Hadamard の構成が歴史的に先に登場し、Walsh がその上に構築されたことが、この慣例の理由の 1 つである
  • 実用上は Hadamard 行列を作る方法のほうがよく文書化されており、単純で効率的なビット操作方式もある
  • 教科書的な構成は、1×1 配列から始め、前の行列 H_{n-1} を 4 つのタイルにコピーする方式である
    • 左上、右上、左下はそのままコピーする
    • 右下はすべての符号を反転する
    • この拡張には Kronecker product 表記 が使われるが、実際の動作はコピーと符号反転である
  • 構成ステップを n 回経ると、Hadamard 行列のサイズは常に 2^n × 2^n になる
  • 特定セルの Hadamard 値は、x & y を計算した後、その結果で立っている bit 数が偶数か奇数かで求められる
    • 立っている bit 数が奇数なら -1、偶数なら +1
    • C コードでは __builtin_popcount(x & y) % 2 として実装する

Walsh-Hadamard 変換の実装

  • Hadamard 行列を直感的な Walsh 順に変えるには、行を sequency 基準でソートする必要がある
  • 最も単純な方法は、各行で符号変化の回数を数えることである
  • ほかのビット操作方式も可能である
    • Walsh の行番号を、それ自身を 1bit 右にシフトした値と XOR して Gray code を作る
    • 最後の nbit の順序を反転し、Hadamard の行マッピングを計算する
  • こうして作った Walsh 配列で DCT 実装の基底を置き換えれば、“discrete square transform” と逆変換を作れる
  • 技術的には、この変換は Walsh–Hadamard transform(WHT) である
  • サンプル入力 1 1 1 1 5 5 5 5 を DCT で処理すると、複数の周波数 bin に高調波成分が広がる
    • DCT : +24.00 -10.25 -0.00 +3.60 +0.00 -2.41 -0.00 +2.04
  • 同じ入力を方形波変換で処理すると、F_0F_1 にだけ非ゼロ成分が現れる
    • SQFT : +24.00 -16.00 +0.00 +0.00 +0.00 +0.00 +0.00 +0.00
  • 逆変換 isqft() は元の入力を復元する
    • ISQFT : +1.00 +1.00 +1.00 +1.00 +5.00 +5.00 +5.00 +5.00

スペクトログラム比較と実用上の位置づけ

  • Gorillaz の “DARE” から取得した 11 秒のオーディオクリップ で、DCT スペクトログラムと Walsh-Hadamard スペクトログラムを比較する
  • Walsh-Hadamard 変換は低性能なコンピュータでも計算効率がよく、特定種類のデータに適しているため、いくつかのニッチな用途で使われている
  • 結論は WHT をもっと使うべきだということではなく、離散 Fourier 変換が真実を独占しているわけではない という点である
  • スペクトログラムは 44.1kHz mono オーディオファイルから DCT と WHT で計算されている
    • 入力サンプル window は 512
    • transform stepover は 1
    • 出力配列サイズは約 512 × 485k
    • ピクセル強度は、正規化した絶対値に gamma 約 0.4 を適用している
    • 画像は Lanczos resampling でリサイズし、黒–水色–白の線形 colormap でレンダリングしている
  • Walsh-Hadamard を画像圧縮に試した事例として http://rotormind.com/blog/2019/hadamard-days-night/ も併せて紹介されている

1件のコメント

 
GN⁺ 2024-04-08
Hacker News の意見
  • 数学的には、フーリエ変換は時間信号を特定の直交ベクトル基底で表現する方法にすぎない
    地表面の変位ベクトルも北/東方向の基底で表現できるし、ある道路の方向とそれに垂直な方向でも表現できる
    時間依存の信号や「きれいな」関数は無限次元ベクトル空間にあるため想像しにくいが、核となる数学は似たように働く
    フーリエ変換では基底ベクトルが調和関数であり、周波数領域は、無限に多くの調和関数の組み合わせとして信号を示す一つの「地図」である
    Walsh–Hadamard 変換のような別の基底による地図も同じように実在し、時間領域での表現も、私たちに馴染みがあるだけで複数ある地図の一つにすぎない

    • 正しい回答で、さらに付け加えるなら、フーリエ変換は時間信号だけでなく、区分的に連続・微分可能で Dirichlet 積分可能な関数にも適用される
      画像処理、微分方程式の解法、高速な乗算など応用は多い
      数学的にはこうした変換は無損失なので、変換後の関数は元の関数と正確に同じ情報を持ち、変換後のものだけがあっても元に戻せる
      工学的には、特定の周波数成分のような不要な情報を捨てるために変換することが多く、この点がしばしば曖昧になる
      結局は、一つの関数を見る複数の視点の一つである
    • 以前は私も「別の基底」だと考えていたが、最近ではその比喩は少し危うい、少なくともすべてではないと思っている
      特に多次元空間では、一般的な多次元フーリエ変換は、その空間に平坦な計量がある場合にだけ適切に機能する
      宇宙そのものが曲がっていることを考えると、警告サインのように見える
      最近、特定の双曲格子でフーリエ級数を一般化した興味深い研究があり、その結果、位置空間よりフーリエ空間の次元のほうが高くなり得る
      しかもこの「フーリエ空間」の次元は格子の離散化方法によって変わるため、ある2次元格子は4次元の周波数類似領域を持ち、別の2次元格子は8次元の類似領域を持つことがある
      https://arxiv.org/abs/2108.09314 または https://www.pnas.org/doi/full/10.1073/pnas.2116869119
    • 昔の天文学者は周転円を付け加えた地球中心の宇宙モデルを信じており、さらに精度が必要になると周転円をさらに追加した
      完全に間違ったモデルだったが、実質的にはフーリエ級数を関数近似器として使っていたことになる
    • 全体としては同意するが、すべての正規直交基底は周波数スペクトルを分割する
      多項式のような基底を使っても、結局は周波数成分で関数を組み立てていることになる
      フーリエ基底は、各要素が特定の周波数に対応するという点で特別だ
      ただし各基底は目的に合わせて設計されたものに近く、基底変換はスペクトルを、分析しにくい形で並べ替えることがある
      その場合は滑らかさのような別の性質を分析することになる
      関心のあるほとんどの関数には特徴的なスペクトルがあるが、フーリエ基底がすべての問いに答えてくれるわけではない
    • 基底ベクトルは必ずしも垂直である必要はないのでは?
      北と北東のように、ある程度垂直成分があれば [n, e] も別の座標で表現できる
      具体的な係数は酔っているせいで間違っているかもしれないが、要点は可能だということだ
  • 修士のとき、力学系グループでホワイトボードの前で交わした会話を思い出す
    「左側からシステムにエネルギーが注入され、右側のここで散逸します」
    「でもシステムは回転不変だから、左も右もないでしょう」
    「周波数空間の話です」
    「ああ、実空間のことだと思っていました」
    「ばかですか? 誰が実空間で考えるんですか?」

    • ちょっと待って、周波数空間にも左と右があるのか?
      抽象的な表現なのだから、空間次元の左右上下とは直接関係ないのでは
    • これまでの人生で象牙の塔の学界に向けて吐き出してきた軽蔑だけを集めても、小さな島国一つを10年は動かせそうだ
    • 修士課程では互いをばか呼ばわりするものなのか、それともかなり脚色しているのか?
  • 複素指数基底関数が線形時不変(LTI)システムの固有ベクトルであるという点で、フーリエ基底は独特だ
    他の変換にはこの性質がない
    回路、通信チャネル、アンテナなど多くの現実のシステムは LTI であり、この性質のおかげで異なる周波数で送られた信号は干渉しない
    だからフーリエ変換は他の変換より広く使われている
    量子物理でも位置と運動量の波動関数としてフーリエ対を使うつながりがあり、他の変換にはこうした性質はない

    • この話を持ち出す人がほとんどいないことに驚いた
      電気工学の背景では、解析のために多くのシステムを線形または非常に弱い非線形だと仮定し、信号もおおむね周期的なのでフーリエ変換が自然になる
      畳み込みは乗算になり、複素指数の時間微分は j*omega を掛けることになる
      畳み込みや時間微分を扱うより、乗算をするほうがずっとよい
      「よくある特定の状況で便利だからフーリエ表現を使う」と受け入れれば、別の問題に別の数学的変換を使うことも驚くことではない
    • 技術的にはラプラス基底の特殊な場合ではないのか?
      多くの講義が、最も一般的な両側ラプラス変換をきちんと扱わず、両側フーリエ変換から片側ラプラス変換へすぐ移る点をいつも不思議に思っていた
      https://en.wikipedia.org/wiki/Two-sided_Laplace_transform
  • 「実在する場所」なのかと聞かれると、昔の光学実験を思い出す
    絵をいくつかのレンズに通すと周波数の平面ができ、さらにレンズを通ってスクリーンに投影される
    その周波数平面の一部を遮ると画像が変わる
    扱うのがものすごく難しく、St Andrew’s の Dr Bruce Sinclair にはとても感謝している
    物理実験室での作業は物事がどう動くのかを見せてくれるが、実験から数か月後に理論を見直すとかなり迷子になる

    • 周波数の平面は常に生じるものではないのか?
      だから開口が解像度を制限し、反射望遠鏡で回折スパイクが生じる、といった話につながるのだと思う
    • ソフトウェアでもできる: https://imagemagick.org/Usage/fourier/#noise_removal
    • 光学で周波数領域表現を得られる性質は、光源ベースのさまざまな顕微鏡法や分光法でかなり便利
    • フーリエ光学である
      シュリーレン法もこのように動作する
  • DFT のもう一つの興味深い一般化として Lomb-Scargle 変換がある
    時間領域で固定された測定間隔を必要としない
    天体物理のように測定間隔が一定でない場合に、周期信号の周波数を見つけるためによく使われる
    https://iopscience.iop.org/article/10.3847/1538-4365/aab766 は一般的な紹介で、https://docs.astropy.org/en/stable/timeseries/lombscargle.ht... は Python の astropy ライブラリでの使い方をよく説明している

    • 最近、Prometheus のデータを周波数領域で表現すると、容量計画における週次・日次・年次のアクセスパターンを可視化するのに役立つのではないかとよく考えている
      オートスケーリングは障害につながる性能低下を避ける助けにはなるが、年間予算がいくらであるべきか、その理由までは教えてくれない
      ただし Prometheus のデータは、本当のサンプリング間隔とは言いにくい
      クラスター内の各マシンが一定間隔で報告していても、互いに同期しているわけではない
  • 別の見方をすると、蝸牛はフーリエ変換の「実際の」実装のように見なせる
    https://www.britannica.com/science/sound-physics/The-ear-as-...

    • 蝸牛はむしろ記事の要旨を裏付けている
      周波数領域には変換するが、フーリエ変換をしているわけでも近似しているわけでもない
      蝸牛が「実装」している時間→周波数領域変換は、ウェーブレット変換により近い
      蝸牛をフーリエ変換として解釈するのは、目の錐体細胞が赤・緑・青の光にだけ反応すると考える誤りに似ている
      実際には各細胞は一定範囲の周波数にわたって異なる反応を示す
      錐体細胞は低・中・高周波数領域でピークを持ち、その両側で減衰する。蝸牛の有毛細胞はピーク周波数の倍音に二次ピークがある、よりウェーブレット的な応答曲線を持つ
      専門家ではなく熱心なアマチュアなので、もっと詳しい人が訂正してくれることを期待している
    • 生物学が物理学にどれほど深く根ざしているかは忘れがちだ
      大学時代、私たちの幹細胞株が骨へ分化してしまう問題があったが、実は環境の硬さが幹細胞に感知できるシグナルだった
      硬い培養皿が細胞に、骨細胞になるべきだと伝えていたようなものだった
  • 記事では「Hadamard 行列の行を sequency の基準で並べるには、ゼロ交差の回数を数えるよりエレガントなアルゴリズムを知らない」と書いていたが、行列を見てパターンを推測したところ、すでに知られた方法だった
    https://en.wikipedia.org/wiki/Walsh_matrix によると、Walsh 行列の sequency 順序は、Hadamard 行列にまずビット反転順列を適用し、続いて Gray-code 順列を適用することで得られる

  • 記事は非常に一般的で哲学的な問いを投げかけているが、その後は他の直交基底や変換も見つけられるのだから、周波数領域はそれほど特別ではないと言っている
    それでも周波数領域とフーリエ変換は、他の多くの変換より特別だと思う
    自然の中で直接観察できるからだ
    例えばレンズは、平行光に載った入力画像の2次元フーリエ変換を行い、それをスクリーンで見ることができる
    また、格子やプリズムの出力を CCD に投影して光の波長や周波数を測定できるが、これも周波数領域の直接測定である
    RF 波でも似た測定が可能だ

  • 正弦波は Helmholtz 波動方程式の自然な解であるという点で特別だ
    矩形波には無限エネルギーのような別の問題もある
    この記事は数学者やコンピューター科学者には筋が通るかもしれないが、音と波の根本的な物理を見落としている

    • 正弦波は微分演算子の固有関数であるという点でも特別だ
      物理的な結果は、おそらくその性質の結果である可能性が高い
      結局、現代数学の中心的な教訓は、対象を複数の観点から見ることが有用だということにある
    • その通りだ
      非常に多くの物理的物体が調和振動子であり、これは物理学にかなり根本的な基盤を持っている
      フーリエ解析を使える別の場面もたくさん思い浮かぶが、正弦波は物理的にはより「実在」しており、どんな基底集合でも表現できるという話はより「妥当」な側に近い
      「実在」という言葉は、現象の背後に実際の振動子があるという感覚を与えるようだ
      矩形波は信号と導関数の両方に不連続があるため物理的ではなく、自然は不連続を本当に好まない
    • 周波数領域は、自然に存在する多くのシステムを近似的に説明する線形時不変演算において、数学を非常に簡単にしてくれる
      例えば Gibbs 現象は、あるカットオフ周波数より上のすべての周波数を 0 にした周波数応答の逆フーリエ変換から自然に現れる
      矩形波の周波数領域が Gibbs 現象をどのように説明するのかは気になる
      おそらくシステムが非線形であるかのように、基本矩形波周波数の高調波が現れるのだと思う
  • 学部で物理・数学を学ぶ中で、関数 f(x) の値を無限に多くの x で知っていることと、f の周波数成分を無限に多くの周波数で知っていることは等価だ、という結論に至った
    哲学的には、2つの表現は同じように「実在」している
    ある問題が、一方の表現ではもう一方より解きやすいだけ

    • 完全に同意する
      時間領域から周波数領域へ変えるのは、座標系の変換のようなもの
      時間領域に狭いピークが1つある信号は、ピーク位置のデルタ1つで非常に小さく疎に表現できるが、周波数領域ではそのように圧縮された表現にはならない
      逆に、時間領域の正弦波信号はそちらでは compact ではないが、周波数領域ではデルタが数個あればよい
      時間と周波数は同じものを表す2つの方法であり、場合によっては一方の領域がより簡単で、別の場合にはその逆になる
      時間領域で有界なものは周波数領域では非有界になり、その逆も成り立つことを証明できる
      したがって、一方の領域で compact なものは、もう一方の領域へ移ると必ず広がる
      量子力学では、位置と運動量が上の時間・周波数と同じように共役変数なので、位置が有界なら運動量は非有界になり、その逆も成り立つ
      これが Heisenberg の不確定性原理の核となる考え方である