Verifying Genetic Algorithm on Simple Functions¶
- Goal: find the minimum of $f(\boldsymbol{x}) = \boldsymbol{x}_1^2 + \boldsymbol{x}_2^2 + \boldsymbol{x}_3^2 + \boldsymbol{x}_4^2 + \boldsymbol{x}_5^2$
More hints on Jupyter notebook/programming norms¶
- Typesetting math equations: \$f(x) = x_1^2 + x_2^2 + x_3^2 + x_4^2 + x_5^2\\$
- Rules to name variables/functions: camel case, snake case, PEP 8
Naming conventions according to PEP 8:
| Type | Naming Convention | Examples |
|---|---|---|
| Function | Use a lowercase word or words. Separate words by underscores to improve readability. | function, my_function |
| Variable | Use a lowercase single letter, word, or words. Separate words with underscores to improve readability. | x, var, my_variable |
| Class | Start each word with a capital letter. Do not separate words with underscores. This style is called camel case. | Model, MyClass |
| Method | Use a lowercase word or words. Separate words with underscores to improve readability. | class_method, method |
| Constant | Use an uppercase single letter, word, or words. Separate words with underscores to improve readability. | CONSTANT, MY_CONSTANT, MY_LONG_CONSTANT |
| Module | Use a short, lowercase word or words. Separate words with underscores to improve readability. | module.py, my_module.py |
| Package | Use a short, lowercase word or words. Do not separate words with underscores. | package, mypackage |
In [1]:
import numpy as np
import matplotlib.pyplot as plt
from itertools import product
# make sure the result is reproducible
np.random.seed(0x1a2b3c4d)
In [2]:
INIT_SIZE = 50
POP_SIZE = 100
CROSS_RATE = 0.95
MUTATION_RATE = 0.005
N_CROSS = 500
N_GEN = 200
N_DIM = 5
Generating initial population¶
An important trick of Python programming is to avoid using loops (because they are slow!). Check out how I generate the initial population without explicit loops.
- We want to avoid negative numbers as they are encoded by two's complement, but we also want to include negative numbers in our genetic algorithm. What do we do? We limit the range of a variable $x$ between 0 and 255, and then when we put $x - 127$ into the equation instead. In this way, we are able to generate both positive and negative integers smoothly.
- Notice how I specify the data type with
dtype=np.uint8. This will be useful afterwards.
In [3]:
init_pop = np.random.randint(0, 100, size=(INIT_SIZE, N_DIM), dtype=np.uint8)
print('shape:', init_pop.shape)
print('excerpt:\n', init_pop[:10, :])
Generating chromosome¶
I am using np.uint8 because it corresponds to an unsigned 8-bit integer. We can use np.unpackbits to convert such arrays to bits.
In [4]:
test_uint8_arr = np.array([0, 1, 2], dtype=np.uint8)
test_uint8_arr_bits = np.unpackbits(test_uint8_arr)
print(test_uint8_arr_bits)
# use np.packbits to reverse the operation
test_uint8_arr_repack = np.packbits(test_uint8_arr_bits)
print(test_uint8_arr_repack)
In [5]:
# chromosome can be acquired rather easily
init_pop_ch = np.unpackbits(init_pop, axis=1)
print('shape:', init_pop_ch.shape)
print('excerpt:\n', init_pop_ch[:10, :])
Computing fitness score¶
We can let the fitness score be $\frac{1}{f(\boldsymbol{x}) + 1}$. When $f(\boldsymbol{x})$ increases, the fitness score decreases, which is exactly what we want.
In [6]:
def fitness_score(bitarr):
# get numerical array, convert to np.int for calculation
arr = np.packbits(bitarr).astype(np.int)
# allow both positive and negative values
arr -= 127
# computing f(x) without loops!
fx = np.sum(np.square(arr))
return 1.0 / (fx + 1.0)
Selection¶
In [7]:
def selection(pop_ch):
#if pop_ch.shape[0] < POP_SIZE: return pop_ch
# apply fitness score function on each sample
fitnesses = np.apply_along_axis(fitness_score, axis=1, arr=pop_ch)
total_fitness = np.sum(fitnesses)
# compute fitness proportions
fitness_prop = fitnesses / total_fitness
#print(fitness_prop)
# generate partial sums with np.cumsum
# no need to compute m repetitively
prop_cumsum = np.cumsum(fitness_prop)
#print(prop_cumsum)
samples = []
# begin sampling
for _ in range(POP_SIZE):
t = np.random.uniform(0.0, 1.0)
for i in range(prop_cumsum.size):
if t < prop_cumsum[i]:
samples.append(pop_ch[i, :])
break
return np.asarray(samples)
Crossover¶
In [8]:
def crossover(pop_ch):
# partition the population into two groups
males, females = np.array_split(pop_ch, 2, axis=0)
# get indices of each group
male_ind = np.arange(males.shape[0])
female_ind = np.arange(females.shape[0])
assert male_ind.size * female_ind.size >= POP_SIZE
# keep track of combinations that are used
used_ind = set()
result = []
# here we assume the size of Cartesian product is far greater than
# our desired sample size. Therefore, the probability of used index
# collision is very low
for _ in range(POP_SIZE):
# generate an unused combination
comb = None
while True:
x = np.random.randint(0, male_ind.size - 1)
y = np.random.randint(0, female_ind.size - 1)
if (x, y) not in used_ind:
comb = (x, y)
used_ind.add(comb)
break
male = males[comb[0], :]
female = females[comb[1], :]
t = np.random.uniform(0.0, 1.0)
if t < CROSS_RATE:
# partition male and female sample (at the middle)
male1, male2 = np.array_split(male, 2)
female1, female2 = np.array_split(female, 2)
# crossover two samples
new1 = np.concatenate([male1, female2])
new2 = np.concatenate([female1, male2])
result.append(new1)
result.append(new2)
else:
result.append(male)
result.append(female)
return np.asarray(result)
Mutation¶
In [9]:
def flip_bit(bit):
return 1 if bit == 0 else 0
def mutation(pop_ch):
# get a copy of pop_ch to avoid modifying the original array
result = pop_ch.copy()
# get total number of bits in the population
total_bits = pop_ch.size
# generate a random number for each bit
random_probs = np.random.uniform(0.0, 1.0, size=total_bits)
# find out indices that are less than mutation probability
mutation_indices = np.where(random_probs < MUTATION_RATE)[0]
# convert flat index to 2D index
pop_indices = np.unravel_index(mutation_indices, pop_ch.shape)
# flip the bits
for i in range(mutation_indices.size):
idx = (pop_indices[0][i], pop_indices[1][i])
#print(idx)
result[idx] = flip_bit(result[idx])
return result
Putting things together¶
In [10]:
# a function to generate statistics of each generation
def get_gen_stat(pop_ch):
fitnesses = np.apply_along_axis(fitness_score, axis=1, arr=pop_ch)
return fitnesses
In [11]:
pop_ch = init_pop_ch
stats = []
for _ in range(N_GEN):
pop_ch = selection(pop_ch)
pop_ch = crossover(pop_ch)
pop_ch = mutation(pop_ch)
stats.append(get_gen_stat(pop_ch))
# add one final selection to limit the output size
pop_ch = selection(pop_ch)
In [12]:
final_pop = np.packbits(pop_ch, axis=1).astype(np.int)
print(final_pop - 127)
In [13]:
# graph generation stats
max_fitnesses = [np.max(arr) for arr in stats]
min_fitnesses = [np.min(arr) for arr in stats]
avg_fitnesses = [np.mean(arr) for arr in stats]
xs = np.arange(len(max_fitnesses))
plt.plot(xs, max_fitnesses, 'b', label='max fitness')
plt.plot(xs, min_fitnesses, 'r', label='min fitness')
plt.plot(xs, avg_fitnesses, 'g', label='average fitness')
plt.legend()
plt.xlabel('generation')
plt.ylabel('fitness')
plt.show()
In [ ]: