<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/numerical-differentiation-with-finite-differences-in-r/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="http://www.aaronschlegel.com/numerical-differentiation-finite-differences-r/"> R – Aaron Schlegel</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<div>Part 1 of 7 in the series <a rel="nofollow" href="http://www.aaronschlegel.com/series/numerical-analysis/" title="Numerical Analysis" target="_blank">Numerical Analysis</a></div>
<p>Numerical differentiation is a method of approximating the derivative of a function <img alt="f" title="f" src="https://s0.wp.com/latex.php?latex=f&bg=ffffff&fg=000&s=0" class="latex"> at particular value <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex">. Often, particularly in physics and engineering, a function may be too complicated to merit the work necessary to find the exact derivative, or the function itself is unknown, and all that is available are some points <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex"> and the function evaluated at those points. Numerical differentiation, of which finite differences is just one approach, allows one to avoid these complications by approximating the derivative.</p>
<p>The most straightforward and simple approximation of the first derivative is defined as:</p>
<img alt=" f^\prime (x) \approx \frac{f(x + h) - f(x)}{h} \qquad h &gt; 0 " title=" f^\prime (x) \approx \frac{f(x + h) - f(x)}{h} \qquad h &gt; 0 " src="https://s0.wp.com/latex.php?latex=+f%5E%5Cprime+%28x%29+%5Capprox+%5Cfrac%7Bf%28x+%2B+h%29+-+f%28x%29%7D%7Bh%7D+%5Cqquad+h+%3E+0+&bg=ffffff&fg=000&s=0" class="latex">
<p>The approximation error of this equation can be found by performing a Taylor expansion of <img alt="f(x + h)" title="f(x + h)" src="https://s0.wp.com/latex.php?latex=f%28x+%2B+h%29&bg=ffffff&fg=000&s=0" class="latex"> about <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex">, which gives:</p>
<img alt=" f(x + h) = f(x) + hf^\prime (x) + \frac{h^2}{2} f^{\prime\prime}(\epsilon) " title=" f(x + h) = f(x) + hf^\prime (x) + \frac{h^2}{2} f^{\prime\prime}(\epsilon) " src="https://s0.wp.com/latex.php?latex=+f%28x+%2B+h%29+%3D+f%28x%29+%2B+hf%5E%5Cprime+%28x%29+%2B+%5Cfrac%7Bh%5E2%7D%7B2%7D+f%5E%7B%5Cprime%5Cprime%7D%28%5Cepsilon%29+&bg=ffffff&fg=000&s=0" class="latex">
<p>We then rearrange the expansion and substitute <img alt="f^\prime" title="f^\prime" src="https://s0.wp.com/latex.php?latex=f%5E%5Cprime&bg=ffffff&fg=000&s=0" class="latex"> with the earlier approximation to arrive at an exact form of the approximation equation:</p>
<img alt=" f^\prime (x) = \frac{f(x + h) - f(x)}{h} - \frac{h^2}{2} f^{\prime\prime}(\epsilon) " title=" f^\prime (x) = \frac{f(x + h) - f(x)}{h} - \frac{h^2}{2} f^{\prime\prime}(\epsilon) " src="https://s0.wp.com/latex.php?latex=+f%5E%5Cprime+%28x%29+%3D+%5Cfrac%7Bf%28x+%2B+h%29+-+f%28x%29%7D%7Bh%7D+-+%5Cfrac%7Bh%5E2%7D%7B2%7D+f%5E%7B%5Cprime%5Cprime%7D%28%5Cepsilon%29+&bg=ffffff&fg=000&s=0" class="latex">
<h6>Finite Differences Methods for Numerical Differentiation</h6>
<p>There are three primary differencing techniques, forward, backward, and central, for approximating the derivative of a function at a particular point. Note there are more accurate methods and algorithms for approximating derivatives; however, due to their relative ease of computation and accuracy, differencing methods are still a viable tool for performing numerical differentiation.</p>
<p>The forward difference approximation is the same as the earlier approximation:</p>
<img alt=" f^\prime (x) \approx \frac{f(x + h) - f(x)}{h} " title=" f^\prime (x) \approx \frac{f(x + h) - f(x)}{h} " src="https://s0.wp.com/latex.php?latex=+f%5E%5Cprime+%28x%29+%5Capprox+%5Cfrac%7Bf%28x+%2B+h%29+-+f%28x%29%7D%7Bh%7D+&bg=ffffff&fg=000&s=0" class="latex">
<p>The backward difference approximation is based on the values of the function evaluated at <img alt="x - h" title="x - h" src="https://s0.wp.com/latex.php?latex=x+-+h&bg=ffffff&fg=000&s=0" class="latex"> and <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex">, defined as:</p>
<img alt=" f^\prime \approx \frac{f(x) - f(x - h)}{h} " title=" f^\prime \approx \frac{f(x) - f(x - h)}{h} " src="https://s0.wp.com/latex.php?latex=+f%5E%5Cprime+%5Capprox+%5Cfrac%7Bf%28x%29+-+f%28x+-+h%29%7D%7Bh%7D+&bg=ffffff&fg=000&s=0" class="latex">
<p>The central (or centered) difference approximation is essentially an average of the forward and backward differencing methods:</p>
<img alt=" f^\prime \approx \frac{f(x + h) - f(x - h)}{2h} " title=" f^\prime \approx \frac{f(x + h) - f(x - h)}{2h} " src="https://s0.wp.com/latex.php?latex=+f%5E%5Cprime+%5Capprox+%5Cfrac%7Bf%28x+%2B+h%29+-+f%28x+-+h%29%7D%7B2h%7D+&bg=ffffff&fg=000&s=0" class="latex">
<h6>Error of Approximations</h6>
<p>It was shown earlier the forward and backward difference approximations have an error term of <img alt="\frac{h^2}{2} f^{\prime\prime}(\epsilon)" title="\frac{h^2}{2} f^{\prime\prime}(\epsilon)" src="https://s0.wp.com/latex.php?latex=%5Cfrac%7Bh%5E2%7D%7B2%7D+f%5E%7B%5Cprime%5Cprime%7D%28%5Cepsilon%29&bg=ffffff&fg=000&s=0" class="latex"> where <img alt="\epsilon \in (x, x + h)" title="\epsilon \in (x, x + h)" src="https://s0.wp.com/latex.php?latex=%5Cepsilon+%5Cin+%28x%2C+x+%2B+h%29&bg=ffffff&fg=000&s=0" class="latex">. The error term of the central difference method can be found by performing Taylor expansions on <img alt="f(x + h)" title="f(x + h)" src="https://s0.wp.com/latex.php?latex=f%28x+%2B+h%29&bg=ffffff&fg=000&s=0" class="latex"> and <img alt="f(x - h)" title="f(x - h)" src="https://s0.wp.com/latex.php?latex=f%28x+-+h%29&bg=ffffff&fg=000&s=0" class="latex"> about <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex">.</p>
<p><img alt=" f(x + h) = f(x) + hf^\prime (x) + \frac{h^2}{2} f^{\prime\prime}(x) + \frac{h^3}{6} f^{\prime\prime\prime}(\epsilon) \qquad \epsilon \in (x, x + h) " title=" f(x + h) = f(x) + hf^\prime (x) + \frac{h^2}{2} f^{\prime\prime}(x) + \frac{h^3}{6} f^{\prime\prime\prime}(\epsilon) \qquad \epsilon \in (x, x + h) " src="https://s0.wp.com/latex.php?latex=+f%28x+%2B+h%29+%3D+f%28x%29+%2B+hf%5E%5Cprime+%28x%29+%2B+%5Cfrac%7Bh%5E2%7D%7B2%7D+f%5E%7B%5Cprime%5Cprime%7D%28x%29+%2B+%5Cfrac%7Bh%5E3%7D%7B6%7D+f%5E%7B%5Cprime%5Cprime%5Cprime%7D%28%5Cepsilon%29+%5Cqquad+%5Cepsilon+%5Cin+%28x%2C+x+%2B+h%29+&bg=ffffff&fg=000&s=0" class="latex"><br>
<img alt=" f(x - h) = f(x) - hf^\prime (x) + \frac{h^2}{2} f^{\prime\prime}(x) - \frac{h^3}{6} f^{\prime\prime\prime}(\epsilon) \qquad \epsilon \in (x - h, x) " title=" f(x - h) = f(x) - hf^\prime (x) + \frac{h^2}{2} f^{\prime\prime}(x) - \frac{h^3}{6} f^{\prime\prime\prime}(\epsilon) \qquad \epsilon \in (x - h, x) " src="https://s0.wp.com/latex.php?latex=+f%28x+-+h%29+%3D+f%28x%29+-+hf%5E%5Cprime+%28x%29+%2B+%5Cfrac%7Bh%5E2%7D%7B2%7D+f%5E%7B%5Cprime%5Cprime%7D%28x%29+-+%5Cfrac%7Bh%5E3%7D%7B6%7D+f%5E%7B%5Cprime%5Cprime%5Cprime%7D%28%5Cepsilon%29+%5Cqquad+%5Cepsilon+%5Cin+%28x+-+h%2C+x%29+&bg=ffffff&fg=000&s=0" class="latex"></p>
<p>These two equations are then subtracted, and then solving for <img alt="f^\prime (x)" title="f^\prime (x)" src="https://s0.wp.com/latex.php?latex=f%5E%5Cprime+%28x%29&bg=ffffff&fg=000&s=0" class="latex"> leads to the central difference method:</p>
<img alt=" f^\prime \approx \frac{f(x + h) - f(x - h)}{2h} - \frac{h^2}{12}\bigg({f^{\prime\prime\prime}(\epsilon_1) + f^{\prime\prime\prime}(\epsilon_2)\bigg)} " title=" f^\prime \approx \frac{f(x + h) - f(x - h)}{2h} - \frac{h^2}{12}\bigg({f^{\prime\prime\prime}(\epsilon_1) + f^{\prime\prime\prime}(\epsilon_2)\bigg)} " src="https://s0.wp.com/latex.php?latex=+f%5E%5Cprime+%5Capprox+%5Cfrac%7Bf%28x+%2B+h%29+-+f%28x+-+h%29%7D%7B2h%7D+-+%5Cfrac%7Bh%5E2%7D%7B12%7D%5Cbigg%28%7Bf%5E%7B%5Cprime%5Cprime%5Cprime%7D%28%5Cepsilon_1%29+%2B+f%5E%7B%5Cprime%5Cprime%5Cprime%7D%28%5Cepsilon_2%29%5Cbigg%29%7D+&bg=ffffff&fg=000&s=0" class="latex">
<p>The central difference approximation is more accurate than forward and backward differences and should be used whenever possible. The three difference methods report the same approximations of the following example as the function and its derivative are rather simple; however, it is still best to apply the central difference approximation in actual practice.</p>
<h6>Differencing Approximations Example in R</h6>
<p>Consider the following set of data points:</p>
<table>
<thead>
<tr>
<th><img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex"></th>
<th><img alt="f(x)" title="f(x)" src="https://s0.wp.com/latex.php?latex=f%28x%29&bg=ffffff&fg=000&s=0" class="latex"></th>
</tr>
</thead>
<tbody>
<tr>
<td>0.0</td>
<td>0.00000</td>
</tr>
<tr>
<td>0.2</td>
<td>0.74140</td>
</tr>
<tr>
<td>0.4</td>
<td>1.37180</td>
</tr>
</tbody>
</table>
<pre>
x &lt;- c(0.0, 0.2, 0.4)
fx &lt;- c(0.00000, 0.74140, 1.3718)
</pre>
<p>The function is unknown; however, we would still like to approximate the function’s derivative at the values of <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex">. The following function combines the forward and backward differencing methods to approximate the derivative of the function at each value of <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex">.</p>
<pre>
finite.differences &lt;- function(x, y) {
  if (length(x) != length(y)) {
    stop('x and y vectors must have equal length')
  }
  
  n &lt;- length(x)
  
  # Initialize a vector of length n to enter the derivative approximations
  fdx &lt;- vector(length = n)
  
  # Iterate through the values using the forward differencing method
  for (i in 2:n) {
    fdx[i-1] &lt;- (y[i-1] - y[i]) / (x[i-1] - x[i])
  }
  
  # For the last value, since we are unable to perform the forward differencing method 
  # as only the first n values are known, we use the backward differencing approach
  # instead. Note this will essentially give the same value as the last iteration 
  # in the forward differencing method, but it is used as an approximation as we 
  # don't have any more information
  fdx[n] &lt;- (y[n] - y[n - 1]) / (x[n] - x[n - 1])
  
  return(fdx)
}
</pre>
<p>Using the function, perform each differencing method on the <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex"> values and the evaluated points <img alt="f(x)" title="f(x)" src="https://s0.wp.com/latex.php?latex=f%28x%29&bg=ffffff&fg=000&s=0" class="latex">.</p>
<pre>
finite &lt;- finite.differences(x, fx)
finite
## [1] 3.707 3.152 3.152
</pre>
<p>Let’s assume the data from earlier was taken from the function <img alt="e^x - 2x^2 + 3x - 1" title="e^x - 2x^2 + 3x - 1" src="https://s0.wp.com/latex.php?latex=e%5Ex+-+2x%5E2+%2B+3x+-+1&bg=ffffff&fg=000&s=0" class="latex">.</p>
<pre>
f &lt;- function(x) {
  return(exp(x) - 2 * x^2 + 3 * x - 1)
}
</pre>
<p>Using the central differencing method and the given function, we can approximate the derivative at each point. As <img alt="h" title="h" src="https://s0.wp.com/latex.php?latex=h&bg=ffffff&fg=000&s=0" class="latex"> approaches 0, the approximation becomes more accurate. The following function approximates the values of the derivative at the given points at several different values of <img alt="h" title="h" src="https://s0.wp.com/latex.php?latex=h&bg=ffffff&fg=000&s=0" class="latex"> to demonstrate how the approximations converge.</p>
<pre>
central.difference &lt;- function(f, x) {
  steps &lt;- c(0.1, 0.01, 0.001, 0.0001, 0.00001, 0.000001, 0.0000001)
  n &lt;- length(steps)
  
  fdx &lt;- vector(length = n)
  
  for (h in 1:n) {
    fdx[h] &lt;- (f(x + 0.5 * steps[h]) - f(x - 0.5 * steps[h])) / steps[h]
  }
  
  return(fdx)
}
</pre>
<p>Print the approximations of the derivative.</p>
<pre>
for (i in x) {
  print(central.difference(f, i))
}
## [1] 4.000417 4.000004 4.000000 4.000000 4.000000 4.000000 4.000000
## [1] 3.421912 3.421408 3.421403 3.421403 3.421403 3.421403 3.421403
## [1] 2.892446 2.891831 2.891825 2.891825 2.891825 2.891825 2.891825
</pre>
<p>We see the approximations converge rather quickly as <img alt="h" title="h" src="https://s0.wp.com/latex.php?latex=h&bg=ffffff&fg=000&s=0" class="latex"> approaches 0, giving the following approximations:</p>
<img alt=" f^\prime (0.0) \approx 4.00000 \qquad f^\prime (0.2) \approx 3.421403 \qquad f^\prime (0.4) \approx 2.891825 " title=" f^\prime (0.0) \approx 4.00000 \qquad f^\prime (0.2) \approx 3.421403 \qquad f^\prime (0.4) \approx 2.891825 " src="https://s0.wp.com/latex.php?latex=+f%5E%5Cprime+%280.0%29+%5Capprox+4.00000+%5Cqquad+f%5E%5Cprime+%280.2%29+%5Capprox+3.421403+%5Cqquad+f%5E%5Cprime+%280.4%29+%5Capprox+2.891825+&bg=ffffff&fg=000&s=0" class="latex">
<p>The derivative of the function is <img alt="e^x - 4x + 3" title="e^x - 4x + 3" src="https://s0.wp.com/latex.php?latex=e%5Ex+-+4x+%2B+3&bg=ffffff&fg=000&s=0" class="latex">. We shall compare our approximated values of with the actual values of the derivative at <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex"> to see how the central differencing method faired.</p>
<pre>
fdx &lt;- function(x) {
  return(exp(x) - 4 * x + 3)
}
</pre>
<p>For fun, plot both the function and its derivative to get a visualization of where the function’s derivative is at the values of <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex">.</p>
<pre>
ggplot(data = data.frame(x = 0), mapping = aes(x = x)) + 
  stat_function(fun = f, size = 1.25, alpha = 0.75) +
  stat_function(fun = fdx, size = 1.25, color = 'blue', alpha = 0.75) +
  xlim(-3,3)
