simgeneticalgorithm.m

X(t)=Asin(2*pi *f *t+ q)+n(t) 估计其中的参数为A

% 基于遗传算法的单正弦信号参数估计
% 函数模型
% x(t)=A*sin(2*pi*f*t+q)+noi(t)
% 估计参数 A,f,q 
%在这里 A=2, f=2.5, q=1.5

% 该程序中的目标函数
%  F(A ,f,q ) = 叠加abs(Y(i)- A*sin(2*pi*f*i*T+q))  i=1...SimpleSize
%  0 <=  x1 , x2  <= XRang    


clear all
close all
%--------------------------------
%函数参数设定
T = 0.1; %采样周期
pi = 3.14159;
SimpleSize = 100; %采样数
%--------------------------------

% Some parameters of the Simple Genetic Algorithm
PopulationSize  = 40 ;   %The population of the colony
MaxGeneration = 400 ; %The maximal generation 
ProCro = 0.9 ;               %The probability of crossover
ProMut = 0.001 ;            %The probability of mutation 

XRang =2.5  ;

Vmax = 1000 ; %固定最大值
Vmin  = 0 ; %固定最小值
% FunctionMode >= 0 find the maximum of the function
% FunctionMode < 0   find the minimum of the function
FunctionMode = -1 ; % 设定该题为求最小值

ChromosomeLen1 = 10  ; % The lenght of the first chromosome
ChromosomeLen2 = 10  ; % The lenght of the second chromosome
ChromosomeLen3 = 10  ; % The lenght of the second chromosome
ChromosomeLen  = 10  ;

ChromosomeLength = ChromosomeLen1 + ChromosomeLen2 + ChromosomeLen3 ; % The total lenght of chromosome

X1 = zeros ( 1 ,  PopulationSize );%随机生成零向量，用0初始化
X2 = zeros ( 1 ,  PopulationSize );
X3 = zeros ( 1 ,  PopulationSize );

%采集数据-----------------------
Noise = randn(1,SimpleSize); %生成随机白噪声序列
Y = zeros(1,SimpleSize); %函数采样初值

for j = 1 : SimpleSize
Y(j) = 2*sin(2*pi*2.5*j*T+1.5)+Noise(j);
end % end of for
%-------------------------------

% Generate the first generation of the population
rand('state',sum(100*clock)); % Resets the generator to a different state each time？？？？

for j = 1 : PopulationSize
    for k = 1 : ChromosomeLength
        Temp = rand(1); % 随机生成1×1的0－1的数
        if Temp < 0.5
            Population( j ) . Chromosome( k ) = 0;
        else
            Population( j ) . Chromosome( k ) = 1;
        end % End of if
    end % End of for - k     
    
    %X1( j ) = DecodeChromosome(Population( j ) . Chromosome, 1 , ChromosomeLen );
    %X2( j ) = DecodeChromosome(Population( j ) . Chromosome, ChromosomeLen+1 , ChromosomeLen );
    %X3( j ) = DecodeChromosome(Population( j ) . Chromosome, 2*ChromosomeLen+1 , ChromosomeLen );
    
    Population( j ).Value = 0 ;
    Population( j ).Fitness = 0 ;
end % End of for - j

%X1 = 2 * XRang * X1/1023 - XRang; %解码公式 将二进制代码对应的十进制数转换为定义域（－XRange, XRange）内的一系列值
%X2 = 2 * XRang * X2/1023 - XRang;
%X3 = 2 * XRang * X2/1023 - XRang;

%figure
%plot( X1 , X2 ,X3, '*') , grid; % 初始群体

BestValue = zeros( 1 , MaxGeneration ) ;
AverageValue = zeros( 1 , MaxGeneration ) ;

