%AK Specify parameters in a vector that can be read in separately to avoid
%updating code every optimization
%Can enter this vector differently if pulling in from other files, etc.

clock = clock

% %First, load in the pmin from
load('20151027_ps_20151021_Revised10_PMCA_PeakHeight_10000_1','pminp');
%Results_psearch_20150330_All3_noMito_10000_code.mat
% %Second, set p = pminp 
p = pminp;

%Bounds as specified by code slight alteration to p(16) and p(20) to see if
%the drop out of the mito can be reduced
lowbounds = [p(1),0.2,p(3),0.01,p(5),0.0005,p(7),p(8),p(9),p(10),...      %10 values
    0.001,p(12),p(13),p(14),p(15),0.01,p(17),p(18),p(19),50,...                    %10 values
    p(21),p(22),p(23),5,1,50,0.001,p(28),0.5,p(30),...                          %10 values
    p(31) 280];                            
upbounds = [p(1),250,p(3),50,p(5),0.05,p(7),p(8),p(9),p(10),...                %10 values
    1,p(12),p(13),p(14),p(15),0.1,p(17),p(18),p(19),500,...                              %10 values
    p(21),p(22),p(23),450,10,1000,1,p(28),1.5,p(30),...                               %10 values
    p(31) 420];  

%specify the size of the population
popsize = 1400;
%Develop the initial population
M = zeros(popsize,length(lowbounds));
r = rand(popsize,length(lowbounds));
for i = 1:popsize
    for j = 1:length(lowbounds)-1
        if lowbounds(j) == upbounds(j)
            M(i,j) = lowbounds(j);
        else
        M(i,j) = 0.8*p(j)+0.4.*r(i,j).*(p(j)); 
        end
    end
    for j = length(lowbounds)
        %The original value was 350, but this wasn't saved in p the same
        %way so we are coding it outright
        M(i,j) = 0.8.*350+0.4.*350.*r(i,j);
    end
end


%Need to perform minimization of parameter sets

%*******Sets options for running the genetic algorithm
%--------------------------------------------------------------------------
%Added output into the output of ga on March 2, 2015
%PopulationSize specifies the size of the population
%IniialPopulation is used to seed the genetic algorithm
%PlotFcns is an array of handles to functions that plot data computed by
% % the algorithm
options = gaoptimset('PopulationSize',popsize,'Generations',100, ...
    'InitialPopulation',M,'PlotFcns', ...
    {@gaplotbestf,@gaplotscores,@gaplotstopping,@gaplotdistance});

[pmin,Smin,exitflag,output] = ga(@AK_Optimization_ToMinimize_20151119_TCell_DiffS,...
    length(lowbounds),[],[],[],[],lowbounds,upbounds,[],options);


%*******Sets options for running the fmincon
%--------------------------------------------------------------------------
%This enables plotting of the results
% options = optimset('MaxFunEvals',6000,'PlotFcns',{@optimplotfval,@optimplotstepsize});
% % % Try running fmincon: find minimum of constrained nonlinear multivariable
% % % function
% % 
% % %Here, the pminp is loaded in from
% %Results_psearch_20150330_All3_noMito_10000_code.mat
% 
% %Load in AAK_GA_20150715_Young_2000_1.mat for
% %the pmin value to load into fmincon
% 
% [pminf,Smin,exitflag,output]=fmincon(@AK_Optimization_ToMinimize_20150708_TCell_DiffS,...
%     pmin,[],[],[],[],lowbounds,upbounds,[],options);



%********Sets options for running the patternsearch
%--------------------------------------------------------------------------
%This enables plotting of the results
% options = psoptimset('PlotFcns',{@psplotbestf,@psplotmeshsize});
% 
% %loaded in AK_fmincon_20150412_TCell_oldbounds_10000 for pminf
% 
% [pminp,Smin,exitflag,output] = ...
%     patternsearch(@AK_Optimization_ToMinimize_20150708_TCell_DiffS, ...
%     pminf,[],[],[],[],lowbounds,upbounds,[],options);


%  
%*******Sets options for running the simulannealbnd
%--------------------------------------------------------------------------
%Find minimum of function using simulated annealing algorithm
%This enables plotting of the results
% options = saoptimset('PlotFcns',...
%     {@saplotbestf,@saplotbestx,@saplotf,@saplottemperature});
% 
% [pmin,Smin,exitflag,output] = simulannealbnd(@AK_Optimization_ToMinimize_20150331_TCell,...
%     pmin,lowbounds,upbounds,options);