用R语言实现牛顿迭代求最优点

牛顿迭代的基本原理见百科http://baike.sogou.com/v1846821.htm?fromTitle=%E7%89%9B%E9%A1%BF%E8%BF%AD%E4%BB%A3%E6%B3%95
基本公式：

泰勒展开公式：f(x) = f(x0)+f'(x0)*(x-x0)+1/2*f''(x0)*(x-x0)^2

                   →f'(x) = f'(x0)+f''(x0)*(x-x0)  #一阶导为零是极值

                   →x = x0 - f'(x0)/f''(x0)         #公式

下面是R语言的代码实现 函数为y = x^2+2

x = 5  #初始值，自己定
diedai   = 40 #迭代次数
g        = function(x){x^2+2} #目标函数
g.prime  = function(x){2*x}    #一阶导
g.2prime = function(x){2}      #二阶导

#公式
for(i in 1:diedai){x = x - g.prime(x)/g.2prime(x)} 
x
 

下面是结果

> x
[1] 5
>  
> 
> x = 5  #初始值，自己定
> diedai   = 40 #迭代次数
> g        = function(x){x^2+2} #目标函数
> g.prime  = function(x){2*x}    #一阶导
> g.2prime = function(x){2}      #二阶导
> 
> #公式
> for(i in 1:diedai){x = x - g.prime(x)/g.2prime(x)} 
> x
[1] 0

---------------------------------------------------------------------------------------------------------

下面是复杂函数的牛顿迭代求最优点

x = 5  #初始值，自己定
diedai   = 40 #迭代次数
g        = function(x){sin(x)^3+log(cos(x))+x^(-1)} #目标函数
x=5
g.prime  = function(x){
fun = expression(sin(x)^3+log(cos(x))+x^(-1))
D(fun,'x')
f = deriv(D(fun,'x'),'x',func=T)
f(x)}    #一阶导
g.2prime = function(x){
D(D(fun,'x'),'x')
f = deriv(D(D(fun,'x'),'x'),'x',func=T)
f(x)}      #二阶导

#公式
for(i in 1:diedai){x = x - g.prime(x)/g.2prime(x)} 
x

--------------------------------------------------------------------------------------------------------

 

转载于:https://my.oschina.net/u/3272626/blog/856173