%%%%%%%%%×××××××××××Genetic Algorithm
for Gen = 1 : MaxGeneration
    % Calculate the value according to the function which we want to find its minimum 
    % In this  program the function is
    %  F(A ,f,q ) = 叠加(Y(i)- A*sin(2*pi*f*i*T+q))^2  i=1...SimpleSize
    %  s.t          -XRang  <=  x1 , x2  <= XRang    
    for j = 1 : PopulationSize
        X1( j ) = DecodeChromosome(Population( j ) . Chromosome, 1 , ChromosomeLen ); % 对应A
        X2( j ) = DecodeChromosome(Population( j ) . Chromosome, ChromosomeLen+1 , ChromosomeLen ); % 对应f
        X3( j ) = DecodeChromosome(Population( j ) . Chromosome, 2*ChromosomeLen+1 , ChromosomeLen ); % 对应q
        
       A = XRang * X1( j )/1023 ; 
       f = XRang * X2( j )/1023 ;
       q = XRang * X3( j )/1023 ;
       
       sum = 0;
       for i = 1 : SimpleSize
         sum = sum + abs(Y(i)- A*sin(2*pi*f*i*T+q))  ;
       end
       TempValue = sum ;
       Population( j ).Value = TempValue;
       
    end % End of for - j
    
    %Calculate fitness value
    for j = 1 : PopulationSize
       TempValue = CalculateFitnessValue( Population( j ).Value ,  Vmax , Vmin , FunctionMode); 
       Population( j ).Fitness = TempValue;
    end % End of for - j    
    
    % Find out the best and worst individual of this generation
    BestValue( Gen ) = Population( 1 ).Value ;
    for j = 1 : PopulationSize
        if Population( j ).Value < BestValue( Gen )  
            BestValue( Gen ) = Population( j ).Value;
        end % End of if
    end % End of for - j 
    
    % Calculate the AverageValue
    Sum = 0 ;
    for j = 1 : PopulationSize
        Sum = Sum + Population( j ).Value;
    end % End of for - j    
    AverageValue( Gen ) = Sum / PopulationSize ;

    % Reproduce a chromosome by proportional selection
    Sum = 0 ;
    CFitness = zeros ( 1 , PopulationSize ); % The cumulative fitness value
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % Calculate cumulative fitness 
    for j = 1 : PopulationSize
       Sum = Sum + Population( j ).Fitness ;
    end % End of for - j       
    
    for j = 1 : PopulationSize
       CFitness( j ) = Population( j ).Fitness / Sum ; % 相对适应度
    end % End of for - j       
    for j = 2 : PopulationSize
       CFitness( j ) = CFitness( j - 1 ) + CFitness( j ) ; %?????1, 1+2, 1+2+3, 1+2+3+...
    end % End of for - j           
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%选择
    % Select operation
    for j = 1 : PopulationSize
       RandP = rand( 1 ) ;  % 随机生成0－1之间的数
       Index = 1;
       while RandP  > CFitness( Index )
           Index =  Index + 1 ;
       end % End of  while       
       NewPopulation( j )  = Population( Index ) ;        
    end % End of for - j           
    
    for j = 1 : PopulationSize
       Population( j )  = NewPopulation( j ) ;        
    end % End of for - j      
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%交叉
    % Crossover two chromosome
    IndexPair = zeros( 1, PopulationSize);
    for j = 1 : PopulationSize
       IndexPair( j ) = j;
    end % End of for - j   
    
    % 随机配对
    for j = 1 : PopulationSize
       Point = floor ( (PopulationSize - j - 1)  *  rand( 1 ) ); % B = floor(A) ,B取小于或等于A的整数
       Temp = IndexPair( j ) ;
       IndexPair( j ) = IndexPair( j + Point );
       IndexPair( j + Point ) = Temp;
    end % End of for - j       
    
    for j = 1 : 2 : ( PopulationSize - 1 )
        RandP = rand( 1 ) ;
        if RandP < ProCro
            
            Point = floor ( rand( 1 ) *  ChromosomeLen ) ; % 随机交叉点
            if Point == 0
                Point = 1;
            end
            
            for k = Point : (ChromosomeLen + Point - 1)
                Temp = Population ( IndexPair( j ) ) . Chromosome( k ) ;
                Population ( IndexPair( j ) ) . Chromosome( k )  = Population ( IndexPair( j + 1) ) . Chromosome( k ) ;
                Population ( IndexPair( j + 1) ) . Chromosome( k )  = Temp ;
            end % End of for - k
            
        end % End of if
    end % End of for - j
        
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%变异
    % Mutation of a chromosome    
    for j = 1 : PopulationSize
        for k = 1 : 2*ChromosomeLen
            RandP = rand( 1 ) ;            
            if RandP < ProMut
                
                if Population( j ) . Chromosome( k ) == 1
                    Population( j ) . Chromosome( k ) = 0 ; 
                else
                    Population( j ) . Chromosome( k ) = 1 ;
                end% End of if            
                    
            end% End of if            
        end % End of for - k        
    end % End of for - j   
    
