//************************************************************** // File:GA_v15.sce // Authors: Gerd Doeben-Henisch // Version Start: January-18, 2010 // Version Last: May-23, 2011 //****************************************************************** // Function: Implements a simple GA according to Goldberg (1989) // //***************************************************************** // Apri-l5, 2011 : Adding a mutation operator // //*************************************************** // April-8, 2011 : Adding an automatic generation of a population with random values // //**************************************************** // April-21, 2011: Adding a frequency counter for binary strings // //************************************************** // // April-22-27, 2011: Adding an evaluation of the quality of sets of solutions // //************************************************** // // April-28, 2011: Computing the probability of certain patterns // //************************************************* // // April-29, 2011: GA without fitness // //************************************************* // // May-01, 2011: Virtual fitness function for comparison // //************************************************* // // May-05, 2011: Conditioned Mutation // //************************************************* // // May-12, 2011: Fitness ratios, convergence times // //************************************************* // // May-18, 2011: Purified version of ga() // and : more converting functions vec22dec() // //************************************************* // // May-22, 2011: Improving usage of parameter p in case // of crossoverPrep with helpfunctions: the flag is fixed to 'l+7' and // freed from p // Corrections to crossoverPrep(); mismatch with while - conditions // //************************************************* // // May-23, 2011: Added a counter for the number of genomes in a POP which have // reached 100\% // //************************************************* // Helpful scilab-functions from the scilab libraries // // nancumsum ? This function returns the cumulative sum of the values of a matrix // tabul ? frequency of values of a matrix or vector // variance ? variance of the values of a vector or matrix // strange ? range // stdevf ? standard deviation // bin2dec := translate a binary string of '1','0' into a decimal number // // SW-STRUCTURE // // All data are organized in a dynamical table. Left hand the actual genes, in the middle // supporting parameters and to the right intermediate modifications of the genes to the left. // // POP FORMAT // l := length of strings // p := number of cells between string <1...l> and <l+p+1, ..., l+p+l> // Pos 1-l := String // Pos l+1 := Decimal Value of binary string // Pos l+2 := Fitness for l+1 and l+p // Pos l+3 := Percentage of overall fitness // Pos l+4 := Expected count according to fitness // Pos l+5 := Realized count // Pos l+6 := 2nd Decimal value of a compund fitness // Pos l+7 := crossover flag // p := 7 // Variable for overall Fitness FITNESS_ALL = 0 // Variable for average Fitness AFITNESS = 0 //*************************************************** //Function to generate automatically a population with random values // Input: // l := length of strings // n := size of population // Output: // A population according to the above defined format // function[POPX]=popgen(n,l) //Generate a matrix with only '0's POPX = zeros(n,l) MaxCells = n*l // Distribute randomly '1's for i=1:n for j=1:l if( rand() < 0.5 ) then POPX(i,j) = 0, else POPX(i,j) = 1, end end end endfunction //*************************************************** //Function to generate random numbers between [1,up] // // n := number of choices function[Xs]=randInt(n, up) Xs = [] for i=1:n x=0 while(x==0) x = round(up*rand()) end Xs(i) = x end endfunction //*************************************************** //Function to compute descriptive parameters of a population // Input: // l := length of strings // n := size of population // Output: // Pc := Percentage of possible space (PSpace) // Pb := Probability of gitting a certain individual of PSpace function[Pc, Pb]=popDescr(n,l) // Cardinality of PSpace PSPC = 2^l // Percentage Pc Pc = n/(PSPC/100) // Probability Pb Pb = n/2^l endfunction //*************************************************** //Function to compute the frequencies of the elements of a population // // Input: // l := length of strings // n := size of population // N := number of events // show := if '1' then display the POPX-matrix // // Output: // FX := Vector of freuencies ordered along the decimal values of a string (1st column) // FXDif := Vector of Differences compared to theoretical expectation in % (2nd column) // MEAN := The mean value of the differences // STD := The standard deviation of the differences // // Idea: Generate a sequence of populations, count the occurences, and check the degree of agreement between // empirical values and theoretica lexpectations function[DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,FX, POPX]=frequCheck(l,n,N,show) // Install the counters r= (2^l) c= 3 FX = zeros(r,c) // Generate an internal index for display for j=1:r FX(j,1)=j-1 end // Generate the populations for i=1:N [POPX]=popgen(n,l) if show==1 then printf("POPX first\n"), disp(POPX), end //Compute decimal value ('0' is placed at index '1' !!!) // Count occurences for j=1:n v=POPX(j,:), d = vec2dec(v,l), [POPX(j,l+1)]= d, FX(d+1,2)=FX(d+1,2)+1 end//n if show==1 then printf("POPX second\n"), disp(POPX), end end //N for j=1:r FX(j,3)=(FX(j,2)/n)/N //Normal frequency FX(j,4)=(FX(j,2)/n)/(((1/2^l)*N)/100) // Frequency as percentage end//m disp(FX) FREQ1=tabul(FX(:,3)) MEAN1=mean(FREQ1(:,1)) STD1= stdevf(FREQ1(:,1), FREQ1(:,2)) FREQ=tabul(FX(:,4)) MEAN=mean(FREQ(:,1)) STD= stdevf(FREQ(:,1), FREQ(:,2)) MAX=max(FREQ(:,1)) MIN=min(FREQ(:,1)) DIST=MAX-MIN DIST2=DIST/2 printf("Number of Events n * N = %d\n",n*N) // Show Plot // clf(), xdel, // plot2d(FREQ(:,1), FREQ(:,2)) endfunction //*************************************************** //Function to compute the frequencies of the elements of a population // and using a fitness function to compute 'virtual' fitness to compare // the populations with others by fitness. But the function does not make use of these values. // // Input: // l := length of strings // n := size of population // N := number of events // show := if '1' then display the POPX-matrix // // Output: // FX := Vector of freuencies ordered along the decimal values of a string (1st column) // FXDif := Vector of Differences compared to theoretical expectation in % (2nd column) // MEAN := The mean value of the differences // STD := The standard deviation of the differences // // Idea: Generate a sequence of populations, count the occurences, and check the degree of agreement between // empirical values and theoretica lexpectations function[FITNESS_ALL_PERC1,DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,FX, POPX]=frequCheck1(l,n,N,show) // Install the counters r= (2^l) c= 3 FX = zeros(r,c) FITNESS_ALLLOG1=[] // Generate an internal index for display for j=1:r FX(j,1)=j-1 end // Generate the populations for i=1:N [POPX]=popgen(n,l) if show==1 then printf("POPX first\n"), disp(POPX), end //Compute decimal value ('0' is placed at index '1' !!!) // Count occurences for j=1:n v=POPX(j,:), d = vec2dec(v,l), [POPX(j,l+1)]= d, //decimal values FX(d+1,2)=FX(d+1,2)+1, //occurences [POPX(j,l+2)]= fitness1(POPX(j,l+1)) //fitness values end//n [FITNESS_ALL]=fitness(POPX,l) FITNESS_ALLLOG1(i)=FITNESS_ALL if show==1 then printf("POPX second\n"), disp(POPX), end end //N for j=1:r FX(j,3)=(FX(j,2)/n)/N //Normal frequency FX(j,4)=(FX(j,2)/n)/(((1/2^l)*N)/100) // Frequency as percentage end//m disp(FX) FREQ1=tabul(FX(:,3)) MEAN1=mean(FREQ1(:,1)) STD1= stdevf(FREQ1(:,1), FREQ1(:,2)) FREQ=tabul(FX(:,4)) MEAN=mean(FREQ(:,1)) STD= stdevf(FREQ(:,1), FREQ(:,2)) MAX=max(FREQ(:,1)) MIN=min(FREQ(:,1)) DIST=MAX-MIN DIST2=DIST/2 printf("Number of Events n * N = %d\n",n*N) // Compute the maximal value and the percentage of success [MAXFIT]=maxfit02(l,n) FITNESS_ALL_PERC1=[] for i=1:N, FITNESS_ALL_PERC1(i) = FITNESS_ALLLOG1(i)/(MAXFIT/100), end // Show graphical results if show==2 then clf(), xdel, plot2d([1:1:N], [FITNESS_ALL_PERC1]), end endfunction //*************************************************** //Function to repeat frequency checks // // Input: // l := length of strings // n := size of population // N := Number of events // K := how often steps should be applied // // Output: // FX := Vector of freuencies ordered along the decimal values of a string (1st column) // FREQ := Vector of frequencies and their occurences // MEAN := The mean value of the differences // STD := The standard deviation of the difference // DIST2 := Half distances between MAX and MIN function[MEANX, DIST2X, MEAN1X, STD1X,DIST20, MEAN0]=freqCheckRep(l,n,N,K,show) // Counter for MEANX, DIST2X... MEAN1X=[] STD1X =[] MEANX=[] DIST2X =[] for j=1:K [DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,FX, POPX]=frequCheck(l,n,N,show) MEAN1X(j,1) = MEAN1 STD1X(j,1) = STD1 MEANX(j,1) = MEAN DIST2X(j,1) = DIST2 end //j MEAN10 = mean(MEAN1X(:,1)) STD10 = mean(STD1X(:,1)) MEAN0 = mean(MEANX(:,1)) DIST20 = mean(DIST2(:,1)) // Show Plot // clf(), xdel, // plot2d([DIST2SUM(:,1)], [DIST2SUM(:,4)]) endfunction //*************************************************** //Function to repeat frequency checks with gassimple0() // // Input: // l := length of strings // n := size of population // N := Number of events // K := how often steps should be applied // // Output: // FX := Vector of freuencies ordered along the decimal values of a string (1st column) // FREQ := Vector of frequencies and their occurences // MEAN := The mean value of the differences // STD := The standard deviation of the difference // DIST2 := Half distances between MAX and MIN function[MEANX, DIST2X, MEAN1X, STD1X,DIST20, MEAN0]=freqCheckRep0(POP,l,p,n,N,K, MThreshold,show) // Counter for MEANX, DIST2X... MCount = 0 MEAN1X=[] STD1X =[] MEANX=[] DIST2X =[] XDIM=[] for j=1:K XDIM(j,1)=j*N [MCount,DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,POP,FX]=gasimple0(POP,l,p,n,N, j,MThreshold,MCount,show) disp(MCount) MEAN1X(j,1) = MEAN1 STD1X(j,1) = STD1 MEANX(j,1) = MEAN DIST2X(j,1) = DIST2 end //j MEAN10 = mean(MEAN1X(:,1)) STD10 = mean(STD1X(:,1)) MEAN0 = mean(MEANX(:,1)) DIST20 = mean(DIST2(:,1)) // Show Plot if show==2 then clf(), xdel, plot2d([XDIM(:,1)], [DIST2X(:,1)]), end endfunction //*************************************************** //Function to compute the changing quality of solution sets // // Input: // l := length of strings // n := size of population // N := Number of events // K := how often should be incremented // // Output: // MEAN=:= Mean value of events // DIST20 := Mean value of 1/2-Distance of events between MAX and MIN // DIST2X := Series of DIST0s // MEANX := Series of MEANOs // FX := Summary of all values function[MDX]=solutionSetApproach(l,n,N,K,show) // Counter MDX=zeros(K,6) for j=1:K disp(j) [DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,FX, POPX]=frequCheck(l,n,j*N,show) MDX(j,1) = j*N MDX(j,2) = MEAN MDX(j,3) = DIST2 MDX(j,4) = STD MDX(j,5) = MEAN1 MDX(j,6) = STD1 end //j // Show Plot clf(), xdel, x=[MDX(:,1)], y1=[MDX(:,3)], subplot(3,1,1) plot2d(x,y1), xtitle=("DIST2[%]") y2=[MDX(:,4)],subplot(3,1,2) plot2d(x,y2), xtitle=("STD[%]") y3=[MDX(:,6)], subplot(3,1,3) plot2d(x,y3), xtitle=("STD") endfunction //*************************************************** //Function to integrate multiple frequency simulations // // Input: // l := length of strings // n := size of population // step := Inrement for number of events // K := how often steps should be applied // // Output: // DIST2SUM := Summary of distances DIST2 with regard to number of events function[DIST2SUM]=distSum(l,n,step,K,show) // Counter for DIST2 DIST2SUM=zeros(K,4) for j=1:K disp(j) [DIST2,STD, MEAN, FREQ,FX, POPX]=frequCheck(l,n,j*step,show) DIST2SUM(j,1) = j*step DIST2SUM(j,2) = MEAN DIST2SUM(j,3) = DIST2 DIST2SUM(j,4) = DIST2SUM(j,3)/(MEAN/100) end //j // Show Plot clf(), xdel, plot2d([DIST2SUM(:,1)], [DIST2SUM(:,4)]) endfunction //*************************************************** //Function to compute the overall fitness of a matrix POP // Sum up all the l+2-th positions // l := length of fitness strings function[FITNESS_ALL]=fitness(POP,l) FITNESS_ALL=0; [r,c]=size(POP); for j=1:r FITNESS_ALL=FITNESS_ALL+POP(j,l+2) end endfunction //*************************************************** //Function to compute the maximal individual fitness of a matrix POP // Compute all the l+2-th positions // l := length of fitness strings function[MFITNESS]=maxfitness(POP,l) // Extracts column l+2 and selects the maximal values MFITNESS = max(POP([:],l+2)) endfunction //*************************************************** //Function to compute the relative fitness of a string // along with the average fitness // Compute all the l+3-th positions // l := length of fitness strings function[POP,AFITNESS]=rfitness(POP,l, FITNESS_ALL) [r,c]=size(POP); for j=1:r POP(j,l+3)=POP(j,l+2)/FITNESS_ALL end AFITNESS=FITNESS_ALL/r endfunction //*************************************************** //Function to compute the fitness ratio r=f1/f0 // // l := length of fitness strings // f1 := fitness of all individuals with a allele value '1' at a position // f0 := fitness of all individuals with a allele value '0' at a position // r = f1/f0 // // Assumption: POP has the relativ fitness available at l+3 function[POP,ratio,f1,f0]=fratio10(POP,l,show) [r,c]=size(POP); f1=0; f0=0; for i=1:r for j=1:l if (POP(i,j)== 1) then f1=f1+POP(i,l+3), else f0=f0+POP(i,l+3) end end //j end//i ratio=f1/f0 if show == 2 then disp(POP), disp(f1), disp(f0), disp(ratio), end endfunction //************************************************************* // Computing the worst case convergence time tcw // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) // POP := a predefined population (can be done automatically) // // f1 := fitness of all individuals with a allele value '1' at a position // f0 := fitness of all individuals with a allele value '0' at a position // r = f1/f0 // tcw := worst case convergence time // // Assumption: POP has the relativ fitness available at l+3 function[POP,tcw,tcabin,ratio,f1,f0]=ftcw(POP,l,p,n,show) for i=1:n, v=POP(i,:), d= vec2dec(v,l), //decimal value [POP(i,l+1)]= d, //decimal value end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [FITNESS_ALL]=fitness(POP,l) [MFITNESS]=maxfitness(POP,l) [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) [POP,ratio,f1,f0]=fratio10(POP,l,show) tcw = log((n-1)^2)/log(ratio) tcabin = log(n-1)/log(ratio) endfunction //*************************************************** //Function to compute the relative count of a string // along with the realized new count by 'rounding up' // Compute all the l+4 and l+5-th positions // // Problem: if the sum of the new proposals is bigger than n // Solution: // l := length of fitness strings // n := number of individuals in POP function[POP]=newMember(POP,l,n) [r,c]=size(POP); for j=1:r POP(j,l+4)=POP(j,l+3)*n POP(j,l+5)=round(POP(j,l+4)) end //Check for sum of proposals r2 = 0 for j=1:r, r2=r2+POP(j,l+5), end while r2>n z=round(rand()*r) if z==0 then z=1,end //if POP(z,l+5) = POP(z,l+5)-1 r2=0, for j=1:r, r2=r2+POP(j,l+5), end//for end//while r2 endfunction //*************************************************** // Dummy function to replace newMember in simplegas0 // // along with the realized new count by 'rounding up' // Compute all the l+4 and l+5-th positions // l := length of fitness strings // n := number of individuals in POP function[POP]=newMember0(POP,l,n) [r,c]=size(POP); for j=1:r POP(j,l+4)=1 POP(j,l+5)=1 end endfunction //*************************************************** // Make a copy of a string from row_old to row_new: L --> R // where the copy starts at position l+p+1 // l := length of string // j := row_old // k := row_new // p := number of parameters (P=5) function[POP]=strcpyLR(POP,l,p,j,k) show=0 [r,c]=size(POP); // Testing the boundaries if (j < 1) then if show == 1 then printf('ERROR: Position of old String outside of Matrix! r = %d, j = %d', r,j),end, j=1, elseif (j > r) then if show == 1 then printf('ERROR: Position of old String outside of Matrix!, r = %d, j = %d', r,j),end, j=r, elseif (k < 1) then if show == 1 then printf('ERROR: Position of new String outside of Matrix!, r = %d, k = %d', r,k),end, k=1 elseif (k > r) then if show == 1 then printf('ERROR: Position of new String outside of Matrix!, r = %d, k = %d', r,k), end, k=r end // Making a copy for i=1:l POP(k,l+p+i) = POP(j,i) end endfunction //*************************************************** // Make a copy of a string from row_old to row_new: R --> L // where the copy starts at position l+p+1 // l := length of string // j := row_old // k := row_new // p := number of parameters (p=5) function[POP]=strcpyRL(POP,l,p,j,k) show=0 [r,c]=size(POP); // Testing the boundaries if (j < 1) then if show == 1 then printf('ERROR: Position of old String outside of Matrix! r = %d, j = %d', r,j),end, j=1 elseif (j > r) then if show == 1 then printf('ERROR: Position of old String outside of Matrix!, r = %d, j = %d', r,j),end, j=r elseif (k < 1) then if show == 1 then printf('ERROR: Position of new String outside of Matrix!, r = %d, k = %d', r,k), end, k=1 elseif (k > r) then if show == 1 then printf('ERROR: Position of new String outside of Matrix!, r = %d, k = %d', r,k),end, k=r end // Making a copy for i=1:l POP(j,i) = POP(k,l+p+i) end endfunction //*************************************************** // Make a copy of all strings from the old rows to the new ones // where the new copies will start at position l+p+1 // l := length of string // j := row_old // k := row_new // p := number of parameters // r2 := memory of position for new strings function[POP]= newPop(POP,l,p,n) [r,c]=size(POP); r2 = 1 for j=1:r if POP(j,l+7) > 0 then r3 = POP(j,l+7)-1, for k=r2:(r2+r3), //printf('\n k = %d\n', k) [POP]=strcpyLR(POP,l,p,j,k) end r2=k+1 end end endfunction //*************************************************** // Make a copy of all strings from the old rows to the new ones // where the new copies will start at position l+p+1 // // Problem: If the sum of new copies is greater than n then this causes a fault! // Solution: At least check the sum before operation // // l := length of string // j := row_old // k := row_new // p := number of parameters // r2 := memory of position for new strings function[POP]= newPop2(POP,l,p,n) [r,c]=size(POP); //Check for sum of proposals r2 = 0 for j=1:r, r2=r2+POP(j,l+5), end if r2 >n then error('In newPop2() sum of proposals >n!!!'),end r2 = 1 for j=1:r if POP(j,l+5) > 0 then r3 = POP(j,l+5)-1, for k=r2:(r2+r3), //printf('\n k = %d\n', k) [POP]=strcpyLR(POP,l,p,j,k) end r2=k+1 end end endfunction //*************************************************** // Prepare the crossover operations within a population POP // by randomly selecting the new strings and copy them onto the // POP base area [1,l] // // General Assumption: the number of members n is even!!! // // The strings are assumed to be in the area starting at l+p+1 // The parameter at l+p has been changed to a flag '1' := not yet mated // j := target position in the base area for copy action function[POP]= crossoverPrep(POP,l,p,n) [r,c]=size(POP); // Set a flag at position l+7 with '1' for j=1:r POP(j,l+7)=1 end // Select randomly the strings for transfer j = 1 //:= Baseline for filling up with strings while(j < n+1) S = 1 //:= Flag for searching a string to be copied while(S == 1) k1 = round(n * rand()) if k1==0 then k1=1 end //if if POP(k1,l+7) <> 0 then POP=strcpyRL(POP,l,p,j,k1) //mshow(POP,n,l+p+l) POP(k1,l+7) = 0, S=0 else S=1 end //if end //while == S // printf('k1 = %d\n',k1) // printf('j = %d\n',j) j = j+1 // mshow(POP,n,l+p+l) end // while == j endfunction //*************************************************** // Apply a crossover operation onto two adjacent strings at intersection x // The strings are assumed to be randomly paired in the area at position 1 // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) function[POP]= crossover(POP,l,p,n) // Take a first row r=1 while(r < n) // printf('r = %d\n',j) // Select randomly a point in this row x = round(l * rand()), if(x == 0) then x=1,end, // printf('x = %d\n',x) for i=x:l, m=POP(r+1,i) POP(r+1,i) = POP(r,i) POP(r,i) = m end r=r+2 end endfunction //*************************************************** // Apply a mutation operation onto one randomly selected string // Select randomly a point in the string and replace it's // value by it's inverse '1' --> '0', '0' --> '1' function[POP]= mutation(POP,l,p,n) // Find a point in a string c = round(l * rand()) if c == 0 then c=1 // This does not change anything end // Find a row in the table r = round(n * rand()) if r == 0 then r=1 // This does not change anything end // replace the value if(POP(r,c) == 1) then POP(r,c) = 0, else POP(r,c) = 1 end printf('mutation point at row = %d col = %d\n',r,c) endfunction //*************************************************** // Print the content of a matrix // // d := depth of matrix 'downwards' // w := width of matrix from left to right function[M]= mshow(M,d,w) for j=1:d,// printf('\n j= %3.1d : ',j) for i=1:w, //printf(' %3.1d ',M(j,i)) end //printf('\n') end endfunction //*************************************************** // Translate strings with binaries '0', '1' as decimal numbers // // v = vector of binaries from a POP-matrix (left part) // l := length of string // D := decimal computed out of binaries function[D]= vec2dec(v,l) str=string(v(1:l)) for i=2:l, str(1)=str(1)+str(i) end D=bin2dec(str(1)) endfunction //*************************************************** // Translate strings with binaries '0', '1' as decimal numbers // // v = vector of binaries from a POP-matrix (left part) // l := length of string // D := decimal computed out of binaries function[D]= vec22dec(v,l1,l2) str=string(v(l1:l2)) for i=2:l2-l1+1, str(1)=str(1)+str(i) end D=bin2dec(str(1)) endfunction //*************************************************** // Translate strings with binaries '0', '1' and n-many compartments as decimal numbers // // v = vector of binaries from a POP-matrix (left part) // l := length of string // b := number of bins for one number // D := Array of decimals function[D]= vec222dec(v,l,b,show) D=[] k=1 for j=1:b:l+1-b if show==1 then printf('j = %d ',j),end str=string(v(j:j+b-1)) for i=2:b, str(1)=str(1)+str(i) end //i if show==1 then printf('k = %d ',k),end D(k)=bin2dec(str(1)) if show==1 then printf('D(k) = %d\n',D(k)),end k=k+1 end //j endfunction //*************************************************** // maxfit01 builds a string of length l and computes the decimal value // // l := length of string // MF := decimal computed out of binaries function[MF]= maxfit01(l) FitString="" v=ones(1,l) str=string(v) for i=1:l, FitString = FitString + str(i) end MF=bin2dec(FitString) endfunction //*************************************************** // maxfi02 computes with maxfit01 the decimal value of a string with length l and // then compares the maximal value according to the fitness function y=x^2 // multiplied by the number n of strings in a population // // l := length of string // MF := decimal computed out of binaries function[MAXFIT]= maxfit02(l,n) MAXFIT=n*(maxfit01(l)^2) endfunction //*************************************************** // Simple fitness-function f=(x^2) // // D := decimal computed out of binaries function[F]= fitness1(D) F=D^2 endfunction //*************************************************** // Fitness-function f for a world W which is realized as a table t:D --> R // IF (d,r) in t then f(d1,d2)=100, ELSE f(d1,d2)=0 // // d1,d2 := decimal values of population POP // W := a world organized as a table // n := number of elements in W // F := fitness value function[F]= fitness2(d1,d2,W) [r,c]=size(W) upper =10 goal = upper*2 i=1 while i < r+1, if W(i,1) == d1 & W(i,2) == d2 then F=goal, i=r+1, else F=round(upper * rand()), i=i+1, end //if end//while endfunction //*************************************************** // Fitness-function f for a world W which is realized as a table t:D --> R // IF (d,r) in t then f(d1,d2)=goal, ELSE f(d1,d2)=random(upper) // // D := Array of pairs of decimal values // W := a world organized as a table // F := fitness value function[F,goal]= fitness3(D,W) [r,c]=size(W) [r2,c2]=size(D) F2=[] upper =10 // range of numbers for non-goals [0,upper] goal = upper*2 j=1 //Index for D for i=1:r //Look to every pair in W if W(i,1) == D(j) & W(i,2) == D(j+1) then F2(i)=goal, else F2(i)=round(upper * rand()), end //if j=j+2 end // for i F=sum(F2) endfunction //************************************************************* // Function to compute only the fitness of a population // // Input: // p := number of parameters between strings left and right (actually p=5) // l := length of strings // n := number of elements in POP // Output // POP := a population // MFITNESS := Maximal Fitness of each individual function[POP,MFITNESS, FITNESS_ALL]=popFit(l,p,n) [POP]=popgen(n,l+p) MFITNESS=0 for j=1:n, v=POP(j,:), [POP(j,l+1)]= vec2dec(v,l) end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [FITNESS_ALL]=fitness(POP,l) [MFITNESS]=maxfitness(POP,l) endfunction //************************************************************* // Function to compute the maximal fitness of an individuum // for several populations // // Input: // p := number of parameters between strings left and right (actually p=5) // l := length of strings // n := number of elements in POP // run := number of choices // Output // POP := a population // MFITNESS := Maximal Fitness of each individual // MFITNESSX := Array with maximal fitness function[POP,MEAN, MFITNESSX]=popFitX(l,p,n,run) MAXFITUP = ((2^l)^2) for k=1:run [POP]=popgen(n,l+p) for j=1:n, v=POP(j,:), [POP(j,l+1)]= vec2dec(v,l) end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [MFITNESSX(k)]=maxfitness(POP,l) end for k=1:run MFITNESSX(k)=MFITNESSX(k)/(MAXFITUP/100) end MEAN = mean(MFITNESSX) // Show graphical results //clf(), xdel, //plot2d([1:1:run], MFITNESSX) endfunction //************************************************************* // Function to compute the maximal fitness of an individuum // for several populations // for increasing n // // Input: // p := number of parameters between strings left and right (actually p=5) // l := length of strings // n := number of elements in POP // run := number of choices // Output // POP := a population // MFITNESS := Maximal Fitness of each individual // MFITNESSX := Array with maximal fitness function[MEAN, MEANX]=popFitXN(l,p,m,run) for i=1:m [POP,MEAN, MFITNESSX]=popFitX(l,p,i,run) MEANX(i) = MEAN end MEAN = mean(MEANX) // Show graphical results clf(), xdel, plot2d([1:1:m], MEANX) endfunction //************************************************************* // Function to compute the maximal fitness of a population // for several generations // // Input: // p := number of parameters between strings left and right (actually p=5) // l := length of strings // n := number of elements in POP // run := number of generations // Output // POP := a population // FITNESS_ALL := Maximal Fitness of a population // FITNESS_ALLX := Array with maximal fitness // MEAN := Mean value of array function[STD,MEAN, FITNESS_ALLX]=popMaxFitX(l,p,n,run) MAXFITUP = ((2^l)^2)*n for k=1:run [POP]=popgen(n,l+p) for j=1:n, v=POP(j,:), [POP(j,l+1)]= vec2dec(v,l) end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [FITNESS_ALLX(k)]=fitness(POP,l) end for k=1:run FITNESS_ALLX(k)=FITNESS_ALLX(k)/(MAXFITUP/100) end MEAN = mean(FITNESS_ALLX) mf=tabul(FITNESS_ALLX) STD=stdevf(mf([:],1), mf([:],2)) // Show graphical results //clf(), xdel, //plot2d([1:1:run], MFITNESSX) endfunction //************************************************************* // All Functions Unified // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) // POP := a predefined population (can be done automatically) // run := number of cycles // MThreshold := number of cycles which have to be waited until the mutation operators will be applied function[POP,FITNESS_ALLLOG]=gasimple(POP,l,p,n,run, MThreshold) MCount = 0 FITNESS_ALLLOG=[] for cyc = 1:run [M]= mshow(POP,n,l+p) for j=1:n, v=POP(j,:), [POP(j,l+1)]= vec2dec(v,l) end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [FITNESS_ALL]=fitness(POP,l) FITNESS_ALLLOG(cyc)=FITNESS_ALL [MFITNESS]=maxfitness(POP,l) [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) [POP]=newMember(POP,l,n) [POP]=newPop(POP,l,p,n) [POP]= crossoverPrep(POP,l,p,n) [POP]= crossover(POP,l,p,n) printf('Mutationcount = %d\n', MCount) if(MCount > MThreshold) then [M]= mshow(POP,n,l+p), [POP]= mutation(POP,l,p,n), [M]= mshow(POP,n,l+p), MCount=0, end MCount = MCount + 1 end //for == run // Compute the maximal value and the percentage of success [MAXFIT]=maxfit02(l,n) FITNESS_ALL_PERC=[] for i=1:run, FITNESS_ALL_PERC(i) = FITNESS_ALLLOG(i)/(MAXFIT/100), end // Show graphical results clf(), xdel, //plot2d([1:1:run], FITNESS_ALLLOG) plot2d([1:1:run], FITNESS_ALL_PERC) endfunction //************************************************************* // A system with GA without any additional parameters // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) // POP := a predefined population (can be done automatically) // N := number of cycles // MThreshold := number of cycles which have to be waited until the mutation operators will be applied function[FITNESS_ALL_PERC,POP]=ga(POP,l,p,n,N, MThreshold,show) MCount = 0 FITNESS_ALLLOG=[] for cyc = 1:N for j=1:n, v=POP(j,:), d= vec2dec(v,l), //decimal value [POP(j,l+1)]= d, //decimal value end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [FITNESS_ALL]=fitness(POP,l) FITNESS_ALLLOG(cyc)=FITNESS_ALL [MFITNESS]=maxfitness(POP,l) [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) [POP]=newMember(POP,l,n) [POP]=newPop(POP,l,p,n) [POP]= crossoverPrep(POP,l,p,n) [POP]= crossover(POP,l,p,n) if show==1 then printf('Mutationcount = %d\n', MCount), end if(MCount > MThreshold) then [POP]= mutation(POP,l,p,n), if show==1 then disp(POP), end MCount=0, end MCount = MCount + 1 end //for == N printf("Number of Events n * N = %d\n",n*N) // Compute the maximal value and the percentage of success [MAXFIT]=maxfit02(l,n) FITNESS_ALL_PERC=[] for i=1:N, FITNESS_ALL_PERC(i) = FITNESS_ALLLOG(i)/(MAXFIT/100), end // Show graphical results if show==2 then clf(), xdel, plot2d([1:1:N], FITNESS_ALL_PERC), end endfunction //************************************************************* // A system with GA // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) // POP := a predefined population (can be done automatically) // N := number of cycles // MThreshold := number of cycles which have to be waited until the mutation operators will be applied function[FITNESS_ALL_PERC,DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,FX, POP]=ga0(POP,l,p,n,N, MThreshold,show) // Install the counters r= (2^l) c= 6 FX = zeros(r,c) MCount = 0 FITNESS_ALLLOG=[] // Generate an internal index for display for j=1:r FX(j,1)=j-1 end for cyc = 1:N for j=1:n, v=POP(j,:), d= vec2dec(v,l), //decimal value [POP(j,l+1)]= d, //decimal value FX(d+1,2)=FX(d+1,2)+1, //occurences end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [FITNESS_ALL]=fitness(POP,l) FITNESS_ALLLOG(cyc)=FITNESS_ALL [MFITNESS]=maxfitness(POP,l) [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) [POP]=newMember(POP,l,n) [POP]=newPop(POP,l,p,n) [POP]= crossoverPrep(POP,l,p,n) [POP]= crossover(POP,l,p,n) if show==1 then printf('Mutationcount = %d\n', MCount), end if(MCount > MThreshold) then [POP]= mutation(POP,l,p,n), if show==1 then disp(POP), end MCount=0, end MCount = MCount + 1 end //for == N for j=1:r FX(j,3)=(FX(j,2)/n)/N //Normal frequency FX(j,4)=(FX(j,2)/n)/(((1/2^l)*N)/100) // Frequency as percentage end//m disp(FX) FREQ1=tabul(FX(:,3)) MEAN1=mean(FREQ1(:,1)) STD1= stdevf(FREQ1(:,1), FREQ1(:,2)) FREQ=tabul(FX(:,4)) MEAN=mean(FREQ(:,1)) STD= stdevf(FREQ(:,1), FREQ(:,2)) MAX=max(FREQ(:,1)) MIN=min(FREQ(:,1)) DIST=MAX-MIN DIST2=DIST/2 printf("Number of Events n * N = %d\n",n*N) // Compute the maximal value and the percentage of success [MAXFIT]=maxfit02(l,n) FITNESS_ALL_PERC=[] for i=1:N, FITNESS_ALL_PERC(i) = FITNESS_ALLLOG(i)/(MAXFIT/100), end // Show graphical results if show==2 then clf(), xdel, plot2d([1:1:N], FITNESS_ALL_PERC), end endfunction //************************************************************* // GA algorithm without fitness // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) // POP := a predefined population (can be done automatically) // N := number of events // K := Number of repetitions from calling function // MThreshold := number of events which have to be waited until the mutation operators will be applied function[MCount,DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,POP,FX]=gasimple0(POP,l,p,n,N, K,MThreshold,MCount,show) // Install the counters r= (2^l) c= 5 FX = zeros(r,c) // Generate an internal index for display for j=1:r FX(j,1)=j-1 end for cyc = 1:N*K //Start Cyle for events if show ==1 then disp(cyc), end for j=1:n v=POP(j,:), d = vec2dec(v,l), [POP(j,l+1)]= d, FX(d+1,2)=FX(d+1,2)+1 end//n if show ==1 then printf("\nFX after Count\n\n"), disp(FX), end [POP]=newMember0(POP,l,n) //Prepares copies by inserting '1' at position l+4 and l+5 [POP]=newPop(POP,l,p,n) //Copies old strings to the right [POP]= crossoverPrep(POP,l,p,n) [POP]= crossover(POP,l,p,n) if show ==1 then printf('Mutationcount = %d\n', MCount),disp(POP), end if(MCount > MThreshold) then [POP]= mutation(POP,l,p,n), if show ==2 then printf('Mutationcount (MCount, cyc)= (%d,%d)\n', MCount, cyc), end MCount = 0, end MCount = MCount + 1 if show==1 then printf("POP after Operators\n"), disp(POP), end end //cyc=N for j=1:r FX(j,3)=(FX(j,2)/n)/(N*K) //Normal frequency FX(j,4)=(FX(j,2)/n)/(((1/2^l)*(N*K))/100) // Frequency as percentage end//m disp(FX) FREQ1=tabul(FX(:,3)) MEAN1=mean(FREQ1(:,1)) STD1= stdevf(FREQ1(:,1), FREQ1(:,2)) FREQ=tabul(FX(:,4)) MEAN=mean(FREQ(:,1)) STD= stdevf(FREQ(:,1), FREQ(:,2)) MAX=max(FREQ(:,1)) MIN=min(FREQ(:,1)) DIST=MAX-MIN DIST2=DIST/2 printf("Number of Events n * N * K = %d\n",n*N*K) endfunction //************************************************************* // GA algorithm with fitness and conditioned mutation // // Works like 'normal' GA but uses mutation only if the actual maximal fitness // stays 'constant' over MT-many cycles and is below te theoretical maximum // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) // POP := a predefined population (can be done automatically) // N := number of events // MT := mutation trigger; number of values to monitor for to trigger mutation function[FITNESS_ALL_PERC,DIST2,STD, MEAN, FREQ,STD1, MEAN1, FREQ1,FX, POP]=ga1(POP,l,p,n,N, MT,show) //Theoretical Maximum for Fitness according to y=x^2 TMaxFit = ((2^l)-1)^2*n // Install the counters r= (2^l) c= 6 FX = zeros(r,c) MCount = 0 FITNESS_ALLLOG=[] MTriggerEstimator=[] // Generate an internal index for display for j=1:r FX(j,1)=j-1 end for cyc = 1:N for j=1:n, v=POP(j,:), d= vec2dec(v,l), //decimal value [POP(j,l+1)]= d, //decimal value FX(d+1,2)=FX(d+1,2)+1, //occurences end for i=1:n,[POP(i,l+2)]= fitness1(POP(i,l+1)) end [FITNESS_ALL]=fitness(POP,l) FITNESS_ALLLOG(cyc)=FITNESS_ALL [MFITNESS]=maxfitness(POP,l) if show ==1 then disp(FITNESS_ALL), end MTriggerEstimator(cyc)=FITNESS_ALL if show ==1 then disp(MTriggerEstimator(cyc)), end [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) [POP]=newMember(POP,l,n) [POP]=newPop(POP,l,p,n) [POP]= crossoverPrep(POP,l,p,n) [POP]= crossover(POP,l,p,n) if show ==1 then disp(MCount), end if(MCount > MT) then if show ==1 then disp(MTriggerEstimator(cyc)), disp(MTriggerEstimator(cyc-MT)), disp(TMaxFit), end if (MTriggerEstimator(cyc) == MTriggerEstimator(cyc-MT) & MTriggerEstimator(cyc) < TMaxFit) then [POP]= mutation(POP,l,p,n), if show==1 then disp(POP), end end MCount=0, end MCount = MCount + 1 end //for == N for j=1:r FX(j,3)=(FX(j,2)/n)/N //Normal frequency FX(j,4)=(FX(j,2)/n)/(((1/2^l)*N)/100) // Frequency as percentage end//m disp(FX) FREQ1=tabul(FX(:,3)) MEAN1=mean(FREQ1(:,1)) STD1= stdevf(FREQ1(:,1), FREQ1(:,2)) FREQ=tabul(FX(:,4)) MEAN=mean(FREQ(:,1)) STD= stdevf(FREQ(:,1), FREQ(:,2)) MAX=max(FREQ(:,1)) MIN=min(FREQ(:,1)) DIST=MAX-MIN DIST2=DIST/2 printf("Number of Events n * N = %d\n",n*N) // Compute the maximal value and the percentage of success [MAXFIT]=maxfit02(l,n) FITNESS_ALL_PERC=[] for i=1:N, FITNESS_ALL_PERC(i) = FITNESS_ALLLOG(i)/(MAXFIT/100), end // Show graphical results if show==2 then clf(), xdel, plot2d([1:1:N], FITNESS_ALL_PERC), end endfunction //************************************************************* // A system with GA for a world with a table // // A world W organized as a table maps a set D into a set R // l := length of strings // p := number of cells between string <1...l> and <l+p+1, ..., l+p+l> // Pos 1-l := String // Pos l+1 := Decimal Value of binary string // Pos l+2 := Fitness for l+1 and l+p // Pos l+3 := Percentage of overall fitness // Pos l+4 := Expected count according to fitness // Pos l+5 := Realized count // Pos l+6 := 2nd Decimal value of a compund fitness // Pos l+p := flag for crossover // p := 7 // POP := a predefined population (can be done automatically) // n := number of elements in POP // N := number of cycles // MT := number of cycles which have to be waited until the mutation operators will be applied function[FITNESS_ALL_PERC,POP]=gaw1(W,POP,l,p,n,N, MT,l1,show) MCount = 0 FITNESS_ALLLOG=[] for cyc = 1:N if show == 1 then disp('cyc ='),disp(cyc),disp(' '),end for j=1:n, v=POP(j,:), d1= vec22dec(v,1,l1), //decimal value1 d2= vec22dec(v,l1+1,l), //decimal value2 [POP(j,l+1)]= d1, //decimal value [POP(j,l+6)]= d2, //decimal value end if show == 1 then printf('New Numbers in POP '), disp(POP), end for i=1:n,[POP(i,l+2)]= fitness2(POP(i,l+1),POP(i,l+6),W) end [FITNESS_ALL]=fitness(POP,l) if show == 1 then printf('FITNESS_ALL = %f\n',FITNESS_ALL), end FITNESS_ALLLOG(cyc)=FITNESS_ALL [MFITNESS]=maxfitness(POP,l) [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) if show == 1 then disp(POP), end [POP]=newMember(POP,l,n) [POP]=newPop2(POP,l,p,n) if show == 1 then disp(POP), end [POP]= crossoverPrep(POP,l,p,n) [POP]= crossover(POP,l,p,n) if show == 1 then disp(POP),printf('MCount = %d\n',MCount), end if(MCount > MT) then [POP]= mutation(POP,l,p,n), printf('Mutation at cycle = %d\n',cyc), if show==1 then disp(POP), end MCount=0, end MCount = MCount + 1 end //for == N printf("Number of Events n * N = %d\n",n*N) // Compute the maximal value and the percentage of success MAXFIT=20*n //The value '20' is set in the function fitness2() FITNESS_ALL_PERC=[] for i=1:N, FITNESS_ALL_PERC(i) = FITNESS_ALLLOG(i)/(MAXFIT/100), end // Show graphical results if show==2 | show == 1 then clf(), xdel, plot2d([1:1:N], FITNESS_ALL_PERC), end endfunction //************************************************************* // A system with GA for a world with a table and with // multi-compartment genomes // // A world W organized as a table maps a set D into a set R // l := length of strings // p := number of cells between string <1...l> and <l+p+1, ..., l+p+l> // Pos 1-l := String // Pos l+1 := Decimal Value of binary string // Pos l+2 := Fitness for l+1 and l+p // Pos l+3 := Percentage of overall fitness // Pos l+4 := Expected count according to fitness // Pos l+5 := Realized count // Pos l+6 := 2nd Decimal value of a compund fitness // Pos l+7 := flag for crossover // p := 7 // POP := a predefined population (can be done automatically) // n := number of elements in POP // N := number of cycles // MT := number of cycles which have to be waited until the mutation operators will be applied // V := version function[FITNESS_ALL_PERC,MAXGENOMECOUNT,POP,V]=gaw2(W,POP,l,n,N, MT,b,show,V) p = 7 MCount = 0 FITNESS_ALLLOG=[] for cyc = 1:N if show == 1 then disp('cyc ='),disp(cyc),disp(' '),end for i=1:n , v=POP(i,:), [D]= vec222dec(v,l,b,show), [POP(i,l+2),goal]= fitness3(D,W), end //for i if show == 1 then printf('New Numbers in POP '), disp(POP), end [FITNESS_ALL]=fitness(POP,l) if show == 1 then printf('FITNESS_ALL = %f\n',FITNESS_ALL), end FITNESS_ALLLOG(cyc)=FITNESS_ALL [MFITNESS]=maxfitness(POP,l) [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) // if show == 1 then disp(POP), end [POP]=newMember(POP,l,n) if show == 1 then disp(POP), end [POP]=newPop2(POP,l,p,n) if show == 1 then disp(POP), end [POP]= crossoverPrep(POP,l,p,n) [POP]= crossover(POP,l,p,n) if show == 1 then disp(POP),printf('MCount = %d\n',MCount), end if(MCount > MT) then [POP]= mutation(POP,l,p,n), printf('Mutation at cycle = %d\n',cyc), if show==1 then disp(POP), end MCount=0, end MCount = MCount + 1 end //for == N printf("Number of Events n * N = %d\n",n*N) // Compute the maximal value and the percentage of success [r,c]=size(POP) [r2,c2]=size(W) MAXFIT=goal*r2*r //The value 'goal' is set in the function fitness3() MAXFITGENOME = goal*r2 // To count the complete genomes in a population FITNESS_ALL_PERC=[] MAXGENOMECOUNT = 0 for i=1:r, if POP(i,l+2) == MAXFITGENOME then MAXGENOMECOUNT =MAXGENOMECOUNT+1, end, end for i=1:N, FITNESS_ALL_PERC(i,1)=i, FITNESS_ALL_PERC(i,2) = FITNESS_ALLLOG(i)/(MAXFIT/100), end // Show graphical results if show==2 | show == 1 then clf(), xdel, scf(V), plot2d([1:1:N], FITNESS_ALL_PERC(:,2)), end endfunction //************************************************************* // Computing the worst case convergence time tcw // // p := number of parameters between strings left and right (actually p=5) // l := length of strings ( actually l=5) // n := number of elements in POP (actually n=4) // POP := a predefined population (can be done automatically) // // f1 := fitness of all individuals with a allele value '1' at a position // f0 := fitness of all individuals with a allele value '0' at a position // r = f1/f0 // tcw := worst case convergence time // // Assumption: POP has the relativ fitness available at l+3 function[POP,tcw,tcabin,ratio,f1,f0]=ftcw2(W,POP,l,p,n,l1,show) for j=1:n, v=POP(j,:), d1= vec22dec(v,1,l1), //decimal value1 d2= vec22dec(v,l1+1,l), //decimal value2 [POP(j,l+1)]= d1, //decimal value [POP(j,l+6)]= d2, //decimal value end for i=1:n,[POP(i,l+2)]= fitness2(POP(i,l+1),POP(i,l+6),W) end [FITNESS_ALL]=fitness(POP,l) [MFITNESS]=maxfitness(POP,l) [POP,AFITNESS]=rfitness(POP,l,FITNESS_ALL) [POP,ratio,f1,f0]=fratio10(POP,l,show) tcw = log((n-1)^2)/log(ratio) tcabin = log(n-1)/log(ratio) endfunction //*********************************************************************** // Examples // POP = [0 1 1 0 1 13 169 0 0 0; 1 1 0 0 0 24 576 0 0 0; 0 1 0 0 0 8 64 0 0 0; 1 0 0 1 1 19 361 0 0 0;] W=[1 3; 2 1; 3 2] W2=[1 3; 2 1; 3 2; 4 4; 5 7; 6 5]