SciPy, a Python library for mathematics, science and engineering

SciPy is a library in Python, or more precisely a collection of mathematical algorithms and other functions particularly used in science. SciPy is built on NumPy, a library that extends the Python language to better manage mathematical calculations.

Thanks to SciPy, a Python session has nothing to envy to other environments like Matlab, IDL, Octave and SciLab. SciPy will prove to be an irreplaceable tool for anyone dealing with the world of scientific calculations.

SciPy, a Python library for mathematics, science and engineering

How the SciPy library is organized

SciPy is organized in different sections, each of which covers different topics of scientific calculation. The structure of the library is therefore composed of different modules (sub-packages), each inherent to a specific application topic

  • cluster – Clustering algorithms
  • constants – Mathematical and physical constants
  • fftpack – Fourier Transform functions
  • integrate – Integration and ordinary differential equations
  • interpolate – Interpolation and smoothing spline
  • io – Input and Output
  • linalg – Linear algebra
  • ndimage – N-dimensional image processing
  • odr – Regression of orthogonal distances
  • optimize – Optimization and root-finding functions
  • signal – Signal Processing
  • sparse – Sparse matrices and associated functions
  • spatial – Spatial data structures and algorithms
  • special – Other special features
  • stats – Statistical distributions and related functions

Generally, when working with SciPy, only the necessary modules are imported individually. For example:

>>> from scipy import some_module
>>> some_module.some_function() 

Before starting

SciPy is a library built on NumPy, so every time we need to work with arrays we need to import it first.

>>> import numpy as np

Moreover, when working with scientific data it is often necessary to use graphical representations. So another very important (and useful) library is matplotlib.

>>> import matplotlib as mpl
>>> import matplotlib.pyplot as plt

Online documentation

Both SciPy and NumPy and Matplotlib have a lot of online documentation. As for the official documentation, many HTML and PDF files are available at docs.scipy.org.

If instead you want to go into more detail on the subject and have many examples and tutorials on SciPy, I recommend this book. It is not really very recent, but it is very exhaustive on the subject.

Elegant SciPy book - O'Reilly

Elegant SciPy: The Art of Scientific Python

J.Nunez-Iglesias et al. – O’ Reilly 2017

As for NumPy, I recommend this book. There is an introductory chapter to this very simple and intuitive library, full of examples. In a few hours, you will become familiar with this library by understanding all the basic features. Furthermore there are other chapters on the use of Matplotlib applied to scientific libraries and practical applications.

Python Data Analytics (2nd Edition)

Python Data Analytics: With Pandas, NumPy, and Matplotlib (2nd Edition)

Fabio Nelli – Apress 2018

Example – the Voronoi Diagram

It is impossible in a single article to speak entirely of all the features of the SciPy library, as well as going into detail for each individual module.

A good introduction to SciPy may be to choose a sample application, such as the Voronoi diagrams. With a few lines of code, you can guess the similarities with Matlab, SciLab and Octave, and how it is possible to work interactively with a mathematical or scientific problem. Great merit however is that SciPy is completely integrated with Python, with the possibility of integrating it with any other library or application.

Before starting with the code, however, we see a brief introduction to the Voronoi diagrams, also known as Dirichlet tessellation. The operation consists in dividing an area into n convex polygons, determined by n points defined in this area. The polygons must be generated in such a way that every point within this polygon is closer to its point of generation (one of the n points) than to any other.

Although they started with the introduction of computer science, the Voronoi diagrams were first considered in 1644 by René Descartes. Dirichlet used them in 1850 for the study of positive quadratic forms. Finally they were definitively studied and developed by Voronoi in 1907, extending them even to cases with higher dimensions. Voronoi diagrams find many applications in various fields such as computer graphics, epidemiology, geophysics and meteorology.

Starting to program

To develop the following example I recommend using the Jupyter notebook or opening an IPython session, so as to work interactively, entering and evaluating the snippet code for snippets.

Note: To gradually see the evolution of the example, add the command plt.show() at the end and then delete it and continue on.

First we import the necessary SciPy libraries and modules.

import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial import KDTree, Voronoi

There are two approaches to the Voronoi problem. One that can be solved with the KDTree module with which we will determine the areas of the diagram. Each area will be assigned a different color, and the other with the Voronoi module, with which we will calculate the separation segments between the various areas.

Now insert the 9 points between a rectangular area of size 7×8

points = np.array([[0, 0], [4, 4], [5, 6.5], 
                   [4.2,5.8], [1.7, 1.4], [3.1, 7.6], 
                   [3, 6], [3.4, 5.2], [1.8, 2.2], 
                   [2.7, 0.4], [6.5, 1.8], [2.3, 2.8]])

Let’s now pass to the two objects corresponding to the two methods: KDTree () and Voronoi ().

tree = KDTree(points)
vor = Voronoi(points)

We create a grid by dividing the area in question into 40000 particles, (200×200). The number is a compromise between quality and speed of calculation. If you want you can change these values to your liking. Divide the area into smaller particles. For convenience we will extend the area of representation a (-1, -1) in order to represent the point of origin (0,0) as well.

x = np.linspace(-1, 7, 200)
y = np.linspace(-1, 8, 200)
xx, yy = np.meshgrid(x, y)
xy = np.c_[xx.ravel(), yy.ravel()]

Now all we have to do is represent everything using the matplotlib functions. All the values necessary for the representation have already been calculated, once the tree and vor objects have been generated.

First we insert the 9 points we entered into the representation. These will be represented as black circles. Then we delimit the area of representation with the X axis going from -1 to 7 and the Y axis going from -1 to 8.

plt.plot(points[:, 0], points[:, 1], 'ko')
plt.xlim(-1, 7); plt.ylim(-1, 8)

Now we will represent the segments that separate the various areas of the Voronoi diagram, joining the Voronoi vertices calculated through the object vor.

for simplex in vor.ridge_vertices:
    simplex = np.asarray(simplex)
    if np.all(simplex >= 0):
       plt.plot(vor.vertices[simplex, 0], vor.vertices[simplex, 1], 'k-')

What is much more difficult to represent are the segments that start at the vertices of Voronoi and go towards infinity. Here’s how to represent them.

center = points.mean(axis=0)
for pointidx, simplex in zip(vor.ridge_points, vor.ridge_vertices):
     simplex = np.asarray(simplex)
     if np.any(simplex < 0):
         i = simplex[simplex >= 0][0] # finite end Voronoi vertex
         t = points[pointidx[1]] - points[pointidx[0]]  # tangent
         t = t / np.linalg.norm(t)
         n = np.array([-t[1], t[0]]) # normal
         midpoint = points[pointidx].mean(axis=0)
         far_point = vor.vertices[i] + np.sign(np.dot(midpoint - center, n)) * n * 100
         plt.plot([vor.vertices[i,0], far_point[0]],
              [vor.vertices[i,1], far_point[1]], 'k-')

Finally we represent the various areas of the Vorodoi diagram with different colors, using the tree object.

plt.pcolor(x, y, tree.query(xy)[1].reshape(200, 200))

And finally we just have to send the command to represent the diagram.

plt.show()

Here is the Voronoi diagram corresponding to the 9 points entered

Leave a Reply