Rで近似曲線（直線ではない曲線）
Rで回帰分析のため、回帰直線ではなく、近似曲線を求めてみます。
今回使うのは以下のようなデータ。
No,データ 1,121696 2,122012 3,122469 4,123033 5,123756 6,124454 7,125019 （略） 46,125983 47,126249 48,126643 49,127263 50,127816
このデータをモデル化してplotしておきます。
x <- data[,1] y <- data[,2] model = (y~x)
明らかに５次式で表されるグラフのようですが、簡易的に三次式で近似します。
近似する場合には、まず未定係数の三次式を作成します。
y = a*x^3 + b*x^2 + c*x + d
この数式の係数を求めるために、nls関数を用います。
nlsを用いるときに、a,b,c,dのおおよその初期値を与えておくようです。（この初期値が適切でないと実行できないようです）
> result <- nls(y ~ (a * x^3 + b * x^2 + c * x + d), start = c(a=0, b=0, c=0,d=121696))
> result
Nonlinear regression model
model: y ~ (a * x^3 + b * x^2 + c * x + d)
data: parent.frame()
a b c d
6.766e-01 -5.159e+01 1.058e+03 1.203e+05
residual sum-of-squares: 35571879
Number of iterations to convergence: 1
Achieved convergence tolerance: 6.783e-08
> summary(result)
Formula: y ~ (a * x^3 + b * x^2 + c * x + d)
Parameters:
Estimate Std. Error t value Pr(>|t|)
a 6.766e-01 5.279e-02 12.82 < 2e-16 ***
b -5.159e+01 4.094e+00 -12.60 < 2e-16 ***
c 1.058e+03 9.032e+01 11.71 2.11e-15 ***
d 1.203e+05 5.372e+02 224.00 < 2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 879.4 on 46 degrees of freedom
Number of iterations to convergence: 1
Achieved convergence tolerance: 6.783e-08 resultの結果から
a b c d 6.766e-01 -5.159e+01 1.058e+03 1.203e+05
y = 0.6766x^3 + 51.586x^2 + 1058x + 12030 と求められます。
また、このグラフをplotするには
predict.c <- predict(result) plot(x, predict.c, type="l")
とすればよいです。
もとのグラフと合わせると
plot(model) par(new=T) plot(x, predict.c, type="l")
上記のようにおおよそ近似した曲線が描かれます。
なお、同様に五次式で近似しようとして、下記のように試みましたがエラーになるようです。この初期値について適切な値というのはどうやって求めるものなのだろうか・・・
> result <- nls(y ~ (a * x^5 + b * x^4 + c * x^3 + d * x^2 * e * x + f), start=c(a=0.1, b=0.2, c=0.3, d=40, e=50, f=10000)) 以下にエラー nlsModel(formula, mf, start, wts) : パラメータの初期値で勾配行列が特異です
このやりかたを「曲線あてはめ（カーブフィッティング）」というらしく、数式が三次曲線であったりとか、成長曲線であるロジスティック曲線（指数関数）に対しての未定係数（ロジスティック関数の場合はa,b,cの三つ）を求めるのに利用することができます。Excelを使って行う場合には、ソルバーを利用して求めるものになりますね。
