Newton-Raphson算法简介及其R实现

本文简要介绍了Newton-Raphson方法及其R语言实现并给出几道练习题供参考使用。 下载PDF格式文档(Academia.edu)

  • Newton-Raphson Method
    Let $f(x)$ be a differentiable function and let $a_0$ be a guess for a solution to the equation $$f(x)=0$$ We can product a sequence of points $x=a_0, a_1, a_2, dots $ via the recursive formula $$a_{n+1}=a_n-frac{f(a_n)}{f'(a_n)}$$ that are successively better approximation of a solution to the equation $f(x)=0$.
  • R codes

    There are 4 parameters in this function:
    • f is the function you input.
    • tol is the tolerance (default $1e-7$).
    • x0 is the initial guess.
    • N is the default number (100) of iterations.
    The process will be end up until either the absolute difference between two adjacent approximations is less than tol, or the number of iterations reaches N.
  • Examples
    Generally speaking, the "guess" is important. More precisely, according to Intermediate Value Theorem we can find two values of which function value are larger and less than 0, respectively. Then choosing the one, which first derivative is larger than another, as the initial guess value in the iterative formula. This process will guarantee the convergence of roots. Let's see some examples.
    • Example 1
      Approximate the fifth root of 7.
      Solution:
      Denote $f(x)=x^5-7$. It is easily to know that $f(1)=-6 < 0$ and $f(2)=25 > 0$. Additionally, $f'(1)=5 < f'(2)=80$, so we set the initial guess value $x_0=2$. By Newton-Raphson method we get the result is 1.47577316159. And $$f(1.47577316159)approx 1.7763568394e-15$$ which is very close to 0. R codes is below:
      # Example 1
      f = function(x){x^5 - 7}
      h = 1e - 7
      df.dx = function(x){(f(x + h) - f(x)) / h}
      df.dx(1); df.dx(2)
      # [1] 5.0000009999
      # [1] 80.0000078272
      app = newton(f, x0 = 2)
      app
      # [1] 1.68750003057 1.52264459615 1.47857137506 1.47578373325 1.47577316175
      # [6] 1.47577316159
      f(app[length(app)])
      # [1] 1.7763568394e-15
    • Example 2
      The function $f(x)=x^5-5x^4+5x^2-6$ has a root between 1 and 5. Approximate it by Newton-Raphson method.
      Solution:
      We try to calculate some values first. $f(1)=-5, f(2)=-34, f(3)=-123, f(4)=-182, f(5)=119$, so there should be a root between 4 and 5. Since $f'(4)=40 < f'(5)=675$, hence $x_0=5$ is a proper initial guess value. By Newton-Raphson method we get the result is 4.79378454069 and $$f(4.79378454069)approx -2.84217094304e-14$$ which is a desired approximation. R codes is below:
      # Example 2
      f = function(x){x^5 - 5 * x^4 + 5 * x^2 - 6}
      x = c(1 : 5)
      f(x)
      # [1]   -5  -34 -123 -182  119
      h = 1e-7
      df.dx = function(x){(f(x + h) - f(x)) / h}
      df.dx(4); df.dx(5)
      # [1] 40.0000163836
      # [1] 675.000053008
      app = newton(f, x0 = 5)
      app
      # [1] 4.82370371755 4.79453028339 4.79378501861 4.79378454069 4.79378454069
      f(app[length(app)])
      # [1] -2.84217094304e-14
    • Example 3
      A rectangular piece of cardboard of dimensions $8 imes 17$ is used to make an open-top box by cutting out a small square of side $x$ from each corner and bending up the sides. Find a value of $x$ for which the box has volume 100.
      Solution:
      Firstly, building the model. $V(x)=x(8-2x)(17-2x)=100$, that is, we want to find the root of equation $$f(x)=x(8-2x)(17-2x)-100=0Leftrightarrow f(x)=4x^3-50x^2+136x-100=0$$ We know that $0 < x < 4$ and hence try to calculate some non-negative integers: $$f(0)=-100, f(1)=-10, f(2)=4, f(3)=-34, f(4)=-100$$ Note that there are two intervals may have roots: $(1, 2)cup (2,3)$. Since $$f'(1)=48 > f'(2)=-16 > f'(3)=-56$$ so we set the initial guess values $x_0=1$ and $x'_0=2$ (i.e. there are two separate iteration procedures). By using Newton-Raphson method we obtain the result are 11.26063715644 and 2.19191572127 respectively. Both of them are quite accurate. R codes is below:
      # Example 3
      f = function(x){4 * x^3 - 50 * x^2 + 136 * x - 100}
      x = c(0 : 4)
      f(x)
      # [1] -100  -10    4  -34 -100
      h = 1e-7
      df.dx = function(x){(f(x + h) - f(x)) / h}
      df.dx(1); df.dx(2); df.dx(3)
      # [1] 47.9999962977
      # [1] -16.0000024607
      # [1] -56.0000012229
      app1 = newton(f, x0 = 1)
      app2 = newton(f, x0 = 2)
      app1; app2
      # [1] 1.20833334940 1.25768359879 1.26062673622 1.26063715631 1.26063715644
      # [1] 2.24999996155 2.19469026652 2.19192282154 2.19191572132 2.19191572127
      f(app1[length(app1)]); f(app2[length(app2)])
      # [1] 2.84217094304e-14
      # [1] -2.84217094304e-14
      


作者:赵胤
出处:http://www.cnblogs.com/zhaoyin/
本文版权归作者和博客园共有,欢迎转载,但未经作者同意必须保留此段声明,且在文章页面明显位置给出原文连接,否则保留追究法律责任的权利。

原文地址:https://www.cnblogs.com/zhaoyin/p/4125450.html