|
[推定されたセン]
上の R コード例の sample (ねつ造曲線+正規乱数) を用いた推定． q の真の値は 0.7 …… あ， 推定値は 0.75 だから値としては そんなにズレてないぞ． f と C はほぼ正しく推定できてる． |
mlfitting
のことを書いたら，
これのユーザーである牧さんから質問メイルいただく
……
悪夢のごとく非線型なる非直角双曲線，
すなわち光合成曲線の推定に関して，
である．
mlfitting
と
R
の両方で推定を工夫してみる．
推定すべきパラメーターに与える初期値に
スルどく反応してしまうあたりが，
なんとも面倒なところで．
R
では下のように nls という
関数を用いる．
xy <- scan("sample.txt", list(x = 0, y = 0))
require(nls)
Pmax = 10.2 # きめうち
psynthesis <- nls(
y ~ (f * x + Pmax - sqrt((f * x + Pmax)^2 - 4 * (f * x * q * Pmax)))
/ (2 * q - C), # この式は間違いと翌日判明
data = xy,
start = list(f = 0.1, q = 0.5, C = 0.0)
)
summary(psynthesis)
(結果の出力例)
Parameters:
Estimate Std. Error t value Pr(>|t|)
f 0.04712 0.00218 21.612 < 2e-16 ***
q 0.75033 0.05345 14.038 < 2e-16 ***
C 0.10854 0.01582 6.862 6.63e-10 ***
q
の推定が悪いためで，
R
では下の図のように過大推定となり
mlfitting
では過小推定ぎみになる．
|
|
[推定されたセン]
上の R コード例の sample (ねつ造曲線+正規乱数) を用いた推定． q の真の値は 0.7 …… あ， 推定値は 0.75 だから値としては そんなにズレてないぞ． f と C はほぼ正しく推定できてる． |
R の nls
では Gauss-Newton 法でパラメーターを推定値に
収束させている
……
いっぽうで牧さんにご教示いただいたところによると，
先方でよく用いられてる「でるたぐらふ」
では Levenberg-Marquardt 法を使ってるそーだ．
このアルゴリズムは知らなかったので
手元の数値計算本をながめてみると，
なかなかホメられている．
R と言いますか
S 言語 いやいや統計学全般に関する
一言居士たる英国人 D. Ripley 教授が
「Levenberg-Marquardt 法なんてたいしたこたねーよ」
と
いちゃもんつけてる
と，
これはかとー先生に教えてもらった．
mlfitting
はヘボいわけだが
……
今回の体験から，
自前の汎用最尤推定プログラムがあれば
他のソフトウェアが出してきた結果を突き合わせたり
できるとわかった
……
ということで，
使い勝手から推定精度にいたるまで
あれこれ問題あるにせよ，
これを運転可能な状態に維持しておくことは重要だな．
何人かの
mlfitting
ユーザーのおかげで存在価値がわかった，
と言いますか．
R の glm による計算)
がかなりうまくいってるようだ，
という連絡をいただく
……
ひょっとして実際の観測データにたいして
この推定法を適用した
初めての日本人生態学者になるんではなかろーか．