In a previous discussion on optimization, we touched on the importance of optimization within complex systems. To tackle such these non-intuitive problems, we were introduced to man-made and nature-inspired techniques. Topping the list of nature-inspired algorithms was Particle Swarm Optimization (PSO) — easily one of my favorite optimization routines and one that we can get up and running quickly.

In this tutorial, we’ll do just that — start with a blank script and build a PSO routine from scratch using the Julia programming language. Before we start, it’s worth noting that you can easily find a pre-built PSO package and implement it into your algorithm pipeline. Personally, though, I find that building my own code helps me understand the algorithm on a deeper level. With that, let’s get started!

Implementation

Every analysis algorithm is based on some sort of dataset and that’s where we’ll start as well. Since the goal of our optimization routine is to find a global maximum, we’ll want to build a dataset that has both local and global extrema. For our purposes, a combination of sines and cosines will work nicely within a specified window. Let’s begin by building our “toy” data.

For illustration purposes, we’ll define a two-variable function:

[math] f(x,y)=cos(x/5) * sin(y/3) [/math]

In reality, though, many systems will be much more complex and contain more than two variables. In fact, this is one of the strongest traits of optimization routines like PSO — instead of sprinkling particles over two dimensions (ex. [math] f(x,y) [/math]), we would sprinkle them over the N-dimensional space (ex [math] f(x,y,z,N,q,b,l,t) [/math]) of our more complex problem!

import numpy as np
import matplotlib.pyplot as plt

# Plot data
fig = plt.figure()
ax = plt.axes(projection='3d')

# Define 3D toy data
x = np.arange(start=-1,stop=10,step=0.1)
y = np.arange(start=-1,stop=10,step=0.1)
X, Y = np.meshgrid(x, y)
Z = np.sin(X/2)*np.cos(Y/2)

surf = ax.plot_surface(X, Y, Z, cmap = plt.cm.cividis)
ax.set_xlabel('x', labelpad=20)
ax.set_ylabel('y', labelpad=20)
ax.set_zlabel('z', labelpad=20)

fig.colorbar(surf, shrink=0.5, aspect=8)

plt.show()

Generate particles

With our three-dimensional space defined, we can start sprinkling in our particles. We begin by choosing random (x,y) positions within the bounds our initial function space. For each particle position (xpos,ypos), we find the value at that coordinate (zpos). We then store each of these arrays in our current dataframe (df) for bookkeeping.

(in progress)