{{category R}}
R
!!!非線形回帰分析
{{outline}}
----

!!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の法則にしたがうか
{{pre
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))
}}
{{ref_image rank.freq.png}}


!Zipfの法則のモデル式
{{ref_image zipf.png}}

{{pre
fit <- nls(freq ~ a*(rank^b), start=list(a=6257,b=1), data=frag.short, trace=T)
## 1.635864e+15 (2.79e+03): par = (6257 1)
## 9.472708e+09 (6.64e+00): par = (15.44635 0.9995952)
## 1.408801e+09 (2.46e+00): par = (15.48336 0.8352084)
## 2.321331e+08 (7.96e-01): par = (40.77552 0.4017383)
## 2.161933e+08 (3.00e+00): par = (151.2654 -0.09240324)
## 2.119514e+08 (5.33e+00): par = (414.3436 -0.3945464)
## 1.899423e+08 (3.94e+00): par = (1295.112 -0.6965551)
## 8.183550e+07 (2.58e+00): par = (2768.631 -0.4813824)
## 4.029910e+07 (1.45e+00): par = (5244.751 -0.6189928)
## 1.347308e+07 (1.85e-01): par = (7529.115 -0.5842222)
## 1.307798e+07 (2.22e-02): par = (7636.171 -0.6021353)
## 1.307275e+07 (4.16e-03): par = (7585.197 -0.599355)
## 1.307257e+07 (7.86e-04): par = (7593.71 -0.5998829)
## 1.307256e+07 (1.49e-04): par = (7592.11 -0.5997831)
## 1.307256e+07 (2.82e-05): par = (7592.413 -0.599802)
## 1.307256e+07 (5.34e-06): par = (7592.356 -0.5997984)
summary(fit)
## 
## 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: 15 
## Achieved convergence tolerance: 5.342e-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)
}}
!得られたパラメタを使ったモデル式
{{ref_image zipf.tf.png}}

{{ref_image zipf.fit.png}}

*どのくらいフィットしているか、残差分析プロット
 plot(resid(fit))
{{ref_image zipf.resid.png}}