Apr 29, 2009

統計学復習メモ12: 不偏推定量とは

以前にt分布を使うパラメーターの区間推定の方法を覚えたが、この前ちょっと推定の実践のネタを思いついてやってみようとしたら、そのパラメーターが正規分布に従うものではないことに気付き、1手目でつまずいた。これでは勉強した意味が無いと思い、推定というものの基礎を勉強し直すことにした。

区間推定をするにもまず、その区間の基準となる値(推定量)を点推定しないと始まらない。その点推定量が満たすのが好ましい性質として、不偏性と一致性がある。推定量θ^がパラメーターθの不偏推定量であるとは、
E(¥hat{¥theta})=¥theta ...(1)
が成り立つことをいう。

出たE(¥hat{¥theta})=¥theta。学生時代にこのθの上の屋根の意味が理解できなかったことが、私が推定・検定の1手目でつまずいた原因の1つであることは間違いない。そもそもなぜ統計学では推定するパラメーターをθと書くのか。高校の三角関数の授業で習って以来、私にとってはθといえば角度に違いないのであり、パラメーターをθと表すのは、円周率をxと表すとか、y=ax+bのxとyが定数でaとbがグラフの縦軸・横軸であるくらいの混乱の元である。

そもそも、θ^がθの推定値だとしたら、(1)が成り立たないとはどういう状態なのだろうか?それに、θが未知だからθを推定するのに、θ^の平均がθであることをどうやって確認するのか?
…という疑問を持ったのは筆者だけなのだろうか。

これは、θが確率変数のパラメーターであり、その推定量θ^を標本から求める方法を決めた時、仮にθによって決まる母集団からの標本を使った場合、θ^の平均がθに一致するなら、そのθ^は不偏推定量である、という意味らしい。

例題として、X1〜Xnが平均μ、分散σ2の母集団からの標本として、まず標本平均X~=E(X)=1/n ΣXiがμの不偏推定量かどうかを考える。
E(¥bar{X})=E¥left(¥frac{1}{n}¥sum_{i=1}^{n}X_i¥right)=¥frac{1}{n}¥sum_{i=1}^{n}E(X_i)=¥frac{1}{n}n¥mu=¥mu
従ってX~はμの不偏推定量である。

次に、標本分散S^2=¥frac{1}{n}¥sum_{i=1}^{n}(X_i-¥bar{X})^2、不偏分散U^2=¥frac{1}{n}¥sum_{i=1}^{n-1}(X_i-¥bar{X})^2がσ2の不偏推定量かどうかを考える。

ここで、Var(X)=E(X^2)-¥left(E(X)¥right)^2の公式とVar(¥bar{X})=¥frac{1}{n}Var(X)の公式を用いると、
E(X_i^2)=Var(X_i)+¥left(E(X_i)¥right)^2=¥sigma^2+¥mu^2
E(¥bar{X})=Var(¥bar{X})+¥left(E(¥bar{X})¥right)^2=¥frac{¥sigma^2}{n}+¥mu^2
が求まるので、
E(S^2)=¥frac{n-1}{n}¥sigma^2
E(U^2)=¥frac{n}{n-1}E(S^2)=¥sigma^2
となる。従って、S2はσ2の不偏推定量ではなく、U2はσ2の不偏推定量である。

See more ...

Posted at 23:34 in 数学 | WriteBacks (0)
WriteBacks

Apr 26, 2009

続・ポリゴンできれいな球を作る

先週は面積が均一なポリゴンで球を作れなかったので、Java3DのプリミティプであるSphereクラスの球を参考にすることにした。

テストアプリ(Java3D)の起動ページへ
(今回のはカーソルキーの左右とPageUp/PageDownキーで回転します)
ソースコード
画面サンプル

Sphereクラスのコンストラクタのdivisionsの値を変えることにより、ポリゴンの密度を変えることができる。右端から反時計回りに、divisionsの値が4,8,12,16,24,32の球である。
それぞれの三角形の面積が大体等しい、きれいな球面である。どういう仕組みなのだろうか。

