In reply to Ksharp.

Nice.
I've managed to run KSHARP's GA code with 2 variables. The code is attached below.

I think it's cleaner to apply the ## operator (elementwise power operator) to the objective function to be minimized (SSE in KSHARP's code).
I wonder which command in the GA setup forces the group split to be relatively equal in size of the number of observations. Thanks in advance for the answer.

proc iml;
use units nobs nobs;
read all var {Cholesterol AgeAtStart};
close;

start function(x) global(Cholesterol, AgeAtStart, nobs, group, mean_std);
if countunique(x)=group then do;
do i=1 to group;
idx=loc(x=i);
temp=AgeAtStart[idx];
mean_std[i,1]=mean(temp);
mean_std[i,2]=std(temp);

temp1=Cholesterol[idx];
mean_std[i,3]=mean(temp1);
mean_std[i,4]=std(temp1);

end;
sse=sum((mean_std[,]-mean_std[><,])##2) ;
end;
else sse=999999;
return (sse);
finish;

group=5; /* <--Change it(divide into 5 groups)*/
mean_std=repeat(j(group,2,.),1,2);
/*print mean_std;*/

encoding=j(2,nobs,1);

encoding[2,]=group;

id=gasetup(2,nobs,123456789);
call gasetobj(id,0,"function");
call gasetsel(id,10,1,1);
call gainit(id,100,encoding);

niter = 50;
do i = 1 to niter;
call garegen(id);
call gagetval(value, id);
end;
call gagetmem(mem, value, id, 1);

col_mem=t(mem);
create group var {col_mem};
append;
close;

print value[l = "Min Value:"] ;
call gaend(id);
quit;

data want;
merge group units(keep=AgeAtStart Cholesterol);
run;
proc means data=want mean std;
class col_mem;
var AgeAtStart Cholesterol;
run;