摘自http://www.cnblogs.com/hxsyl/p/5240905.html
該文中后面的代碼是我在學校里面編寫的,當時應該是可以用的,里面的圖確實是我畫出來的,但是我現在不確定是否正確,因為現在畢業后沒有matlab工具,我后面會根據現在的matlab的代碼改寫成python的代碼,到時候再把新的代碼傳上來。(20190816補充)
1、遺傳算法介紹
遺傳算法,模擬達爾文進化論的自然選擇和遺產學機理的生物進化構成的計算模型,一種不斷選擇優良個體的算法。談到遺傳,想想自然界動物遺傳是怎么來的,自然主要過程包括染色體的選擇,交叉,變異(不明白這個的可以去看看生物學),這些操作后,保證了以后的個基本上是最優的,那么以后再繼續這樣下去,就可以一直最優了。
2、解決的問題
先說說自己要解決的問題吧,遺傳算法很有名,自然能解決的問題很多了,在原理上不變的情況下,只要改變模型的應用環境和形式,基本上都可以。但是遺傳算法主要還是解決優化類問題,尤其是那種不能直接解出來的很復雜的問題,而實際情況通常也是這樣的。
本部分主要為了了解遺傳算法的應用,選擇一個復雜的二維函數來進行遺傳算法優化,函數顯示為y=10*sin(5*x)+7*abs(x-5)+10,這個函數圖像為:

