LoginSignup
0
1

More than 5 years have passed since last update.

Rで多項式曲線フィッティング

Last updated at Posted at 2013-01-21

PRML演習1.1に従って係数を求め、図1.3 多項式曲線フィッティングと同様の図を描きます。

1.PolynomialFitting.png

set.seed(0)
x<-seq(0.05,1,0.1)
t<-sin(2*pi*x)+rnorm(10,0,0.5)
makeA<-function(m){
    A<-matrix(nrow=m+1,ncol=m+1)
    for (i in 0:m) {
        for (j in 0:m) {
            A[1+i,1+j]<-sum(x^(i+j))
        }
    }
    A
}
makeB<-function(m){
    b<-0:m
    for (i in 0:m) {
        b[i+1]<-sum(x^i*t)
    }
    b
}
xrange<-c(-0.5,1.5)
yrange<-c(-2,2)
plot(x,t,xlim=xrange,ylim=yrange)
mrange=1:9
for (M in mrange) {
    A<-makeA(M)
    b<-makeB(M)
    w<-solve(A,b)
    estimate<-function(x){
        i <- 0:(length(w)-1)
        (sapply(x,function(xn){ (sum(w*xn^i)) })) # (w[1]+w[2]*x+w[3]*x^2+w[4]*x^3...)
    }
    par(new=T)
    curve(estimate,type="l",xlim=xrange,ylim=yrange,col=M)
}
legend("topleft",legend=mrange,col=mrange,lty=1,cex=0.5)

多項式$y(x,\mathbf{w}) = w_0 + w_1 x + w_2 x^2 + \ldots + w_M x^M = \sum_{j=0}^{M} w_j x^j$に観測値をフィッティングする。目標値(観測値)$t_n$、予測値$y(x_n,\mathbf{w})$として、最適な$\mathbf{w}$は、次の二乗和誤差を最小化するもの。

{\displaystyle E(\mathbf{w}) = \frac{1}{2} \sum_{n=1}^N \Bigl(y(x_n,\mathbf{w}) - t_n\Bigr)^2}

二乗和誤差の$\mathbf{w}$についての微分が$\mathbf{0}$となる$\mathbf{w}$を求める。

\begin{aligned}

\frac{\partial}{\partial w_i}E(\mathbf{w}) &= \frac{1}{2} \cdot 2 \sum_{n=1}^N \Bigl\{y(x_n,\mathbf{w}) - t_n\Bigr\}\frac{\partial}{\partial w_i}y(x_n,\mathbf{w})\\

&= \sum_{n=1}^N \biggl( \Bigl(\sum_{j=0}^M w_j x_n^j -t_n\Bigr)\frac{\partial}{\partial w_i}\sum_{j=0}^M w_j x_n^j\biggr)
\\

&= \sum_{n=1}^N \biggl( \Bigl(\sum_{j=0}^M w_j x_n^j -t_n\Bigr) x_n^i \biggr)\\

&= \sum_{n=1}^N \Bigl(x_n^i \sum_{j=0}^M w_j x_n^j - x_n^i t_n\Bigr) = 0\\

\sum_{n=1}^N \Bigl(x_n^i \sum_{j=0}^M w_j x_n^j \Bigr) &= \sum_{n=1}^N x_n^i t_n\\

\sum_{j=0}^M w_j \sum_{n=1}^N x_n^{i+j} &= \sum_{n=1}^N x_n^i t_n

\end{aligned}

$\mathbf{w}$についての連立方程式$\mathbf{A w} = \mathbf{b}$を解く。

\begin{aligned}
\mathbf{A}&=\sum_{n=1}^N \begin{bmatrix} x_n^0 && x_n^1 && \cdots && x_n^M \\ \vdots && \vdots && \ddots && \vdots \\ x_n^M && x_n^{M+1} && \cdots && x_n^{M+M}\end{bmatrix}\\

\mathbf{w}&=\begin{bmatrix}w_0 \\ \vdots \\ w_M\end{bmatrix}\\

\mathbf{b}&=\sum_{n=1}^N \begin{bmatrix} x_n^0 t_n \\ \vdots \\ x_n^M t_n\end{bmatrix}
\end{aligned}

Dropbox にも資料を置いています。

0
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
1