まず、試しにdivisionsの値を変えながら重ねて表示してみると、5,6,7,8は同じ球、9,10,11,12も同じ球であることがわかった。どうやらdivisionsの値は4の倍数に切り上げられて使われるらしい。
また、divisions=4だと正8面体で、divisionsがいくらでも、その正8面体の辺を球面に投影した(球と同じ半径の)3つの円周で球面が8分割されていることがわかる。さらに、その3つの円周の1つと平行な円周で球面が輪切りにされていることがわかる。

より細かい分割方法の特徴を明確にするために、divisions=28の球面でさらに調べてみる。(ここまでで既に出てきた3や8の倍数は避け、divisions=20は正20面体と関係するかも知れないので避ける)

赤い線で囲まれた部分が、球面の1/8である。赤い線は内接正8面体の頂点を結ぶ円周の弦、水色の線は下の方の赤い線と接する円周(以下、赤道)と平行な円周(以下、緯線)の弦である。
赤道と残り2つの経線上に、それぞれ赤い線分が7つある。おそらく、divisionsの値は赤道や残り2つの経線を何等分するかという値であろう。この1/8球面において、赤道は7分割、その上の緯線は6分割、その上は5分割、以下同様に分割されている。1/8球面を三角形に見立てると、三角形を各辺に平行な線で等分割する要領である。

という訳で、その三角形を等分割する点(線形補間)を球面に投影した点でポリゴンを作ってみる。
テストアプリの起動ページへ
ソースコード
画面サンプル

大きな三角形の頂点の辺りは密で、辺の辺りは疎になった。先週のエントリーで書いたのと同じ問題である。やはり線形補間だとうまくいかない。

そこで、Sphereクラスに習って緯線に沿った円周補間に変える。
…と思ったが、クォータニオンを使った球面補間なら上のソースコードのVector3fをQuat4fに変えるだけだし、球面補間の方が1/8球面を均等に分けられるはずなので、球面補間にする。(Sphereクラスも球面補間か?)
テストアプリの起動ページへ
ソースコード
画面サンプル

うむ、丸い。

See more ...

Posted at 18:58 in Java | WriteBacks (0)
WriteBacks

Apr 20, 2009

ポリゴンできれいな球を作るのは難しい

ポリゴンで球を作ろうとすると、真っ先に考えるのは、地球儀の緯度・経度のように2つの角度で球面上の点を表す極座標モデルであろう。しかし、それだと、赤道の辺りのメッシュが荒く、極点の辺りが細かくなり、球面が均一にならない。何とか、なるべく均等な面積のポリゴン片で球を作れないものか…と考えてみた所、次のような方法を思いついた。

球面上の3点からなる三角形の各辺の中点を取り、それらの間を結んで、三角形を4等分する。そのまま、前述の中点を球面に張り付かせるように動かす。

上から、元の三角形、4等分して球面に張り付かせた4つの三角形、そのポリゴンに色をつけたものである。合同な三角形から成る正八面体から始めて、これを繰り返していけば、均一な球面ができるのではないか?と期待して、やってみた。

テストアプリ(Java3D)の起動ページへ
ソースコード
画面サンプル

右から反時計回りに、8, 32, 128, 512, 2048, 8192個の三角形ポリゴンからなる球であり、手前が線で表示したもの、奥が面で表示したものである。

うまくできているようにも見えるが、細かくしても正八面体の面影が残っている。よく見ると、元々頂点だった所が密で、面の中心だった所が疎であり、くっきりと境目がある。これは、上記の方法で三角形を4等分して各頂点を球面に張り付ける時に、真ん中の三角形の面積が他の3つより大きくなるのが原因である。このことは、三角形の辺の中点が球面から遠い、赤道に近い位置の大きな三角形で考えてみるとわかりやすい。
正八面体に代えて正四面体で始めてみると、その影響が顕著になる。

テストアプリ(Java3D)の起動ページへ
ソースコード
画面サンプル

こうなると、却って緯線・経線の四角形ポリゴンで作る球の方が自然に見えてしまう。