end % End of for - Gen

%%%%%%%%×××××××××××end

%-------------------------------输出最终估计的各参数值
for j = 1 : PopulationSize
       X1( j ) = DecodeChromosome(Population( j ) . Chromosome, 1 , ChromosomeLen ); % 对应A
       X2( j ) = DecodeChromosome(Population( j ) . Chromosome, ChromosomeLen+1 , ChromosomeLen ); % 对应f
       X3( j ) = DecodeChromosome(Population( j ) . Chromosome, 2*ChromosomeLen+1 , ChromosomeLen ); % 对应q
        
       A = XRang * X1( j )/1023 ; 
       f = XRang * X2( j )/1023 ;
       q = XRang * X3( j )/1023 ;
       
       sum = 0;
       for i = 1 : SimpleSize
         sum = sum + abs(Y(i)- A*sin(2*pi*f*i*T+q))  ;
       end
       TempValue = sum ;
       Population( j ).Value = TempValue;
end % End of for - j
    
FinalBest= Population( 1 ).Value;
BestNO=1;
for j=1 : PopulationSize
    if FinalBest > Population( j ).Value
       FinalBest = Population( j ).Value;
       BestNO=j;
    end % End of if
end % end of for - j
%-----------------------------------

X1( BestNO ) = DecodeChromosome(Population( BestNO ) . Chromosome, 1 , ChromosomeLen );
X2(BestNO ) = DecodeChromosome(Population( BestNO ) . Chromosome, ChromosomeLen+1 , ChromosomeLen ); 
X3( BestNO) = DecodeChromosome(Population( BestNO ) . Chromosome, 2*ChromosomeLen+1 , ChromosomeLen );

A =  XRang * X1(BestNO)/1023 
f =  XRang * X2(BestNO)/1023
q =  XRang * X3(BestNO)/1023 

%----------------------------------输出最终估计的各参数值

%X1min = -2.5 ;
%X1max = 2.5 ;
%X2min = -3;
%X2max = 3;
%X3min = -2.5 ;
%X3max = 2.5 ;

%figure
%plot( X1 , X2 , X3 ,'*') , grid , axis([X1min X1max X2min X2max X3min X3max]); %最优群体

%figure
%[X , Y] = meshgrid(-XRang : 0.1: XRang);
%Z = 100 * ( X.^2 - Y ).^2 + ( 1 - X ).^2;
%mesh(X , Y , Z);

% AimValue = ones(1 , MaxGeneration) *  (100 * ( XRang^2 + XRang )^2 + ( 1 + XRang )^2) ; %最大值直线


Xmin = 0 ;
Xmax = MaxGeneration ;
Ymin =0;
Ymax =600;
XL = 1 : MaxGeneration;

figure
hold on 
plot( XL ,BestValue , 'r' ) ; 
plot( XL , AverageValue , 'b' ) ; 
%plot( AimValue , 'g' ) ; 
hold off
grid
axis([Xmin Xmax Ymin Ymax]);