Meta-Genetic Optimisation in multidimensional search spaces.

Introduction

The objective of this tutorial is to shed some light on the application of genetic optimisation to systems themselves represented as GA's (Genetic Algorithms). In the following text I will assume you have a basic/intermediate grasp on Genetic Programming. If not there're plenty of helpful resources online, of them my most recommended to beginners would be http://www.ai-junkie.com/ga/intro/gat1.html . A very short article oriented towards readers with little grounding in mathematics.

Concept

As with evolution in biological systems, we see three basic methods essential for the progression of any genetic system, recombination , mutation & fitness evaluation.

I should note that many article & writers will use 'crossover' instead of recombination, in this instance the terms are mostly interchangeable.

Each method and it's intended usage can be roughly defined as such:

Recombination: The production of a successor program where all instructions are sourced from two parents and a new method is constructed from ancestor code.

Mutation: The random shift of functions and/or routines within the program, designed to prevent settling on local minima.

Fitness: The fitness is a general term for a quantifiable rank designed to impose a bias for the intergenerational survival for certain traits and against others. The fitness function is in most cases the most essential element to heuristically direct the curvature of the search space towards the global minima.

Global & Local Maxima/Minima: Although minima and maxima have no explicit mention in the following algorithms they will play a crucial factor in the optimisation potential of each.

Firstly the global extremum, Global Minima and Global Maxima are simply the highest and lowest values within the search domain, whilst the local values refer to any sub-section within the domain separated by an inverse magnitude shift, for example: in an oscillation of constantly rising magnitude, the following peaks maximum will always exceed the pulse prior, the tip of each pulse can be said to be a local maximum whilst only the most recent would be the global maximum.

A point is the local maxima if all surrounding points have a lesser value; a point is the local minima if all surrounding points are of higher value.

The Task

The example task at hand will be a game from the sbs game show, Letters & Numbers and is played as followed.

There are two lists, one of high numbers (25, 50, 75, 100), one of low numbers (1-9).

The player must pick no more than 6 selections of a choice high or low, the numbers are then selected at random from each list, a random number is generated between 1-1000 and from the 6 numbers generated the player must then attempt to create an expression evaluating to the randomly generated number using each number only once, and only the operators [+, -, /, *, ()].

This task has a striking resemblance to the ai-junkie tutorial mentioned earlier. Each answer will consist of up to 12 “genes” or functions each gene is constructed from a string of bits. The genes are to be represented in pairs when the first gene is interpreted as a number, with the first So that the third bit in ever triplet represents whether it’s a high number 0 or low number 1, and the remaining two are interpreted as an integer representing array positioning. For simplicity I’ve restricted the total selections of numbers to a maximum of 4 in each category, each operator to start the parenthesis has been left out, thus only two bits are needed to store every state, however I will still use a triplet for each position to allow for future expansion,

