r - 特定のデータへの曲線の当てはめ

Question

私の論文には次のデータがあります。

28 45
91 14
102 11
393 5
4492 1.77

これに曲線を合わせる必要があります。プロットすると、これが得られます。

ある種の指数曲線がこのデータに適合するはずだと思います。私はGNUplotを使用しています。これに適合する曲線の種類と、使用できる初期パラメーターを誰か教えてもらえますか?

Answer 1

Rがオプションの場合に備えて、使用できる 2 つの方法のスケッチを次に示します。

最初の方法: 一連の候補モデルの適合度を評価する

これは、変数間の関係について既に知っていることや期待することを利用するため、おそらく最良の方法です。

# read in the data
dat <- read.table(text= "x y 
28 45
91 14
102 11
393 5
4492 1.77", header = TRUE)

# quick visual inspection
plot(dat); lines(dat)

# a smattering of possible models... just made up on the spot
# with more effort some better candidates should be added
# a smattering of possible models...
models <- list(lm(y ~ x, data = dat), 
               lm(y ~ I(1 / x), data = dat),
               lm(y ~ log(x), data = dat),
               nls(y ~ I(1 / x * a) + b * x, data = dat, start = list(a = 1, b = 1)), 
               nls(y ~ (a + b * log(x)), data = dat, start = setNames(coef(lm(y ~ log(x), data = dat)), c("a", "b"))),
               nls(y ~ I(exp(1) ^ (a + b * x)), data = dat, start = list(a = 0,b = 0)),
               nls(y ~ I(1 / x * a) + b, data = dat, start = list(a = 1,b = 1))
)

# have a quick look at the visual fit of these models
library(ggplot2)
ggplot(dat, aes(x, y)) + geom_point(size = 5) +
    stat_smooth(method = lm, formula = as.formula(models[[1]]), size = 1, se = FALSE, color = "black") + 
    stat_smooth(method = lm, formula = as.formula(models[[2]]), size = 1, se = FALSE, color = "blue") + 
    stat_smooth(method = lm, formula = as.formula(models[[3]]), size = 1, se = FALSE, color = "yellow") + 
    stat_smooth(method = nls, formula = as.formula(models[[4]]), data = dat, method.args = list(start = list(a = 0,b = 0)), size = 1, se = FALSE, color = "red", linetype = 2) + 
    stat_smooth(method = nls, formula = as.formula(models[[5]]), data = dat, method.args = list(start = setNames(coef(lm(y ~ log(x), data = dat)), c("a", "b"))), size = 1, se = FALSE, color = "green", linetype = 2) +
    stat_smooth(method = nls, formula = as.formula(models[[6]]), data = dat, method.args = list(start = list(a = 0,b = 0)), size = 1, se = FALSE, color = "violet") +
    stat_smooth(method = nls, formula = as.formula(models[[7]]), data = dat, method.args = list(start = list(a = 0,b = 0)), size = 1, se = FALSE, color = "orange", linetype = 2)

オレンジ色の曲線がよく見えます。これらのモデルの相対的な適合度を測定すると、どのようにランク付けされるかを見てみましょう...

# calculate the AIC and AICc (for small samples) for each 
# model to see which one is best, ie has the lowest AIC
library(AICcmodavg); library(plyr); library(stringr)
ldply(models, function(mod){ data.frame(AICc = AICc(mod), AIC = AIC(mod), model = deparse(formula(mod))) })

      AICc      AIC                     model
1 70.23024 46.23024                     y ~ x
2 44.37075 20.37075                y ~ I(1/x)
3 67.00075 43.00075                y ~ log(x)
4 43.82083 19.82083    y ~ I(1/x * a) + b * x
5 67.00075 43.00075      y ~ (a + b * log(x))
6 52.75748 28.75748 y ~ I(exp(1)^(a + b * x))
7 44.37075 20.37075        y ~ I(1/x * a) + b

# y ~ I(1/x * a) + b * x is the best model of those tried here for this curve
# it fits nicely on the plot and has the best goodness of fit statistic
# no doubt with a better understanding of nls and the data a better fitting
# function could be found. Perhaps the optimisation method here might be
# useful also: http://stats.stackexchange.com/a/21098/7744

2 番目の方法: 遺伝的プログラミングを使用して膨大な量のモデルを検索する

これは、カーブ フィッティングへの闇のアプローチにおける一種のワイルド ショットのようです。最初に多くを指定する必要はありませんが、おそらく私のやり方が間違っています...

# symbolic regression using Genetic Programming
# http://rsymbolic.org/projects/rgp/wiki/Symbolic_Regression
library(rgp)
# this will probably take some time and throw
# a lot of warnings...
result1 <- symbolicRegression(y ~ x, 
             data=dat, functionSet=mathFunctionSet,
             stopCondition=makeStepsStopCondition(2000))
# inspect results, they'll be different every time...
(symbreg <- result1$population[[which.min(sapply(result1$population, result1$fitnessFunction))]])

function (x) 
tan((x - x + tan(x)) * x) 
# quite bizarre...

# inspect visual fit
ggplot() + geom_point(data=dat, aes(x,y), size = 3) +
  geom_line(data=data.frame(symbx=dat$x, symby=sapply(dat$x, symbreg)), aes(symbx, symby), colour = "red")

実際、視覚的なフィット感は非常に悪いです。おそらく、遺伝的プログラミングから質の高い結果を得るには、もう少し努力が必要です...

クレジット: G. Grothendieckによるカーブ フィッティングの回答 1、カーブ フィッティングの回答 2。

Answer 2

データが従うべき分析関数を知っていますか？もしそうなら、それはあなたがデータに合うように、あなたが関数の形を選ぶのを助けるかもしれません。

それ以外の場合、データは指数関数的減衰のように見えるため、gnuplotで次のようなものを試してください。ここでは、2つの自由パラメーターを持つ関数がデータに適合しています。

 f(x) = exp(-x*c)*b
 fit f(x) "data.dat" u 1:2 via b,c
 plot "data.dat" w p, f(x)

gnuplotは、最適なものにするために、「via」句にちなんで名付けられたパラメーターを変更します。統計は、現在の作業ディレクトリにある「fit.log」というファイルと同様に、stdoutに出力されます。

c変数は曲率（減衰）を決定し、b変数はすべての値を線形にスケーリングしてデータの正しい大きさを取得します。

詳細については、Gnuplotのドキュメントのカーブフィットのセクションを参照してください。