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)