cd(fileparts(which(mfilename)))
addpath(genpath('./../ICB-DCM-AMICI-35f18fe'))
addpath(genpath('./../ICB-DCM-PESTO-649f3d2'))

amiwrap('model_ErbB','model_ErbB_syms',pwd);

rng(0);

% drug concentrations
conc = [2.5,8,25,80,250,800,2500,8000];
% {'Erlotinib','Lapatinib','PD0325901','PLX_{4720}','CHIR_{265}','Selumetinib','Vandetanib'};
% number of drugs
drug = 6;
% number of cell-lines
ncl = 10;

nexp = ncl*(1+length(conc)*length(drug));

load('Clim')

C_syn = zeros(length(Cmin),nexp);
for ic = 1:size(length(Cmin),1)
    C_syn(ic,:) = Cmin(ic) + (Cmax(ic)-Cmin(ic))*rand(1,nexp);
end

% erlotinib 7
% Selumetinib 11
% PLX-4720 12
% lapatinib 13
% PD0325901 14
% Sorafenib 15
% Vandetanib 16
% chir-265 17

drugidx = [7,13,14,12,17,11,16];
C_syn(drugidx,:) = 0;

np = 4104;

theta = rand(np,1)*4-2;

% generate artificial data
for icl = 1:ncl
    % ref
    I_syn{icl}.ref = (icl-1)*(1+length(conc)*length(drug))+1;
    I_syn{icl}.name = ['TUMOR-Syn' num2str(icl)];
    sol = simulate_model_ErbB(10^12,theta,C_syn(:, I_syn{icl}.ref),[],[]);
    prol_ref = sol.y;
    iexp = 1;
    for idrug = 1:length(drug)
        for iconc = 1:length(conc)
            I_syn{icl}.conc(iconc) = conc(iconc);
            I_syn{icl}.drug(idrug) = drug(idrug); 
            I_syn{icl}.exp(iexp) = (icl-1)*(1+length(conc)*length(drug))+1+iexp;
            C_syn(drugidx(idrug),I_syn{icl}.exp(iexp)) = I_syn{icl}.conc(iconc);
            sol = simulate_model_ErbB(10^12,theta,C_syn(:, I_syn{icl}.exp(iexp)),[],[]);
            I_syn{icl}.D(iexp).Y = sol.y/prol_ref;
            I_syn{icl}.D(iexp).Sigma_Y = 1;
            iexp = iexp+1;
        end
    end
end

parameters.min = -2*ones(np,1);
parameters.max = 2*ones(np,1);
parameters.number = np;
for ip = 1:parameters.number
    parameters.name{ip} = ['p' num2str(ip)];
end


options_MS = PestoOptions();
options_MS.localOptimizerOptions = optimset('algorithm','interior-point',...
                                   'GradObj','on',...
                                   'MaxIter',100,...
                                   'TolFun',0,...
                                   'TolX',1e-8,...
                                   'display','iter',... %        'display','final',...
                                   'MaxFunEvals',1000*parameters.number);
options_MS.comp_type = 'sequential';
options_MS.mode = 'text';
options_MS.save = true;
options_MS.trace = true;
options_MS.tempsave = true;
options_MS.foldername = 'ErbB_synthetic_example';

parameters = getMultiStarts(parameters,@(theta) likelihood_all_ErbB(10^9,theta,C_syn,I_syn,0),options_MS);


