Note

This notebook can be downloaded here: 03_Numerical_Derivation.ipynb

Derivation of numerical data¶

Code author: Emile Roux emile.roux@univ-smb.fr

RISE Slideshow

Scope¶

A finite number \(N\) of data points \((x,y)\) are available, its discribe a function \(f(x)\): compute the derivative \(\dfrac{df}{dx}\)

https://en.wikipedia.org/wiki/Numerical_differentiation

#setup
%load_ext autoreload
%matplotlib nbagg
%autoreload 2
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib as mpl

Create synthetic data¶

For this notebook we use data comming from a known function. In this way we can check the accuracy of the results.

First function

# The function
def f(x):
    return (x**2-3*x+2*np.sin(5*x))

# The derivative
def df(x):
    return 2*x - 3 + 10 * np.cos(5*x)

Seconde function

# The function
def f(x):
    return (x**2)

# The derivative
def df(x):
    return 2*x

Plotting the function and its dericative:

N = 400
xmin, xmax = 0., 4
x = np.linspace(xmin, xmax, N)

fig = plt.figure()
plt.plot(x,f(x),'b',label = "$f(x)$")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>

fig = plt.figure()
plt.plot(x,df(x),'b',label = "$\dfrac{df}{dx}$")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>

Now we assume that only 40 points of this function are known:¶

N = 40
xmin, xmax = 0., 4
xi = np.linspace(xmin, xmax, N)
fi = f(xi)

fig = plt.figure()
plt.plot(xi,fi,'o',label = "the discretize data $(x_i,f_i)$")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>

The one point finite difference formula¶

\(\dfrac{df}{dx} = \dfrac{f(x_{i+1})-f(x_i)}{x_{i+1}-x_i}\)

The numpy diff() function is a fast way to compute this formula:

df_1p = np.diff(fi)/np.diff(xi)

But becarefull the size of the resulting array is equal to \(n-1\).

fig = plt.figure()
plt.plot(x,df(x),'b',label = "Analytical derivative")
plt.plot(xi[:-1],df_1p,'ro',label = "Finite difference derivative (1 point)")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>

The 2 points finite difference formula¶

\(\dfrac{df}{dx} = \dfrac{f(x_{i+1})-f(x_{i-1})}{x_{i+1}-x_{i-1}}\)

Slicing method of numpy array is a good way to do such operation

# Slicing exemple
a=np.linspace(0,6,7)
print(a[2:])
print(a[:-2])

[ 2.  3.  4.  5.  6.]
[ 0.  1.  2.  3.  4.]

With slicing methode of numpy array the 2 points formula can be fast computed:

df_2p = (fi[2:] - fi[:-2])/(xi[2:] - xi[:-2])

fig = plt.figure()
plt.plot(x,df(x),'b',label = "Analytical derivative")
plt.plot(xi[:-1],df_1p,'ro',label = "Finite difference derivative (1 point)")
plt.plot(xi[1:-1],df_2p,'gs',label = "Finite difference derivative (2 points)")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>

Let’s have a look at the error on the derivative computation¶

fig = plt.figure()
plt.plot(xi[:-1],(df_1p-df(xi[:-1])),'ro',label = "Finite difference derivative (1 point)")
plt.plot(xi[1:-1],(df_2p-df(xi[1:-1])),'gs',label = "Finite difference derivative (2 points)")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>

What append if the data are noisy ?¶

fi_noise = f(xi) + 0.2*np.random.randn(fi.size)

fig = plt.figure()
plt.plot(x,f(x),':r',label = "$f(x)$")
plt.plot(xi,fi_noise,'o',label = "noisy data")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>

The 2 points fomrula is the more accurate. Let’s use it to derivate the noisy data:¶

df_2p_n = (fi_noise[2:] - fi_noise[:-2])/(xi[2:] - xi[:-2])

fig = plt.figure()
plt.plot(x,df(x),':r',label = "Analytical derivative")
plt.plot(xi[1:-1],df_2p_n,'-g',label = "Finite difference derivative of noisy data")
plt.legend()
plt.grid()

<IPython.core.display.Javascript object>