function [m,b] = ag_recta(x,y)
%AG_RECTA Straight-line fitting genetic algorithm.
% [M,B] = AG_RECTA(X,Y) finds the the slope M and the y-intercept B of
% the line that better fits to the X,Y data. X and Y are vectors of the
% same length.
%% Ajusta una linea recta usando algoritmos genticos
% Encuentra los parmetros |m| y |b| de la ecuacin de la linea recta
% |Y=m*X+b| que mejor se ajusta a los datos (|x|,|y|) mediante un sencillo
% Algoritmo Gentico (AG). Se compara el resultado obtenido con AG contra
% el obtenido con Mnimos Cuadrados.
% Versin 1.7.1
%% Parmetros iniciales
% Se configua el Algoritmo definiendo parmetros anlogos en la Gentica
generaciones = 800; % iteraciones del algoritmo
n = 200; % tamao de la poblacin de soluciones
probabilidad = 0.6; % probabilidad de que cada gen mute
mutacion = [1 1]; % mximo valor que puede mutar
%% Poblacin inicial
% Se genera aleatoriamente una poblacin inicial de posibles soluciones.
% La primer columna representa el gen que define la pendiente |m| de la
% recta. La segunda columna representa el gen que define la ordenada al
% origen |b|. Cada uno de los n individuos (renglones) representa una
% posible solucin.
poblacion = rand(n,2);
%% Iteraciones
% La poblacin de soluciones evoluciona mediante la cruza y mutacin
s = n;
hw = waitbar(0,'Evolucionando...');
for i = 1:generaciones
%% Orden aleatorio de parejas
poblacion(:,3) = rand(s,1);
poblacion = sortrows(poblacion,+3);
%% Cruza
poblacion=[poblacion; poblacion(1:2:end-1,1) poblacion(2:2:end,2) zeros(fix(s/2),1); poblacion(2:2:end,1) poblacion(1:2:end-1,2) zeros(fix(s/2),1)];
%% Mutacion
s = size(poblacion,1);
poblacion=[poblacion; poblacion+[(rand(s,3)<probabilidad).*(rand(s,3)*2-1).*repmat([mutacion 0],s,1)]];
%% Clculo de la aptitud (o viabilidad)
% La tercer columna de poblacion representa la aptitud (viabilidad) de
% la solucin, menor es mejor.
s = size(poblacion,1);
m = repmat(poblacion(:,1),1,100);
b = repmat(poblacion(:,2),1,100);
X = repmat(x,s,1);
Y = repmat(y,s,1);
Y2=m.*X+b;
poblacion(:,3)=sum(abs(Y-Y2),2);
%% Seleccion
poblacion = unique(poblacion,'rows');
poblacion = sortrows(poblacion,+3);
s = min([size(poblacion,1),n]);
poblacion = poblacion(1:s,:);
waitbar(i/generaciones,hw)
end
close(hw)
%% Resultados
% Se compara el resultado obtenido con AG contra el obtenido con mnimos
% cuadrados
m = poblacion(1,1);
b = poblacion(1,2);
error_ag = poblacion(1,3);
p = polyfit(x,y,1);
y3 = polyval(p,x);
error_mc = sum(abs(y3-y));
plot(x,y,'.',x,polyval(poblacion(1,1:2),x),x,y3,':')
legend('Datos Sintticos','Algoritmos Genticos (AG)','Mnimos Cuadrados (MC)')
disp(['AG: Error = ' num2str(error_ag) ', m = ' num2str(m) ', b = ' num2str(b)])
disp(['MC: Error = ' num2str(error_mc) ', m = ' num2str(p(1)) ', b = ' num2str(p(2))])