怎么樣,還是有一點復雜的吧,當然你還可以任意假設和編寫,只要符合就可以。那么現在問你要你一下求出最大值你能求出來嗎?(這個貌似可以,很容易寫出來----如果再復雜一點估計就不行了)這類問題如果用遺傳算法或者其他優化方法就很簡單了,為什么呢?說白了,其實就是計算機太笨了,同時計算速度又超快,舉個例子吧,我把x等分成100萬份,再一下子都帶值進去算,求出對應的100萬個y的值,再比較他們的大小,找到最大值不就可以了嗎,很笨吧,人算是不可能的,但是計算機可以。而遺傳算法也是很笨的一個個搜索,只不過加了一點什么了,就是人為的給它算的方向和策略,讓它有目的的算,這也就是算法了。
3、如何開始?
不明白遺傳算法的會問怎么開始呢?恩,其實每個算法都有自己的開始方式,遺傳算法也是,首先是選擇個體了。我們知道一個種群中可能只有一個個體嗎?不可能吧,肯定很多才對,這樣相互結合的機會才多,產生的后代才會多種多樣,才會有更好的優良基因,有利於種群的發展。那么算法也是如此,當然個體多少是個問題,一般來說20-100之間我覺得差不多了。那么個體究竟是什么呢?在我們這個問題中自然就是x值了。其他情況下,個體就是所求問題的變量,這里我們假設個體數選100個,也就是開始選100個不同的x值,不明白的話就假設是100個猴子吧。好了,現在有了100個猴子組成的一個種群,那么這個種群應該怎么發展才能越來越好?說到這,我們想想,如何定義這個越來越好呢?這個應該有一個評價指標吧。在我們的這個問題中,好像是對應的y值越大越好是吧。我們甚至可以給他們排個名來決定哪些好哪些不好。我們把這個叫做對於個體的適應度,這應該算是算法的后半部分才對。
4、編碼
首先明白什么是編碼?為什么要編碼?如何編碼?
好,什么是編碼?其實編碼就是把自變量(x)換一下形式而已,在這個形式下,更容易操作其他過程(比如交叉,變異什么的)而已。舉個例子吧,加入我們取x=1, 2, 3,我們可以把x編碼成x=a, b, c,我就說123對應就是abc,為什么要這樣做呢?比如問題里面你能夠獲取的都是abc的組合之類的,那么這樣編碼以后,你就可以再返回去成123來操作了。一般的編碼都是些什么呢?二進制編碼,自然數編碼,矩陣編碼。。等等,很多,不詳細寫了。而用的最多的可以說是二進制編碼,感覺這和人體DNA,基因的排列很相似。想想DNA怎么排的?不就是在兩條長鏈上一對一排的嗎?那么什么是二進制編碼?很簡單,就是1,0,1,0對應的來組合排列而已。比如:1100100010, 0011001001等等,這些都是位數長度為10的二進制編碼。再想想1在計算機的二進制形式是什么?如果八位來表示的話,是不是就是0000 0001;8是不是就是0000 1000;以此類推,那么我們這里也是這樣,把對應的x值換算成這種編碼形式,我們這里可以看到x的范圍是0-5吧,如果按照計算機這樣的方式是不是到0000 0101這里就完事了?想想這樣多短,前面五位都沒有用上多浪費呀,那么要想都用上怎么辦呢?也很簡單,我們把0000 0001不認為是1不就可以了嗎?因為1111 1111是255,那么如果說每一份為1/255的話,那么0000 0001不就是1/255了嗎?這個時候1怎樣表示了?不就是1111 1111嗎?好了我們把范圍擴大一些吧,每一份不是1/255,而是1/255*5,那么這個時候最大值是多少?是不是5,恩,這樣x編碼的范圍不就在0-5之間了嗎?這里就又有問題了,想想這樣做的話,x的最小精度是多少?就是1/255*5.雖然很小,但是在0-1/255*5之間的數你能不能取到?無論如何都取不到吧。那么就又來一個問題,怎樣去擴大這個精度呢?如果要保持0-5不變的話,只能增加位數了,把9位編碼編程10位,20位,100位,哇,夠大了吧,變成100個0,1組合,很恐怖吧,事實上,究竟是多少要視情況而定,一般20位左右感覺就可以了,雖然說越大越好,但是太大了消耗內存,速度慢了,不值。本題中,我們設置它為一個變量,先暫時取為10來實驗,好了,如果還不明白為什么要編碼看下面的吧。知道了交叉與變異,你就知道了。
5、關於交叉與變異
先說變異,什么是變異?簡單,基因發生突變就叫變異,有了編碼的概念,那就在編碼基礎上來說變異。首先就講最簡單的變異,就是個體的變異。現在以10位長的編碼來說,比如把x=3編碼一下,隨便假設為11000 10010吧,好了,在變異操作時,假設第5位變異了(說一下變異就是一位或者多位1或0變成0或1,也只能0,1之間變,沒辦法啊),那么這個時候變成什么了?是不是11001 10010再反編碼回去成x是多少呢?那肯定不是3了,變了呀,是多少是肯定可以反算回去的,這里懶得算了,就假設為3.213吧,發沒發現,這樣一來,x是不是變了?既然變了就好啊,帶到原函數(適應度函數)里面比較這兩個x值對應的哪個y值大一寫,如果后面變異后的大一些是不是就是說產生了好的變異啊,就可以在下一次個體選擇的時候選擇它了。那么想想很多x來一起變異會怎么樣呢?肯定會生成很多很多的解吧,反復這么做會怎么樣呢?只要每次都保留最優解的話,我來個循環100萬次,也總能找到解吧,當然這么多次得花多久,也不太合適,這還只是一個點位在進行變異,若果每次我讓多個點位變異呢?哇,又不可思議了,變化更大了吧。當然,變異不止如此,更多的去看專業的論文吧。知道了變異是干什么的,剩下的都好說了,好了,這還是變異,想想自然界遺傳中除了變異還有什么?交叉吧,那么交叉又是什么?
學過生物的都知道,動物交配時,部分染色體干什么了?是不是交叉了?就是把相應部分的基因交換了,你的給我,我的給你,很有愛吧。再以編碼為例,比如現在隨便從100個x值中選取兩個吧,假設正好選中了x=3和4,對應的編碼假設是11001 10101 和00101 01011,那么怎么交叉呢?我們知道每次交叉的染色體通常是一塊一塊的,恩,這里在算法設計上也來一塊一塊的吧。比如說就把位置在2,3,4號的編碼給整體交叉了,那么x=3對應的位置是100吧,x=4對應的位置是010吧,好,交換以后x=3對應的位置就變成了010,x=4對應的位置就變成100,加回去就變成什么了?x=3是不是就是10101 10101,x=4是不是就是01001 01011了。而現在,把他們再反編碼回去還是x=3和x=4嗎?顯然又不是了吧(當然也有小概率是一樣的吧,很小)。那是什么?不想算,還是假設吧,假設為3.234和4.358吧,好了新的個體是不是又來了?恩,同理,帶到適應度函數里面去吧,再取優秀個體,完事。同樣,有些專門研究這種算法的開發出來各種各樣的交叉方式,什么一個個體的前3個與后一個個體的后三個交叉,中間幾位來交叉等等,總之就是生產新個體。而這樣做的目的在哪呢?無非是三個字,隨機性,充分保證生產新個體具有隨機性,你說你的x=3變異后為3.2,3.2距離3那么近,在一些存在局部最優解的問題上就永遠跳不出局部最優解,相反,你的x=1一下子變異成了x=5,哇,好大的變化,一下從這頭到了那頭,這對於算法的廣闊搜索能力來說是非常好的。
講完了這部分,現在知道了為什么要編碼了吧?如果你不編碼,你說你想要你的x=3怎么去變異?怎么去交叉?當然也不是沒有方法,比如你生成一個小的隨機數加到x=3上,但是你想想這兩種方法哪一個更具有隨機性、普遍性?顯然的。而更多的時候交叉與變異是在一起操作的,先交叉,再變異是普遍遺傳算法的操作步驟。
6、關於選擇的問題
說完了上面的部分,再說說選擇吧,選擇是什么?就是優勝劣汰。好的留下來,差的走人,在自然界中直接gg了是吧。不停地選擇使的種群一直朝着較好的方向行走。
對應到本問題來說,遺傳算法的選擇是什么樣子?在前面說到,每次交叉或者變異是不是產生了新的個體?如果這些個體都保留下來的話,種群是要爆炸的,第一次循環可能有100個x,第二次循環就200個個體,再來那么10萬次,哇哦,多少了,好多。顯然不可能吧。而且在算法里面,我們還規定的是每次循環都必須保證都是100個個體。那么必須在200個個體中剔除100個吧。好了,問題來了,如何剔除呢?有人說很簡單,排名吧,取前100號x不就可以了嗎?排名這個東西真的好嗎?我就不信,憑什么差一點的不能選上,搞不好在下一次變異中一下子沖到了第一呢?這個問題在選擇上也有一些對應的規則,最通用的就是輪盤賭法,簡單來說就是一種概率選擇法(當然還有許多其他的方法,感興趣的自己搜相關的文獻吧,我也沒用過)。什么是輪盤賭法呢?就是把對應所有y值(適應度函數值)加起來,再用各自的y值去除以這個sum值,這樣是不是誰的概率大誰的概率小就很清楚了?然后再隨機生成一個0-1的概率值p,誰在p的范圍里面是不是就選擇誰,比如說x=3時在100個x中y的值最大,那么選擇它的概率是不是就最大,比如說是0.1(0.1小嗎?不小了好吧,想想其他的會是什么,都比0.1小,那么從概率上講,選100次的話,是不是就有10次選擇x=3,其他的都不足10次是吧,那么在下一次100個種群個體中就有10個x=3了,再來一回可能就有20個x=3的個體了。再就是30個,最后就只剩下100個x=3的個體了,它自己在哪里交叉變異已經沒有意義了,如果到了這個時候就意味着這個算法可以結束了)。再詳細點,如下圖所示吧:現在要在下面三個大類中選取100個x個體,輪盤賭轉100次以后,是不是個體數落在s3中的個體多一些,選擇的原理就是這樣,再不明白直接后面的程序吧。

