Thursday, March 14, 2019

SciPy - 7 (Interpolation)

One often has a number of data points, obtained by sampling or experimentation, which represent the values of a function for a limited number of values of the independent variable. It is often required to interpolate, i.e., estimate the value of that function for an intermediate value of the independent variable. Interpolation is a method of constructing new data points within the range of a discrete set of known data points.

Interpolation is simple and convenient in scipy: The interp1d function, when given arrays describing X and Y data, returns an object that behaves like a function that can be called for an arbitrary value of x (in the range covered by X), and it returns the corresponding interpolated y value. See the following program:

from scipy.interpolate import *
from numpy import *
import numpy as N
from matplotlib import pyplot as plt
from scipy import arange , cos , exp
import pylab
import math

def f(x):
return sin(x)

n = arange(0, 10)
x = linspace(0, 9, 100)
y_meas = f(n) + 0.1 * N.random.randn(len(n)) # simulate measurement with noise
y_real = f(x)
linear_interpolation = interp1d(n, y_meas)
y_interp1 = linear_interpolation(x)
cubic_interpolation = interp1d(n, y_meas, kind='cubic')
y_interp2 = cubic_interpolation(x)

fig, ax = plt.subplots(figsize=(10,4))
ax.plot(n, y_meas, 'bs', label='noisy data')
ax.plot(x, y_real, 'k', lw=2, label='true function')
ax.plot(x, y_interp1, 'r', label='linear interp')
ax.plot(x, y_interp2, 'g', label='cubic interp')
ax.legend(loc=3)
plt.show()

The output shows the interpolation-



For a given a set of N points (xi; yi) with i = 1; 2; : : :N, we sometimes need a function which returns yi = f(xi) where x == xi, and which in addition provides some interpolation of the data (xi; yi) for all x.
The function y0 = scipy.interpolate.interp1d(x,y,kind='nearest') does this interpolation based on splines of varying order. Note that the function interp1d returns a function y0 which will then interpolate the x-y data for any given x when called as y0(x). The program below demonstrates this:

from scipy.interpolate import *
from numpy import *
from matplotlib import pyplot as plt
from scipy import arange , cos , exp
import pylab
import math
import numpy as np
import scipy . interpolate

def create_data (n):
    """ Given an integer n, returns n data points
    x and values y as a numpy . array ."""
    xmax = 5.
    x = np. linspace (0, xmax , n)
    y = - x**2
    # make x- data somewhat irregular
    y += 1.5 * np. random . normal ( size =len(x))
    return x, y
   
   
# main program
n = 10
x, y = create_data (n)
#use finer and regular mesh for plot
xfine = np. linspace (0.1 , 4.9 , n * 100)
# interpolate with piecewise constant function (p=0)
y0 = scipy . interpolate . interp1d (x, y, kind ='nearest')
# interpolate with piecewise linear func (p=1)
y1 = scipy . interpolate . interp1d (x, y, kind ='linear')
# interpolate with piecewise constant func (p=2)
y2 = scipy . interpolate . interp1d (x, y, kind ='quadratic')
pylab . plot (x, y, 'o', label ='data point ')
pylab . plot (xfine , y0( xfine ), label ='nearest ')
pylab . plot (xfine , y1( xfine ), label ='linear ')
pylab . plot (xfine , y2( xfine ), label ='cubic ')
pylab . legend ()
pylab . xlabel ('x')

pylab . show ()

The output below shows the different interpolation kinds:




Here I am ending this post regarding Interpolation with SciPy. Next post will discuss about Statistics. Till we meet next keep learning and practicing Python as Python is easy to learn!  
Share:

0 comments:

Post a Comment