Thread by @genkuroki

#Julia言語 『数値計算の常識』という有名な本があって第5章のタイトルが「逆行列よさようなら」です。

Juliaでは計画行列Xによるバックスラッシュ演算

β̂ = X \ y

の一発で最小二乗法も計算できます。

github.com/genkuroki/publ…

https://twitter.com/usazo/status/1437372320700456968

#数楽 Xが縦長の行列で、β, y が縦ベクトルのときの、βの成分に関する連立一次方程式

Xβ = y

は一般に解を持たないのですが、Xβとyのユークリッド距離を最小にするようなβをβ̂と書いて、「解」とみなすのが最小二乗法の考え方です。その「解」を

β̂ = X \ y

と書くことは記号法的に自然です。

#Julia言語 Juliaがバックスラッシュ二項演算子で最小二乗法も可能にしていることの背景には以上のような数学が隠れています。

#数楽 はい、それも正しいです。

しかし、「ムーアペンローズ一般化逆行列」の話は結構複雑で、「部分空間{Xβ|β∈ℝʳ}⊂ℝⁿの点Xβと点yのユークリッド距離が最小になるβ=β̂を求める」の方がより基本的な考え方だと思います。続く

https://twitter.com/tchaikovsky1026/status/1437591438431195140

#数楽 縦長のn×r行列Xがランクがrになる場合にはβ̂が一意的に決まる。

そうでない場合には、β̂が一意的ではなくなるので、数値計算では複数の仕様が考えられる。ノルムが最小のβ̂を選んだり、sparseで計算しているならβ̂は"sparsest"になるようにしたい。後者でβ̂はMP一般化逆行列解ではなくなる。

実際、#Julia言語 のdiscourseでは、sparseの場合の最小二乗法(もしくは直交射影)に関係した質問が定期的に出て来ている感じ。

最近では↓

discourse.julialang.org/t/differences-…

sparse行列AについてA\bはpinv(Matrix(A))*bと一致しない。

github.com/genkuroki/publ…

#数楽 「Moore-Penroseの一般化逆行列」は直観的もしくは視覚的には易しく感じられる「線形写像の像への直交射影」の「実装の詳細」にあたる複雑な話題だと思います。

「実装の詳細」以外にも、「線形写像の像への直交射影」について図を描いて直観的にも理解しておいた方が私は得だと思います。

分岐スレッドへ↓

https://twitter.com/genkuroki/status/1437608571378696193

#数楽 最小二乗法が直交射影の特別な場合であることについては、解説が私のツイログにたくさんあります。

twilog.org/genkuroki/sear…

部分空間のW=Xℝʳの点で点yに最も近いものは、点yの部分空間Wへの直交射影になっています。この手のことは実質的に高校数学でも習っているはずのことだと思います。

#数楽 Xが縦長のn×r行列で、yが長さnの縦ベクトルのとき、長さrの縦ベクトルβに関する一次方程式

Xβ = y

は一般に解を持ちません。しかし、残差eを許して

Xβ + e = y

の形なら解を作れる。そのとき、残差eをユークリッドノルムを最小化するのが最小二乗法です。

#統計 1つ前のツイートのような説明をすると、なぜかより複雑で難しい「Moore-Penroseの一般化逆行列」もしくは「rank X=rの場合のその公式 (X'X)⁻¹X' 」(X'はXの転置)で最小二乗法について勉強してしまった人は、「知っている話と違う！」となってしまう。これはちょっと不幸な感じがします。

#数楽 Moore-Penroseの一般化逆行列の真の理解のためには、直交射影の線形代数をやる必要があり、結果的に私が上で述べたシンプルな理解に戻ることになります。

そして、一般化逆行列に頼らずに直交射影によって理解していれば、解β̂が一意に決まらない場合のsparseでの計算にも対応できます。

もっとはっきり言えば、sparse行列の場合には「#Julia言語 におけるA\bはMoore-Penroseの一般化逆行列の意味でAx=bを解いた結果になる」という主張は誤りです。

それに対して「#Julia言語 において、Axとbの距離が最小になるxはA\bで求められる」ならば、sparse行列でも正しいです。

https://twitter.com/genkuroki/status/1437608587157639169