テイラー展開とボンフェロニ補正、積率母関数、ニュートン法
微分積分６テイラー展開、ボンフェロニ補正、積率母関数、ニュートン法
- 作者: ryamada
- 発売日: 2017/01/29
- メディア: Kindle版
- この商品を含むブログを見る
--- title: "Calculus6 Taylor expansion テイラー展開" author: "ryamada" date: "2017年1月22日" output: html_document: toc: true toc_depth: 6 number_section: true --- ```{r setup, include=FALSE} library(knitr) opts_chunk$set(echo = TRUE) library(rgl) knit_hooks$set(rgl = hook_webgl) ``` # 式 Formula 関数が多項式になっている。 A function is expressed as a polynomial function. $$ \sum_{n=0}^\infty \frac{f^{(n)}(a)}{n!} (x-a)^n $$ # Examples $$ f(x) = \sin{x} $$ ```{r} x <- seq(from=-10,to=10,length=1000) fx <- sin(x) plot(x,fx,type="l") ``` $$ f^{(0)}(x) = \sin{x}\\ f^{(1)}(x) = \cos{x}\\ f^{(2)}(x) = - \sin{x}\\ f^{(3)}(x) = - \cos{x}\\ f^{(4)}(x) = \sin{x} = f^{(0)}(x)\\ $$ $$ f^{(4k)}(x) = \sin{x}\\ f^{(4k+1)}(x) = \cos{x}\\ f^{(4k+2)}(x) = - \sin{x}\\ f^{(4k+3)}(x) = - \cos{x}\\ $$ $a=0$のとき、When $a=0$ $$ f^{(4k)}(0) = \sin{0} = 0\\ f^{(4k+1)}(0) = \cos{0} = 1\\ f^{(4k+2)}(0) = - \sin{0} = 0\\ f^{(4k+3)}(0) = - \cos{0} = -1\\ $$ ```{r} fx.0 <- 0/factorial(0) * x^0 plot(x,fx,type="l") points(x,fx.0,type="l",col=2) ``` ```{r} fx.1 <- fx.0 + 1/factorial(1) * x^1 plot(x,fx,type="l") points(x,fx.1,type="l",col=2) ``` ```{r} fx.2 <- fx.1 + 0/factorial(2) * x^2 plot(x,fx,type="l") points(x,fx.2,type="l",col=2) ``` ```{r} fx.3 <- fx.2 + (-1)/factorial(3) * x^3 plot(x,fx,type="l") points(x,fx.3,type="l",col=2) ``` ```{r} fx.4 <- fx.3 + 0/factorial(4) * x^4 plot(x,fx,type="l") points(x,fx.4,type="l",col=2) ``` ```{r} fx.5 <- fx.4 + 0/factorial(5) * x^5 plot(x,fx,type="l") points(x,fx.5,type="l",col=2) ``` ```{r} fx.6 <- fx.5 fx.7 <- fx.6 + 1/factorial(7) * x^7 plot(x,fx,type="l") points(x,fx.7,type="l",col=2) ``` ```{r} fx.8 <- fx.7 fx.9 <- fx.8 + (-1)/factorial(9) * x^9 plot(x,fx,type="l") points(x,fx.9,type="l",col=2) ``` # Exercise 1 ## Exercise 1-1 Do the same for $a=1$ of $\sin{x}$. ## Exercise 1-2 Do the same for $f(x)= e^x$. ## Exercise 1-3 Do the same of $f(x) = (x-1)(x-2)(x-3)(x-4)$ # Multiple testing correction $N$個の独立な検定を行い$N$個のp-値を得るとする。 帰無仮説が成り立つとき、すべてのp-値が a以上である確率は Assume $N$ independent tests. When null hypothesis is true for all, probability that all p-values are equal or more than a is; $$ f(a) = (1-a)^N $$ この式の$a=0$周辺の$f(a)^{(1)}$までのテイラー展開は Taylor expanstion around $a=0$ upto $f(a)^{(1)}$ is; $$ f(a) \sim \frac{1}{0!}a^0 + (-1)\frac{N}{1!}a^1 = 1- Na $$ この近似を用いたのがボンフェロニ補正。 Bonferroni's correction is based on this approximation. # Exercise 2 ## Exercise 2-1 以下の関数を$a=0$周辺で$f(a)^k; k=0,1,2,...$までテイラー展開したときの$k$と近似値との関係をプロットし、ボンフェロニ補正が「保守的」であることを説明せよ。 Expand $f(a)$ around $a=0$ upto $f(a)^k; k=0,1,2,...$ and plot the relation between k and the approximated values, then describe that Bonferroni's correction is conservative using the plot. $$ f(a) = (1-a)^N $$ # 積率母関数 Moment generating function ## Exercise 3 ### Exercise 3-1 このサイト https://www.probabilitycourse.com/chapter6/6_1_3_moment_functions.php を読み、積率母関数とテイラー展開との関係を説明せよ。 Read the site https://www.probabilitycourse.com/chapter6/6_1_3_moment_functions.php . Describe relation between moment generating function and Taylor expansion. # ニュートン法 Newton's method $$ f(x) = \frac{1}{2}e^{x^3} + 2x - 2 $$ ```{r} x <- seq(from=-1.5,to=1.5,length=100) my.fx <- function(x){ exp(x^3)/2 + 2*x - 2 } fx <- my.fx(x) plot(x,fx,type="l") abline(h=0,col=2) ``` $f(x)=0$の解を近似的に求める。 Solve $f(x)=0$ numerically. ## 単調性の確認 Monotonicity # Exercise 4 ## Exercise 4-1 $f(x)$が単調増加関数であることを示せ。 Show that $f(x)$ is a monotonically increasing function. ## ニュートン法 Newton's method 適当な値$x_0$からスタートし、$\frac{d}{dx}f(x)$を使って、$f(x)=0$の解に次第に近づく。 Take an arbitrary value$x_0$ and get closer to the solution of $f(x)=0$ step-by-step. ```{r} x0 <- 1 my.dfxdx <- function(x){ 3/2 * x^2 * exp(x^3) + 2 } this.x <- x0 this.y <- my.fx(this.x) this.dfxdx <- my.dfxdx(this.x) this.intercept <- this.y -this.x*this.dfxdx new.x <- this.x - this.y/this.dfxdx plot(x,fx,type="l") points(this.x,this.y,pch=20,col=3) points(new.x,0,pch=20,col=4) abline(h=0,col=2) abline(this.intercept,this.dfxdx,col=3) x <- seq(from = new.x -0.1,to=this.x+0.1,length=100) fx <- my.fx(x) plot(x,fx,type="l") points(this.x,this.y,pch=20,col=3) points(new.x,0,pch=20,col=4) abline(h=0,col=2) abline(this.intercept,this.dfxdx,col=3) ``` 接線とx軸の交点は$f(x)=0$の解に近づいている。 The crossing point with the $y=0$ line of the tangent line at $(x_0,f(x_0))$ gets closer to the crossing point of $y=f(x)$ with the $y=0$ line is closer that $x_0$. これを繰り返す。Repeat this procedure. ```{r} this.x <- new.x this.y <- my.fx(this.x) this.dfxdx <- my.dfxdx(this.x) this.intercept <- this.y -this.x*this.dfxdx new.x <- this.x - this.y/this.dfxdx plot(x,fx,type="l") points(this.x,this.y,pch=20,col=3) points(new.x,0,pch=20,col=4) abline(h=0,col=2) abline(this.intercept,this.dfxdx,col=3) x <- seq(from = new.x -0.1,to=this.x+0.1,length=100) fx <- my.fx(x) plot(x,fx,type="l") points(this.x,this.y,pch=20,col=3) points(new.x,0,pch=20,col=4) abline(h=0,col=2) abline(this.intercept,this.dfxdx,col=3) ``` ```{r} this.x <- new.x this.y <- my.fx(this.x) this.dfxdx <- my.dfxdx(this.x) this.intercept <- this.y -this.x*this.dfxdx new.x <- this.x - this.y/this.dfxdx plot(x,fx,type="l") points(this.x,this.y,pch=20,col=3) points(new.x,0,pch=20,col=4) abline(h=0,col=2) abline(this.intercept,this.dfxdx,col=3) x <- seq(from = new.x -0.1,to=this.x+0.1,length=100) fx <- my.fx(x) plot(x,fx,type="l") points(this.x,this.y,pch=20,col=3) points(new.x,0,pch=20,col=4) abline(h=0,col=2) abline(this.intercept,this.dfxdx,col=3) ``` ``` }
