このボードに B-Splines のトピックに関する投稿があることは認識していますが、実際にはそれらの投稿でさらに混乱しているので、誰かが私を助けてくれるかもしれないと思っていました.
0 から 1 の範囲の x 値のシミュレートされたデータがありますdegree = 3。0、0.1、0.2、...、0.9、1 のノットを持つ 3 次スプライン ( ) をデータに当てはめたいと思います。パラメータ推定に B-Spline 基底と OLS を使用します (ペナルティ付きスプラインは探していません)。
bsパッケージの機能が必要だと思いますsplineが、よくわかりませんし、正確に何をフィードするのかもわかりません。
結果の多項式スプラインもプロットしたいと思います。
ありがとう！
## simulate some data - from mgcv::magic
set.seed(1)
n <- 400
x <- 0:(n-1)/(n-1)
f <- 0.2*x^11*(10*(1-x))^6+10*(10*x)^3*(1-x)^10
y <- f + rnorm(n, 0, sd = 2)
## load the splines package - comes with R
require(splines)
OLS の見積もりが必要な場合bs()は、数式で関数を使用します。ノット、多項式の次数などによって与えられる基底関数を提供します。lmbs
mod <- lm(y ~ bs(x, knots = seq(0.1, 0.9, by = 0.1)))
それを線形モデルのように扱うことができます。
> anova(mod)
Analysis of Variance Table
Response: y
Df Sum Sq Mean Sq F value Pr(>F)
bs(x, knots = seq(0.1, 0.9, by = 0.1)) 12 2997.5 249.792 65.477 < 2.2e-16 ***
Residuals 387 1476.4 3.815
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
結び目の配置に関するいくつかの指針。にはデフォルトbsの引数があります。したがって、上記の引数を指定したときに、境界ノットを含めませんでした。Boundary.knotsBoundary.knots = range(x)knots
詳細?bsについては、こちらをお読みください。
コメントでは、適合したスプラインを描画する方法について説明しています。1 つのオプションは、共変量に関してデータを並べ替えることです。これは、共変量が 1 つの場合は問題なく機能しますが、共変量が 2 つ以上の場合は機能する必要はありません。さらなる問題は、観測された値でのみ近似スプラインを評価できることです。xこれは、共変量を密にサンプリングした場合は問題ありませんが、そうでない場合、長い線形セクションでスプラインが奇妙に見える場合があります。
より一般的な解決策はpredict、共変量または共変量の新しい値のモデルから予測を生成するために使用することです。以下のコードでは、上記のモデルに対してこれを行う方法を示し、 の範囲で 100 個の等間隔の値を予測しますx。
pdat <- data.frame(x = seq(min(x), max(x), length = 100))
## predict for new `x`
pdat <- transform(pdat, yhat = predict(mod, newdata = pdat))
## now plot
ylim <- range(pdat$y, y) ## not needed, but may be if plotting CIs too
plot(y ~ x)
lines(yhat ~ x, data = pdat, lwd = 2, col = "red")
それが生み出す
回答の例に基づいて、適合したスプラインをプロットするより簡単な方法は、effectsパッケージを使用することです。
## simulate some data - from mgcv::magic
set.seed(1)
n <- 400
x <- 0:(n-1)/(n-1)
f <- 0.2*x^11*(10*(1-x))^6+10*(10*x)^3*(1-x)^10
y <- f + rnorm(n, 0, sd = 2)
## load the splines package - comes with R
require(splines)
require(car)
require(effects)
## estimate model
mod <- lm(y ~ bs(x, knots = seq(0.1, 0.9, by = 0.1)))
Anova次に、から使用できますcar：
> Anova(mod)
Anova Table (Type II tests)
Response: y
Sum Sq Df F value Pr(>F)
bs(x, knots = seq(0.1, 0.9, by = 0.1)) 2997.5 12 65.477 < 2.2e-16 ***
Residuals 1476.4 387
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
effectsまた、パッケージを使用して、適合したスプラインを簡単にプロットできます。
plot(allEffects(mod))
これは次のように出力されます:
以下も参照してください。