As of such a program might look like this 001100010101001100010101 converted into an array of 1 third the integers it forms the state 1425 1425 which is converted to the value the [2nd small number, first large number, 3rd small number, 2nd large number, -,+,/,- keeping in mind a bit is dropped from the end of each triplet in the interpretation of the operators.

Having stated such, the evaluation is structured as following:

Two nested lists are included under the namespace done to keep a log of which numbers have been called and a comparison performed to ensure only unused numbers can be included in the operation.Code:def evalu(nbits,obits,numbers): num=list(numbers[:]) done=[[],[]] stri='' nu=[] o=[] operators=['+','-','/','*'] for i in range(int(len(nbits)/3)): if not len(num[0])+len(num[1]): break n=(lambda: (1 if nbits[i*3] else 0)+(2 if nbits[i*3+1] else 0))() if (nbits[i*3+2] or not len(num[0])) and len(num[1]): c=0 while num[1][(n+c)%len(num[1])] in done[1] and c < len(num[1]): c=c+1 e=num[1][(n+c)%len(num[1])] if not e in done[1]: nu.append(e) done[1].append(e) else: c=0 while num[0][(n+c)%len(num[0])] in done[0] and c < len(num[0]): c=c+1 e=num[0][(n+c)%len(num[0])] if not e in done[0]: nu.append(e) done[0].append(e) for i in range(int(len(obits)/3)): n=(lambda: (1 if obits[i*3] else 0)+(2 if obits[i*3+1] else 0))() o.append(operators[n]) for i in range(len(nu)-1): if len(o)==0: break stri=stri+str(nu[i])+o[i%len(o)] stri=stri+str(nu[len(nu)-1]) return stri

The initial state of any population is as simple as making an array of arrays of randomly generated Boolean True/False states.

My favoured way of doing this is

This would much more simply be structured as,Code:def randstr(n): return map(lambda x: round(random())==1,range(n))

however, for reasons unknown to me, lambda’s are parsed as an non-importable state.Code:randstr=lambda n:map(lambda x:round(random())==1,range(n))

We’re now ready to tackle the initialisation function,

Don’t let the compaction scare you, to break it down into English the function says:Code:def prodpop(self,size,l,n): return map(lambda x: (lambda y:(lambda z:[y,z,self.getfitness(z,self.goal,self.fitnesstype),0])(evalu(y[0:18],y[18:36],n)))(randstr(self.stringlength)),range(size))

create an array with the length of argument 1, initialise each position as an array as such,

[your randomly generated array of bits, the output of the ‘evalu’ function, the not yet introduced function ‘getfitness’ which is sent the output of the ‘evalu’ function]

and may be more familiar as:

which is a mere character longer in this case however in the case of applications which are not too cpu intensive, I’ve created habit of using maps where applicable.Code:def prodpop(self,size,l,n): out=[] for i in range(size): x=randstr(self.stringlength) y=evalu(x[0:18],[18,36],n) out.append([x,y,self.getfitness(y,self.goal,self.fitnesstype),0]) return out

The fitness function for this game is very simple, yet has the potential to easily become more complex,

We have four categories of ranking answers that represent positive bias and one for everything else, All correct answers are given a value of 10, with every answer only off by 10 or less on either side is given a value of 5 with 3 for a difference of less than 100 and 2 for a difference of less than a difference of 200. All other values are returned 1.

This is simply put,

The inclusion of more types will arise as the layers of metageneration are applied. Now that a generation can be populated and its fitness measured we have all the tools necessary to develop intergenerational development. The splicing will occur to mix a random selection of the fittest genes in the previous generation and a small chance of mutation applied by inversing any bits in the output.Code:def getfitness(self,val,goal,typ): if typ == 'step': dif=abs(eval(val)-goal) if dif == 0: if val not in self.answers: self.answers.append(str(val)) return 10 if dif < 10: return 5 elif dif < 100: return 3 elif dif < 200: return 2 else: return 1

The roulette selection and middle switch from the tutorial linked at the start will be built first, and appear here with the body of the iteration function

And the mutation as such:Code:def splice(self,a,b,typ): if typ=='middle': mid=int(min([len(a),len(b)])/2) return a[:mid]+b[mid:] def iterate(self): self.gen=self.gen+1 roulette=[] for i in self.population: for i2 in range(i[2]): roulette.append(i[0]) pop=map(lambda x: (lambda y:(lambda z:[y,z,self.getfitness(z,self.goal,self.fitnesstype),0])(evalu(y[0:18],y[18:36],self.numbers)))(self.mutate(self.splice( roulette[int(round(random()*(len(roulette)-1)))], roulette[int(round(random()*(len(roulette)-1)))],self.mergetype),self.muttype)) ,range(len(self.population))) self.population=pop

With each iteration, a roulette selection will be made, with a bias towards fitter expressions being their number of appearances hold proportion to fitness, two are randomly selected and merged to form a new function, this function then has a chance to change form and become a function with elements not previously contained within the prior form, if this mutation is an improvement the fitness test will ensure it succeeds.Code:def mutate(self,inp,typ,rate=0.02): out=[] if typ=='bit': for i in inp: out.append((lambda: i if random()<rate else not i)())

With the harder, possibly unfamiliar aspects explained let’s look at the whole code with added types:

Initialised with any input of sufficient mutation rate and population should correctly solve the puzzle whilst having to search only a tiny portion of the problem space. Picture for a second the problem space as an orientable surface and each linear input argument a perpendicular axis, for the means of simplicity let’s temporarily look at the task of a computer trying to draw a circle on the screen to a likeness of another given circle assume the centre position is already known, you have two inputs, r and col, each expression produces a fitness based on its likeness via a Boolean ‘And’ comparison, and the total difference in colour, if we then graph every possible outcome within a range of all colour and a radius of 1-50 if we graphed the difference of the two colours plus the difference in volumes of each circle on z, radius on y, and colours on x, we would end up with a gradual ridged slope looking like this, (z has been replaced with a greyscale).Code:from __future__ import division from time import time from random import random from math import floor def bintoint(inp): o=0 for i in range(len(inp)): o=o+(lambda: int((2**(i+1))/2) if inp[i] else 0)() return o def inttobin(inp): o=map(lambda x:x =='1',reversed(bin(inp)[2:])) return o def randstr(n): return map(lambda x: round(random())==1,range(n)) def nums(bibit): small=map(lambda x:x+1,range(9)) large=[25,50,75,100] out=[[],[]] n=2 for i in [0,1]: if bibit[i]:n=n+(i+1) else:None for i in range(n): r=int(round(random()*9.98-.49)) out[0].append(small[r%len(small)]) del small[r%len(small)] for i in range(6-n): r=int(round(random()*9.98-.49)) out[1].append(large[r%len(large)]) del large[r%len(large)] return out def evalu(nbits,obits,numbers): num=list(numbers[:]) done=[[],[]] stri='' nu=[] o=[] operators=['+','-','/','*'] for i in range(int(len(nbits)/3)): if not len(num[0])+len(num[1]): break n=(lambda: (1 if nbits[i*3] else 0)+(2 if nbits[i*3+1] else 0))() if (nbits[i*3+2] or not len(num[0])) and len(num[1]): c=0 while num[1][(n+c)%len(num[1])] in done[1] and c < len(num[1]): c=c+1 e=num[1][(n+c)%len(num[1])] if not e in done[1]: nu.append(e) done[1].append(e) else: c=0 while num[0][(n+c)%len(num[0])] in done[0] and c < len(num[0]): c=c+1 e=num[0][(n+c)%len(num[0])] if not e in done[0]: nu.append(e) done[0].append(e) for i in range(int(len(obits)/3)): n=(lambda: (1 if obits[i*3] else 0)+(2 if obits[i*3+1] else 0))() o.append(operators[n]) for i in range(len(nu)-1): if len(o)==0: break stri=stri+str(nu[i])+o[i%len(o)] stri=stri+str(nu[len(nu)-1]) return stri class gen(): def __init__(self,popsize,strlen,goal,numbers,mergetype,fitnesstest,muttype,mutationrate): self.mutationrate=mutationrate self.goal=goal self.numbers=numbers self.answers=[] self.fitnesstype=fitnesstest self.stringlength=strlen self.population=self.prodpop(popsize,strlen,self.numbers) self.mergetype=mergetype self.muttype=muttype self.time=time() self.gen=1 def splice(self,a,b,typ): if a > b :return a if b > a :return b if typ=='random':typ=['middle','tri','hex'][int(round(random()*2))] if typ=='middle': mid=int(min([len(a),len(b)])/2) return a[:mid]+b[mid:] if typ=='tri': return map(lambda x:a[x] if floor(x/3)%2 else b[x],range(len(a))) if typ=='hex': return map(lambda x:a[x] if floor(x/6)%2 else b[x],range(len(a))) def mutate(self,inp,typ,rate=0.02): out=[] if typ=='bit': for i in inp: out.append((lambda: i if random()<rate else not i)()) elif typ=='forward': for i in range(int(floor(len(inp)/3))): out=out+(lambda:inttobin((bintoint(inp[i*3:i*3+3])+1)%8) if random()<rate else inp[i*3:i*3+3])() elif typ=='tri': for i in range(int(floor(len(inp)/3))): out=out+(lambda:map(lambda x:not x,inp[i*3:i*3+3]) if random()<rate else inp[i*3:i*3+3])() elif typ=='hex': for i in range(int(floor(len(inp)/6))): out=out+(lambda:map(lambda x:not x,inp[i*6:i*6+6]) if random()<rate else inp[i*6:i*6+6])() return out def prodpop(self,size,l,n): out=map(lambda x: (lambda y:(lambda z:[y,z,self.getfitness(z,self.goal,self.fitnesstype),0])(evalu(y[0:18],y[12:36],n)))(randstr(self.stringlength)),range(size)) return out def iterate(self): self.gen=self.gen+1 roulette=[] for i in self.population: for i2 in range(i[2]): roulette.append(i[0]) pop=map(lambda x: (lambda y:(lambda z:[y,z,self.getfitness(z,self.goal,self.fitnesstype),0])(evalu(y[0:18],y[18:36],self.numbers)))(self.mutate(self.splice( roulette[int(round(random()*(len(roulette)-1)))], roulette[int(round(random()*(len(roulette)-1)))],self.mergetype),self.muttype)) ,range(len(self.population))) self.population=pop def getfitness(self,val,goal,typ): if typ == 'step': dif=abs(eval(val)-goal) if dif == 0: if val not in self.answers: self.answers.append(str(val)) return 10 if dif < 10: return 5 elif dif < 100: return 3 elif dif < 200: return 2 else: return 1 if typ == 'gradient': dif=abs(eval(val)-goal) if dif == 0: if val not in self.answers: self.answers.append(str(val)) return 10 elif dif < 200: return 10-int(round(dif/20))+2 else: return 1

Most, if not every problem in any Euclidian domain or dimension can be visualised as such, and the genetic algorithm is a route of seeking these peaks without iterating over every possible state, take the metaphor of flowing water, it always seeks to go downwards, the expressions will seek to maximise the fitness upon it’s z axis successfully travelling across the other two axis seeking the point of lowest resistance, it does this scanning by points between currently dominant points, as the search location to next evaluate this should work well as many problems exist as finite n-dimensional gradients or curves. Still working with the visualisation, picture a crater on the border on a larger crater, although ideally the water is attracted to pool at the lowest point, any water contained within the higher crater may lack the turbulence to escape past the boundaries of the crater, mutation can be thought of as this turbulence. Based on the distance between the global minima and any local minima around which the entire set is attracted, due to the magnitude of any local maxima separating viable paths of splicing, a function may be unable to evolve to the global minima, by increasing the mutation rate well allow for each new expression to differ more greatly from each ancestor whilst maintaining distribution around the previously most fit, increasing this factor we increase the turbulence of the flow allowing it to break from the confines of any craters and valleys that might not be the global minima. It's very important to find a comfortable range of mutation, too high and it might as well be a random selection, too low and you'll find everything that tries to journey outside your minima may go extinct.

By making the type selection, fitness test type, merge type, mutation rate and population size variable arguments in the class it’s merely the same method to create a bitstring evaluator for handling the inputs of our game solving GA,

This initialisation routine creates 'size' number of strings, all representing genetic algorithms of varied population, mutation rate, and splice/fitness/mutation methods allowing for the next layer of genetic operation to be applied to the initial game routine which will look something like this:Code:def prodpop(self,size,l,n): pop=map(lambda x:randstr(l),range(size)) population=[100,250,350,500,750,1000,1250,1500] merge=['middle','tri','random','hex'] fitness=['step','gradient'] mut=['bit','forward','tri','hex'] out=map(lambda x:[x,gen(population[bintoint(x[0:3])],36,self.goal,n,merge[bintoint(x[3:5])],fitness[bintoint([x[5]])],mut[bintoint(x[6:8])],bintoint(x[8:13])/1.8/100)],pop) return out

With the constructed fitness function there is no bias imposed by the time it takes to compute, this should create an obvious slope towards 111 for the first gene, and should be the next thing to be updated, this search function will navigate the space of directing previous generations how to evolve in order to create a search space of search spaces in order to conclude around a median of the global minima.Code:from random import random from gen import gen, nums, randstr from time import time def bintoint(inp): o=0 for i in range(len(inp) ): o=o+(lambda: (2**(i+1))/2 if inp[i] else 0)() return o class mgen(): def __init__(self,popsize,strlen,mutationrate): self.goal=333#int(round(random()*1000)) self.n=nums([1,0]) self.stringlength=strlen self.mutationrate=mutationrate self.population=self.prodpop(popsize,strlen,self.n) for i in range(maxgen*2): for i2 in self.population: i2[1].iterate() def splice(self,a,b): mid=int(min([len(a),len(b)])/2) return a[:mid]+b[mid:] def mutate(self,inp,rate=0.02): out=[] for i in inp: out.append((lambda: i if random()<rate else not i)()) return out def prodpop(self,size,l,n): pop=map(lambda x:randstr(l),range(size)) population=[100,250,350,500,750,1000,1250,1500] merge=['middle','tri','random','hex'] fitness=['step','gradient'] mut=['bit','forward','tri','hex'] out=map(lambda x:[x,gen(1250,36,self.goal,n,merge[bintoint(x[3:5])],fitness[bintoint([x[5]])],mut[bintoint(x[6:8])],bintoint(x[8:13])/1.8/100)],pop) return out def iterate(self): population=[100,250,350,500,750,1000,1250,1500] pop=[] merge=['middle','tri','random','hex'] fitness=['step','gradient'] mut=['bit','forward','tri','hex'] for i in range(maxgen): for i2 in self.population: i2[1].iterate() roulette=[] for i in self.population: for i2 in range(len(i[1].answers)+1): roulette.append(i[0]) for i in range(len(self.population)): pop.append(self.mutate(self.splice( roulette[int(round(random()*(len(roulette)-2)))], roulette[int(round(random()*(len(roulette)-2)))]))) pop=map(lambda x:[x,gen(population[bintoint(x[0:3])], 36, self.goal, self.n, merge[bintoint(x[3:5])], fitness[bintoint([x[5]])], mut[bintoint(x[6:8])],bintoint(x[8:13])/1.8/100)],pop) self.goal=int(round(random()*1000)) self.n=nums([1,0]) self.population=pop for i in range(maxgen): for i2 in self.population: i2[1].iterate() l=14 w=map(lambda x: 0,range(l)) maxgen=12 b=mgen(80,l,0.03) for i in range(200): print '.' print str(b.goal)+' '+str(b.n) b.iterate() highest=[0,0] for i2 in range(len(b.population)): e=b.population[i2][1].answers if e !=[]: if len(e)>highest[1]: highest=[e,len(e)] print highest for e in b.population: for i2 in range(len(e[0])): w[i2]=w[i2]+(lambda : 1 if e[0][i2] else -1)() print w

Tweet This+ 1 thisPost To Linkedin