7、還差點什么呢
至此,感覺也差不多了吧,選擇完后在重復上述步驟交叉,變異等等,那么什么時候是個頭了?很簡單,辦法就是迭代次數,迭代10次看一下結果,20次,看一下結果,30次,40次,100次,當次數達到一定程度以后,優秀的個體越來越多,大都集中在最優解附近,即使變異或者交叉了也是在這個最優解附近,沒有影響的。在下一次選擇就又變回來了,那么至此就真的結束了。比如說先來結果吧,該問題按我的思路做完后,迭代100次變成什么樣子了?看圖:

下面看代碼:
(1)首先看主函數
function main() clear; clc; %種群大小 popsize=100; %二進制編碼長度 chromlength=10; %交叉概率 pc = 0.6; %變異概率 pm = 0.001; %初始種群 pop = initpop(popsize,chromlength); for i = 1:100 %計算適應度值(函數值) objvalue = cal_objvalue(pop); fitvalue = objvalue; %選擇操作 newpop = selection(pop,fitvalue); %交叉操作 newpop = crossover(newpop,pc); %變異操作 newpop = mutation(newpop,pm); %更新種群 pop = newpop; %尋找最優解 [bestindividual,bestfit] = best(pop,fitvalue); x2 = binary2decimal(bestindividual); x1 = binary2decimal(newpop); y1 = cal_objvalue(newpop); if mod(i,10) == 0 figure; fplot('10*sin(5*x)+7*abs(x-5)+10',[0 10]); hold on; plot(x1,y1,'*'); title(['迭代次數為n=' num2str(i)]); %plot(x1,y1,'*'); end end fprintf('The best X is --->>%5.2f\n',x2); fprintf('The best Y is --->>%5.2f\n',bestfit);
(2)下面看二進制種群生成的方法
%初始化種群大小 %輸入變量: %popsize:種群大小 %chromlength:染色體長度-->>轉化的二進制長度 %輸出變量: %pop:種群 function pop=initpop(popsize,chromlength) pop = round(rand(popsize,chromlength)); %rand(3,4)生成3行4列的0-1之間的隨機數 % rand(3,4) % % ans = % % 0.8147 0.9134 0.2785 0.9649 % 0.9058 0.6324 0.5469 0.1576 % 0.1270 0.0975 0.9575 0.9706 %round就是四舍五入 % round(rand(3,4))= % 1 1 0 1 % 1 1 1 0 % 0 0 1 1 %所以返回的種群就是每行是一個個體,列數是染色體長度
(3)下面看如何把二進制返回對應的十進制
%二進制轉化成十進制函數 %輸入變量: %二進制種群 %輸出變量 %十進制數值 function pop2 = binary2decimal(pop) [px,py]=size(pop); for i = 1:py pop1(:,i) = 2.^(py-i).*pop(:,i); end %sum(.,2)對行求和,得到列向量 temp = sum(pop1,2); pop2 = temp*10/1023;
輸入的是100組0,1編碼的二進制,輸出的是x值,開始取一下種群大小,size(pop),顯然這里py是10了,借着對每一位求和,就是pop1(:,i)=2.^(py-i).*pop(:,i);這里省略用了冒號,,什么依稀呢?就是對所有行都有這個操作,冒號意思就是胸1到100了,那么就其中一個個體來說吧,假設為11001 10010,那么先進性這么一個操作就是什么呢?是不是就是對應的為0或1乘以2的對應次冪,如果1就是管用,是0就不管用。那么這個值就是(2^0)*1+(2^1)*1+0+0+(2^4)*1+....這樣就算出了一個值,因為是10位編碼,所以這個數是結余0-2^9即0-1023.那么最大為多少?1023吧。temp = sum(pop1,2)是對行求和吧,2表示行,1表示列,最后一行是吧它轉化為100組0-10之間的數值了。
(4)下面看計算適應度函數:
%計算函數目標值 %輸入變量:二進制數值 %輸出變量:目標函數值 function [objvalue] = cal_objvalue(pop) x = binary2decimal(pop); %轉化二進制數為x變量的變化域范圍的數值 objvalue=10*sin(5*x)+7*abs(x-5)+10;
(5)如何選擇新的個體
上面所有個體的函數值都計算出來了,存在objvalue中,此時它是不是也是100組y值啊,恩,那么對於現有的隨機生成的100組x,怎么來再選擇100組新的更好的x呢?這里我們把選擇放在了交叉與變異之間了,都可以,如何選擇,就要構造概率的那個輪盤了,誰的概率大,是不是選擇的個體就會多一些?也就是現在的選擇就是100中100個,最后出現的就夠就是以前的100個中最優的x有一個的話,選擇完后,可能就變成5個這個x了,多余的4個是不是相當於頂替了以前的不好的4個x值,這樣才能達到x總數100不變啊。
%如何選擇新的個體 %輸入變量:pop二進制種群,fitvalue:適應度值 %輸出變量:newpop選擇以后的二進制種群 function [newpop] = selection(pop,fitvalue) %構造輪盤 [px,py] = size(pop); totalfit = sum(fitvalue); p_fitvalue = fitvalue/totalfit; p_fitvalue = cumsum(p_fitvalue);%概率求和排序 ms = sort(rand(px,1));%從小到大排列 fitin = 1; newin = 1; while newin<=px if(ms(newin))<p_fitvalue(fitin) newpop(newin,:)=pop(fitin,:); newin = newin+1; else fitin=fitin+1; end end
(6)怎么交叉
%交叉變換 %輸入變量:pop:二進制的父代種群數,pc:交叉的概率 %輸出變量:newpop:交叉后的種群數 function [newpop] = crossover(pop,pc) [px,py] = size(pop); newpop = ones(size(pop)); for i = 1:2:px-1 if(rand<pc) cpoint = round(rand*py); newpop(i,:) = [pop(i,1:cpoint),pop(i+1,cpoint+1:py)]; newpop(i+1,:) = [pop(i+1,1:cpoint),pop(i,cpoint+1:py)]; else newpop(i,:) = pop(i,:); newpop(i+1,:) = pop(i+1,:); end end
(7)怎么變異
%關於編譯 %函數說明 %輸入變量:pop:二進制種群,pm:變異概率 %輸出變量:newpop變異以后的種群 function [newpop] = mutation(pop,pm) [px,py] = size(pop); newpop = ones(size(pop)); for i = 1:px if(rand<pm) mpoint = round(rand*py); if mpoint <= 0; mpoint = 1; end newpop(i,:) = pop(i,:); if newpop(i,mpoint) == 0 newpop(i,mpoint) = 1; else newpop(i,mpoint) == 1 newpop(i,mpoint) = 0; end else newpop(i,:) = pop(i,:); end end
(7)選擇最優個體
%求最優適應度函數 %輸入變量:pop:種群,fitvalue:種群適應度 %輸出變量:bestindividual:最佳個體,bestfit:最佳適應度值 function [bestindividual bestfit] = best(pop,fitvalue) [px,py] = size(pop); bestindividual = pop(1,:); bestfit = fitvalue(1); for i = 2:px if fitvalue(i)>bestfit bestindividual = pop(i,:); bestfit = fitvalue(i); end end
到這里遺傳算法程序就算是結束了。
看部分圖片結果:



