トップ 差分 一覧 Farm ソース 検索 ヘルプ PDF RSS ログイン

nls

*disclaimer
178094

[R]
R

非線形回帰分析



 nls()


nls(formula,start,trace)

nls(式, パラメタの初期値)

残差のプロット

plot(resid(モデル))

References

http://cse.naro.affrc.go.jp/takezawa/r-tips/r/73.html

 例:Syntactic Tree Fragmentsがzipfの法則にしたがうか

fragments.dat <- read.table("fragments.dat.txt")
head(fragments.dat)
##           fragment rank freq
## 1 (PP (IN ) (NP ))    1 6257
## 2      (NP (PRP ))    2 4242
## 3  (S (NP ) (VP ))    3 4079
## 4        (S (VP ))    4 3714
## 5         (DT the)    5 3450
## 6 (NP (DT ) (NN ))    6 3260

frag.short <- head(fragments.dat, 500)
dim(frag.short)
## [1] 500   3
str(frag.short)
## 'data.frame':    500 obs. of  3 variables:
##  $ fragment: chr  "(PP (IN ) (NP ))" "(NP (PRP ))" "(S (NP ) (VP ))" "(S (VP ))" ...
##  $ rank    : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ freq    : int  6257 4242 4079 3714 3450 3260 3101 3085 2850 2569 ...

plot(frag.short$rank, frag.short$freq, xlim=c(0,500), ylim=c(0,7000))


Zipfの法則のモデル式


fit <- nls(freq ~ a/(rank^b), start=list(a=6257,b=1), data=frag.short, trace=T)

1.110042e+08 (4.97e+00): par = (6257 1)
4.156852e+07 (1.59e+00): par = (6095.586 0.6873339)
1.741337e+07 (5.77e-01): par = (6665.465 0.6063987)
1.308021e+07 (2.59e-02): par = (7571.729 0.5972984)
1.307277e+07 (4.32e-03): par = (7599.867 0.6002608)
1.307257e+07 (8.20e-04): par = (7590.96 0.5997113)
1.307256e+07 (1.55e-04): par = (7592.631 0.5998156)
1.307256e+07 (2.94e-05): par = (7592.314 0.5997958)
1.307256e+07 (5.58e-06): par = (7592.374 0.5997996)

Formula: freq ~ a/(rank^b)

Parameters:
   Estimate Std. Error t value Pr(>|t|)    
a 7.592e+03  1.210e+02   62.74   <2e-16 ***
b 5.998e-01  5.694e-03  105.34   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 162 on 498 degrees of freedom

Number of iterations to convergence: 8 
Achieved convergence tolerance: 5.576e-06

plot(frag.short$rank, frag.short$freq, xlim=c(0,500), ylim=c(0,7000))
lines(frag.short$rank, fitted(fit), col="red", lwd=3)


得られたパラメタを使ったモデル式



どのくらいフィットしているか、残差分析プロット

plot(resid(fit))