5

非線形モデルを当てはめようとしていますが、オンラインで良い例を見つけることができません。

この関数に名前はありますか?

線形化できますか?

aパラメータ、b、およびを時間の関数として (グループのように)cランダムな効果で推定しようとしました。ランダム効果を使用せずにモデルを適合させることはできますが、モデルを収束させるのに問題があります。提案を歓迎します (できれば R 内ですが、適切なパッケージであれば何でも構いません)。gtnls

## time, repeated 16 times for 4 replicates from each of 4 groups  
t <- rep(1:20, 16)
## g, group
g <- rep(1:4, each = 80)

## starting to create an example dataset, 
## to see if I can recover known parameters
a <- rep(c(3.5, 4, 4.1, 5), each = 80)
b <- rep(c(1.1, 1.4, 1.8, 2.5), each = 80)
c <- rep(c(0.125, 0.25), each = 160)

## error to add to above parameters
set.seed(1)
e_a <- runif(320, -0.5, 0.5)
e_b <- runif(320, -0.1, -0.1)
e_c <- runif(320, -0.02, 0.02)

## this is my function

f <- function(t, a, b, c) a * (t^b) * exp(-c * t)

## simulate y
y <- f(t = t, a + e_a, b + e_b, c + e_c)

mydata <- data.frame(t = t, y = y, g = g)

library(nlme)
## now fit the model to estimate a, b, c
fm1 <- nlme(y ~ a * (t^b) * exp(-c * t),
            data = mydata,
            fixed = a + b + c~1,
            random = a + b + c ~ 1|g,
            start = c(a = 4, b = 1, c = 0.25),
            method = "REML")
4

1 に答える 1

9

物理学 (および他のいくつかの分野) では、これまたはその変形がHoerl曲線または Hoerl 関数と呼ばれているのを見てきました。cが負で、abが正の場合、スケーリングされたガンマ密度です。

線形化について尋ねるときは、注意が必要です。式y = at ^ b . exp( ct ) は実際にはあなたが意味するものではありません - 観測値 y ( i ) はaと正確に等しくありません。t (i)^ b . exp( ct ( i )) (そうしないと、ほぼすべての 3 つの観測値から正確なパラメーター値が得られます)。

したがって、ノイズは何らかの形でyのモデルに入らなければなりません。添加物ですか?乗法、または何か他のもの?(これも重要ですが、他の理由から: tが変化すると、そのサイズは何らかの形で変化しますか? 異なる観測のノイズ項は独立していますか?)

実際のモデルが y ( i ) = at ( i )^ bの場合。exp( ct ( i ))+ε( i )、これは線形化できません。

実際のモデルがy ( i ) = at ( i )^ bの場合。exp( ct ( i )) . ε(i)、および ε( i )=exp(η( i ))は、線形化可能な (できればゼロ平均の) η( i ) に対してです。

第二形態をとると、

log( y ( i )) = log(a) + b log( t ( i )) + c t ( i ) + log(ε( i ))

また

y *( i ) = a* + b.log( t ( i )) + c. t ( i ) + η( i )

これは、パラメーターa * = log( a )、bおよびc、および誤差項 η( i ) で線形です。したがって、エラーについてそのような仮定を行う準備ができている場合は、そのような線形モデルに適した方法でそれを適合させることができるはずです。その場合、モデル化の方法に影響を与える可能性がある上記の誤差項に関する括弧内の質問について熟考することをお勧めします。

于 2015-10-15T23:13:50.527 に答える