対数尤度を「浮動小数点数で尤度を計算してからlogをとる」で計算するとマイナス無限大になりがち。

尤度を浮動小数点数で計算するとぴったり0.0になりがち。

他にもマイナス無限大が得られるパターンが結構あると思う。
確率がらみの数値計算では、対数尤度は「尤度の対数」で計算しちゃだめで、対数尤度を直接計算するコードを書く必要がある。

例えばgamma(x)やbeta(a, b)の類は使っちゃダメで、loggamma(x)やlogbeta(a, b)を使う。

混合正規分布を使っている場合にはlogsumexp的な工夫が必要になる。
確率密度の値は浮動小数点数では容易に0.0やInfになるので、確率密度を扱うパッケージには、確率密度の値の対数を直接扱う機能が必ず付いているはず。

私はほんの一部の狭い分野を除くあらゆる分野のど素人なので、この手のことも知らずに、結構酷いコードを書いていたことがある。
私は数値計算がらみで疑問が生じたら #Julia言語 のdiscourseで関連のキーワードを検索しています。これ、Juliaを知っているとかなり便利で、GitHubにあるソースコードを確認して理解を深められる。



discourse.julialang.org/search?q=under…
underflowを検索
log(1+x)を計算してくれるlog1p(x)函数が基本初等函数ライブラリの一員であることもJuliaで初めて知った。こんなことも知らなかった。

初等函数は結構難しい。今でも全然自信がない。

• • •

Missing some Tweet in this thread? You can try to force a refresh
 

Keep Current with 黒木玄 Gen Kuroki

黒木玄 Gen Kuroki Profile picture

Stay in touch and get notified when new unrolls are available from this author!

Read all threads

This Thread may be Removed Anytime!

PDF

Twitter may remove this content at anytime! Save it as PDF for later use!

Try unrolling a thread yourself!

how to unroll video
  1. Follow @ThreadReaderApp to mention us!

  2. From a Twitter thread mention us with a keyword "unroll"
@threadreaderapp unroll

Practice here first or read more on our help page!

More from @genkuroki

13 Dec
#統計 添付画像下段はよくあるベイズ版95%信用区間(確信区間)の最高事後密度(HPD)版です。95%の部分は0から1の間の任意の実数に一般化できる。

添付画像上段のベイズ版P値は「信頼区間と検定が表裏一体」という原理を適用して定義されています。

nbviewer.org/github/genkuro…
#統計 ベイズ版95%信用区間もベイズ版のP値も事後分布の情報だけを使って定義されており、通常の信頼区間やP値の作り方とは異なります。

しかし、数学には定義が全然違う量がある条件のもとでほぼ同じ値になることを証明できる場合が多数あって、解析学の基本的な考え方になっています。
#統計 実践的には無視できる違いしかない2つの数学的量を使った分析や推論は平等な扱いをする必要があります。

❌定義の違いは規範の違いから来ている。違いを無視できるほど同じ値になるとしても、それらは異なる使い方をされなければいけない。

などと言うと、単なるトンデモさんになってしまう。
Read 9 tweets
13 Dec
#統計 いやあ、これにはマジでびっくりした。😱

ベイズ版95%信用区間(確信区間)がゼロをまたがなくなるまでNを増やす行為は、P値が5%を切るまでNを増やすp-hacking行為と(小さな誤差を除けば)数学的に同等です。

P値が5%を切るまでNを増やすp-hacking行為の是非も規範の問題で済ませられる?
#統計 2020年3月に出版された豊田秀樹さんの本には、ベイズ統計では【データを取り増しても】よいし【事前登録~も必要ありません】と言っています。(他にも色々!)

これを「ベイズ統計を使えば研究不正扱いされていたことを堂々とできるようになる」と解釈する人達が出ることを危惧していました。
#統計 より正確に言えば、同様の言説を他で見た人達が、「ベイズ統計を使えば研究不正扱いされていたことを堂々とできるようになる」と既に信じている場合もあるのではないかと思いました。

これ、めちゃくちゃまずくないですか?
Read 26 tweets
12 Dec
#Julia言語 いきなり、struct ~ end の中で(内部)コンストラクタをnewを使って定義するのは避けた方が無難。new無しにはできないことをやりたくなったら、初めてそうするのがよいです。

