As many (or very few in the real life haha) people know, today is the Pi Approximation Day ! So it’s time to make a contribution to celebrate this funny day =)
My contribution is to use Python and Pyevolve to approximate Pi number using Genetic Programming approach. I’ve created the functions gp_add(+), gp_sub(-), gp_div(/), gp_mul(*) and gp_sqrt (square root) to use as non-terminals of the GP. The fitness function is very simple too, it simple returns the absolute difference between the Python math.pi and the evaluated individual. I’ve used also a population size of 1k individuals with max tree depth of 8 and the random ephemeral constants as random integers. The best approximation I’ve got while running the GP for about 8 minutes (40 generations) was 3.1416185511, best for 3 digits, you can improve it and let it run for more time to get better approximations.
Here is the formulae I’ve got with the GP (click to enlarge):
And here is the output of the script:
Best (0): 3.1577998365
Error: 0.0162071829
Best (10): 3.1417973679
Error: 0.0002047143
Best (20): 3.1417973679
Error: 0.0002047143
Best (30): 3.1417973679
Error: 0.0002047143
Best (40): 3.1416185511
Error: 0.0000258975
- GenomeBase
Score: 0.000026
Fitness: 15751.020831
Params: {'max_depth': 8, 'method': 'ramped'}
Slot [Evaluator] (Count: 1)
Slot [Initializator] (Count: 1)
Name: GTreeGPInitializator - Weight: 0.50
Doc: This initializator accepts the follow parameters:
*max_depth*
The max depth of the tree
*method*
The method, accepts "grow" or "full"
.. versionadded:: 0.6
The *GTreeGPInitializator* function.
Slot [Mutator] (Count: 1)
Name: GTreeGPMutatorSubtree - Weight: 0.50
Doc: The mutator of GTreeGP, Subtree Mutator
.. versionadded:: 0.6
The *GTreeGPMutatorSubtree* function
Slot [Crossover] (Count: 1)
Name: GTreeGPCrossoverSinglePoint - Weight: 0.50
- GTree
Height: 8
Nodes: 21
GTreeNodeBase [Childs=1] - [gp_sqrt]
GTreeNodeBase [Childs=2] - [gp_div]
GTreeNodeBase [Childs=2] - [gp_add]
GTreeNodeBase [Childs=0] - [26]
GTreeNodeBase [Childs=2] - [gp_div]
GTreeNodeBase [Childs=2] - [gp_mul]
GTreeNodeBase [Childs=2] - [gp_add]
GTreeNodeBase [Childs=2] - [gp_sub]
GTreeNodeBase [Childs=0] - [34]
GTreeNodeBase [Childs=2] - [gp_sub]
GTreeNodeBase [Childs=0] - [44]
GTreeNodeBase [Childs=0] - [1]
GTreeNodeBase [Childs=2] - [gp_mul]
GTreeNodeBase [Childs=0] - [49]
GTreeNodeBase [Childs=0] - [43]
GTreeNodeBase [Childs=1] - [gp_sqrt]
GTreeNodeBase [Childs=0] - [18]
GTreeNodeBase [Childs=0] - [16]
GTreeNodeBase [Childs=2] - [gp_add]
GTreeNodeBase [Childs=0] - [24]
GTreeNodeBase [Childs=0] - [35]
- GTreeGP
Expression: gp_sqrt(gp_div(gp_add(26,
gp_div(gp_mul(gp_add(gp_sub(34,
gp_sub(44, 1)), gp_mul(49, 43)), gp_sqrt(18)),
16)), gp_add(24, 35)))
And finally, here is the source code:
from __future__ import division
from pyevolve import *
import math
def gp_add(a, b): return a+b
def gp_sub(a, b): return a-b
def gp_div(a, b): return 1 if b==0 else a/b
def gp_mul(a, b): return a*b
def gp_sqrt(a): return math.sqrt(abs(a))
def eval_func(chromosome):
code_comp = chromosome.getCompiledCode()
ret = eval(code_comp)
return abs(math.pi - ret)
def step_callback(engine):
gen = engine.getCurrentGeneration()
if gen % 10 == 0:
best = engine.bestIndividual()
best_pi = eval(best.getCompiledCode())
print "Best (%d): %.10f" % (gen, best_pi)
print "\tError: %.10f" % (abs(math.pi - best_pi))
return False
def main_run():
genome = GTree.GTreeGP()
genome.setParams(max_depth=8, method="ramped")
genome.evaluator += eval_func
ga = GSimpleGA.GSimpleGA(genome)
ga.setParams(gp_terminals = ['ephemeral:random.randint(1, 50)'],
gp_function_prefix = "gp")
ga.setMinimax(Consts.minimaxType["minimize"])
ga.setGenerations(50000)
ga.setCrossoverRate(1.0)
ga.setMutationRate(0.09)
ga.setPopulationSize(1000)
ga.stepCallback.set(step_callback)
ga.evolve()
best = ga.bestIndividual()
best.writeDotImage("tree_pi.png")
print best
if __name__ == "__main__":
main_run()
If you are interested why today is the Pi Approximation day, see some resources:
Little Cartoon
Some Background History
Some Pi Approximations