</pre>
<p><a rel="nofollow" href="http://www.aaronschlegel.com/?attachment_id=22592" target="_blank"><img sizes="(max-width: 1060px) 100vw, 1060px" src="https://i2.wp.com/www.aaronschlegel.com/wp-content/uploads/func_deriv-1.png?resize=1060%2C757" alt="Function and its derivative to approximate with finite differences" class="wp-image-22592" data-recalc-dims="1"></a></p>
<p>Print the calculated derivative’s values for each <img alt="x" title="x" src="https://s0.wp.com/latex.php?latex=x&bg=ffffff&fg=000&s=0" class="latex"> and compare them with our previously approximated values.</p>
<pre>
actual &lt;- vector(length = length(x))
central.approx &lt;- c(4.00000, 3.421403, 2.891825)

for (i in 1:length(x)) {
  actual[i] &lt;- fdx(x[i])
}

approx.df &lt;- data.frame(cbind(actual, central.approx, actual - central.approx, finite, actual - finite))
colnames(approx.df) &lt;- c('Actual Values', 'Central Difference Approx', 'Central Differences Error', 'Finite Differences', 'Finite Differences Error')

approx.df
##   Actual Values Central Difference Approx Central Differences Error
## 1      4.000000                  4.000000              0.000000e+00
## 2      3.421403                  3.421403             -2.418398e-07
## 3      2.891825                  2.891825             -3.023587e-07
##   Finite Differences Finite Differences Error
## 1              3.707                0.2930000
## 2              3.152                0.2694028
## 3              3.152               -0.2601753
</pre>
<h6>References</h6>
<p>Burden, R. L., &amp; Faires, J. D. (2011). Numerical analysis (9th ed.). Boston, MA: Brooks/Cole, Cengage Learning.</p>
<p>Weisstein, Eric W. “Numerical Differentiation.” From MathWorld–A Wolfram Web Resource. <a rel="nofollow" href="http://mathworld.wolfram.com/NumericalDifferentiation.html" target="_blank">http://mathworld.wolfram.com/NumericalDifferentiation.html</a></p>
<p>The post <a rel="nofollow" href="http://www.aaronschlegel.com/numerical-differentiation-finite-differences-r/" target="_blank">Numerical Differentiation with Finite Differences in R</a> appeared first on <a rel="nofollow" href="http://www.aaronschlegel.com/" target="_blank">Aaron Schlegel</a>.</p>
	<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="http://www.aaronschlegel.com/numerical-differentiation-finite-differences-r/"> R – Aaron Schlegel</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a rel="nofollow" href="https://www.r-project.org/" title="The R Project for Statistical Computing">R</a> news and <a rel="nofollow" href="https://www.r-bloggers.com/search/tutorial" title="R tutorials">tutorials</a> on topics such as: <a rel="nofollow" href="https://www.r-bloggers.com/search/data%20science" title="Data science">Data science</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Big%20Data" title="Big Data">Big Data, <a rel="nofollow" href="https://www.r-users.com/" title="R jobs">R jobs</a>, visualization (<a rel="nofollow" href="https://www.r-bloggers.com/search/ggplot2" title="ggplot and ggplot2 tutorials">ggplot2</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/boxplot" title="Boxplots using lattice and ggplot2 tutorials">Boxplots</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/map" title="Maps and gis">maps</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/animation" title="Animation in R">animation</a>), programming (<a rel="nofollow" href="https://www.r-bloggers.com/search/RStudio" title="RStudio IDE for R">RStudio</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/sweave" title="Sweave and literate programming">Sweave</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX" title="LaTeX in R">LaTeX</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/SQL" title="SQL and databases">SQL</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/eclipse" title="Eclipse IDE for R">Eclipse</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/git" title="git and github, Version Control System">git</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/hadoop" title="Large data in R using Hadoop">hadoop</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping" title="Web Scraping of google, facebook, yahoo, twitter and more using R">Web Scraping</a>) statistics (<a rel="nofollow" href="https://www.r-bloggers.com/search/regression" title="Regressions and ANOVA analysis tutorials">regression</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/PCA" title="principal component analysis tutorial">PCA</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/time+series" title="Time series">time series</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/trading" title="finance trading">trading</a>) and more...
<img width="1" alt="" src="http://feeds.feedburner.com/~r/RBloggers/~4/On20C9Sm5bU" height="1"></a></div>