公式ドキュメントにも可能な限り内部コンストラクタは避けよと書いてあります。

docs.julialang.org/en/v1/manual/c… Image
#Julia言語 内部コンストラクタを定義せずに、例えば、

struct Foo{A, B}
a::A
b::B
end

と書けば、Foo型オブジェクトのコンストラクタたち

Foo(a::A, b::B) where {A, B}
Foo{A, B}(a, b) where {A, B}

が自動的に定義されます。

ひとまずはコンストラクタとしてこれらを使えば大丈夫。
#Julia言語 クラスがあるプログラミング言語から来た人は、struct ~ end の内側でコンストラクタを定義した方が良いと誤解しがちだと思う。

Juliaでは、クラスベースのOOPを廃して、パラメータ付きの型による多重ディスパッチでエコシステム全体を構築した方が効率的だという話になっている感じ。
Read 8 tweets
12 Dec
に訂正があるので注意!

「ある閾値以上大きくなる」の「ある閾値以上」が文字数のプレッシャに負けて抜け落ちてしまった。

具体的な計算例の方を見れば誤解は起こらないはず。
#統計 誤解しないために見て欲しい計算はこれ↓

nbviewer.org/github/genkuro… Image
Read 5 tweets
12 Dec
#数楽 統計学の勉強と無関係に、弧度法による角度の定義とtanの定義まで戻れば、直線 y = ax とx軸で挟まれた部分の角度が

θ = arc tan a = ∫₀ᵃ dt/(1 + t²)

になることが直接的に導けることを確認すると、高校で習ったことの再構成ができてうれしいかも。続く

genkuroki.github.io/documents/High… ImageImageImage
#数楽 x²+y²=1の右半分とy=txの交点は

(x(t), y(t)) = (1/√(1+t²), t/√(1+t²))

になります。これの軌跡のt=0からt=aまでの長さが弧度法の意味での角度θの定義なので、

θ = arc tan a = ∫₀ᵃ √(x'(t)²+y'(t)²) dt

です。ここでarc tanは角度θをaに対応させる函数tanの逆函数。続く
#数楽 ここまで、高校の数学の教科書に書いてあることしか使っていないことに注意。図を自分で描けば分かり易い。

上の場合にちょっと計算すると

√(x'(t)²+y'(t)²) = 1/(1+t²).

これが arc tan t の導函数であることは、角度とtanの高校数学の教科書における定義からほぼ直接的に導かれる。
Read 7 tweets
12 Dec
かけ算順序固定強制指導を擁護する人達は「算数が苦手な子は理解していないので出て来た数を機械的にかけて済ます。ゆえにかけ算順序の指導はそのような子のために役に立つ」のように言います。

個人的には、これもまた浅はかなウソをついているように感じています。続く
ググって各種の調査結果で児童がどのように問題を解いたり解くことに失敗しているのを見ると、本当に算数が苦手な子が算数を苦手になる理由は出て来た数を機械的にかけているからでは無さそうなことは明らかだと思います。

算数が本当に苦手な子はそういう類の子ではない。続く
例えば、伊藤宏先生の調査結果(添付画像)を見てください。

かけ算順序が問題文に出て来た順序と同じで出題者的に誤りの子は問題文で示された場面を絵で描くことには成功しています。

本当に算数が苦手な子は、そもそもかけ算を使える問題であることが分かっていない。続く

mnavidata.edu-c.pref.miyagi.jp/manage/wp-cont… Image
Read 5 tweets

Did Thread Reader help you today?

Support us! We are indie developers!


This site is made by just two indie developers on a laptop doing marketing, support and development! Read more about the story.

Become a Premium Member ($3/month or $30/year) and get exclusive features!

Become Premium

Too expensive? Make a small donation by buying us coffee ($5) or help with server cost ($10)

Donate via Paypal

Or Donate anonymously using crypto!

Ethereum

0xfe58350B80634f60Fa6Dc149a72b4DFbc17D341E copy

Bitcoin

3ATGMxNzCUFzxpMCHL5sWSt4DVtS8UqXpi copy

Thank you for your support!

Follow Us on Twitter!

:(