テストアプリ(Java3D)の起動ページへ
ソースコード
画面サンプル

何かいい方法がありそうな気がするのだが…

See more ...

Posted at 23:35 in Java | WriteBacks (0)
WriteBacks

Apr 14, 2009

JScrollPaneを使ったミニゲーム

どうしてもJScrollPaneを使ってスクロールする何かを作ってみたくなった。

Javaアプレットの起動ページへ
ソースコード

See more ...

Posted at 21:55 in Java | WriteBacks (0)
WriteBacks

Apr 05, 2009

mod_proxy+mod_perl+MySQLでProxy Error

2月末にDebianのaptitudeが言うまま全パッケージを更新して、Debian 5.xでサポート対象外となったApache 1.3+mod_perlを手作業でインストールした後、mod_perlで動かしているMovableTypeにmod_proxy経由でアクセスすると、時々

502 Proxy Error
Reason: Document contains no data

というエラーが発生するようになった。

しばらく放置した後に起こることから、mod_perl+MySQLでよく問題になる、時間が経つとMySQLへのコネクションが無効になってしまうために起こる、いわゆる"The morning bug"絡みだろうということはすぐに想像がついて、実際MySQLを再起動すると必ず起こったので、それだと確信した。
しかし、通常はmod_perlがDBD::mysqlを介してMySQLに接続すると、"auto_reconnect"の設定が自動的に有効になるので、この問題は起こらないはずである。
(参考:man DBD::mysqlの"mysql_auto_reconnect"の項)

if either the GATEWAY_INTERFACE or MOD_PERL envionment variable is set, DBD::mysql will turn mysql_auto_reconnect on

実際、DBD/mysql.pmでログ出力してみると、mysql_auto_reconnectの値は1になっていた。しかし、mysql.pmの先は.so形式のバイナリなので、実際に再接続が試されているかは調べられなかった。

ここで、"Proxy Error"は、たまたまこの環境で結果としてそうなるだけで、同じ問題に遭遇した人が皆それに行き当たる訳ではないことは落ち着いて考えると明白だったのだが、"mysql proxy error"で検索して情報を探したため、長期に渡って悩むことになってしまった。
それでも、mod_perl+MovableTypeの組み合わせにApache::DBIを加えてコネクション・プーリングすると何かが解決した、という情報を得て、これだ!と思ってApache::DBIをインストールして組み込んでみたが、解決しなかった。

Apache::DBIのドキュメントを読むと、DBDによってはtimeout時にpingメソッドがエラーを返すので、そういう時は自前のpingメソッドを書くと問題を回避できる、と書いてあった。

Most DBI drivers have a working ping() method, but if the driver you're using doesn't have one and the database handle is no longer valid, you will get an error when accessing the database. As a work-around you can try to add your own ping() method using any database command which is cheap and safe

これだ!と思ってDBD/mysql.pmにpingメソッドを加えてみたが、そもそもpingは呼ばれていなかった。

落ち着いてApacheのエラーログをよく見ると、

[notice] child pid 7987 exit signal Segmentation fault (11)

というのが出ていた。すぐさま、そういう時はPHPをコンパイルし直せという情報を見つけて、これだ!と思って作り直すべきPHPを探すと、どこにもインストールされていなかった。

その後の調べで、DBD::mysqlのv4.007で生じたバグであり、v4.006に戻すかv4.009以降に更新すると直るらしいことがわかった。

仕方なく、CPANを使ってDBD::mysqlをコンパイルしようとすると、mysql_configが無い、というエラーが出た。Debian 5.xではmysql_configはlibmysqlclient15-devというパッケージをインストールしないと使えないらしい。これはわかりづらい。

結局、CPANでDBD::mysqlのv4.010をインストールすることによって、うちの"Proxy Error"は解消した。

See more ...

WriteBacks

Hi, Excellent post

Pretty insightful post. Never thought that it was this simple after all. I had spent a good deal of my time looking for someone to explain this subject clearly and you’re the only one that ever did that. Kudos to you! Keep it up

Posted by Faye Ingram at 2010/7/20 19:25:46