## How can i find the optimal parameters that minimize SSE?

### Riccardo Bartoloni (view profile)

on 30 Apr 2019
Latest activity Answered by Alan Weiss

### Alan Weiss (view profile)

on 2 May 2019
Hi, i am working for my university on this code and i have to minimize the last line, the variable "SQ_TOT". The solution of this exercise is the vector "Param" (that give me minimum SSE --> SQ_TOT= 3.1082e+05) but i don't know how to find this vector. How can i find the optimal parameters that minimize SQ_TOT? Thank you all!
clc
clear
t=0.33333:0.33333: 76.33333;
Param=[0.5254010 11.4402640 29861.1830810 2.0010250 -0.051606000 -0.057373000 0.036517000 -0.099394000 1.647586000 81.195773000 1.422753000 -0.078617000 0.036363000 0.106596000 0.128553000]
%Param=ones(1,15);
DUMqA=[Param(5:8)];
expbeta=exp(DUMqA*DUM);
t1=t.^Param(4);
t2=t.^Param(4);
t2(length(t))=[];
t2=[0 t2];
tc=t1-t2;
lt=length(t1);
cohort=zeros(1,lt);
for i=1:length(t1);
summ=[zeros(1,i-1) expbeta(i:length(t1)).*tc(1:length(t)+1-i)];
cohort=[cohort;summ];
end;
cohort(1,:)=[];
B_A=[];
HEY_A=[];
F_A=[];
for k=1:lt;
for i=1:lt;
summ=sum(cohort(k,1:(lt+1-i)));
B_A=[summ B_A];
end;
HEY_A=[HEY_A;B_A];
B_A=[];
end;
F_A=[];
for k=1:lt;
for i=1:lt;
summ=(1-(Param(3)/(Param(3)+HEY_A(k,i)))^Param(2))*(1-Param(1));
B_A=[B_A summ];
end;
F_A=[F_A; B_A];
B_A=[];
end;
DUMqR=[Param(12:15)];
expbeta_R=exp(DUMqR*DUM);
t1R=t.^Param(11);
t2R=t.^Param(11);
t2R(length(t))=[];
t2R=[0 t2R];
tc_R=t1R-t2R;
lt=length(t1);
cohort_R=[];
for i=1:length(t1);
summ=[zeros(1,i-1) expbeta_R(i:length(t1)).*tc_R(1:length(t)+1-i)];
cohort_R=[cohort_R;summ];
end;
cohort_R(1,:)=[];
B_R=[];
HEY_R=[];
for k=1:lt-1;
for i=1:lt;
summ=sum(cohort_R(k,1:(lt+1-i)));
B_R=[summ B_R];
end;
HEY_R=[HEY_R;B_R];
B_R=[];
end;
F_R=[];
for k=1:lt-1;
for i=1:lt;
summ=(1-(Param(10)/(Param(10)+HEY_R(k,i)))^Param(9));
B_R=[B_R summ];
end;
F_R=[F_R; B_R];
B_R=[];
end;
HouseCha1=housed;
housed(229)=[];
HouseCha2=[0 housed];
ChaninPop=HouseCha1-HouseCha2;
F_AZ=[zeros(lt,1) F_A];
EXP_A=[];
EXP_AA=[];
EC_R=zeros(1,lt);
EC=[];
F_RZ=[zeros(228,1) F_R];
EXP_R=[];
EXP_RR=[];
EC_RVERO=[];
EC_A=[zeros(lt,1)];
EXP_AA=[];
EXP_A=[];
EXP_RR=[];
EC=[];
EC_R=[zeros(lt,1)];
EC_AB=[0];
EC_TRAN=[];
EC_RTRAN=[];
for i=1:lt;
for k=1:lt-1;
EXP_AATRAN=(F_A(k,i)-F_AZ(k,i))*(ChaninPop(k)+EC_R(k));
EXP_A=[EXP_A EXP_AATRAN];
ECTRAN=sum(EXP_A(:));
EXP_RRTRAN=(F_R(k,i)-F_RZ(k,i))*EC_A(k);
EXP_R=[EXP_R EXP_RRTRAN];
EC_RVEROTRAN=sum(EXP_R(:));
end;
EXP_AA=[EXP_AA (EXP_A)'];
EXP_A=[];
EC_TRAN=[EC_TRAN ECTRAN];
EC_A=[EC_TRAN zeros(1,lt)];
EC=[EC ECTRAN];
EXP_RR=[EXP_RR (EXP_R)'];
EXP_R=[];
EC_RTRAN=[EC_RTRAN EC_RVEROTRAN];
EC_R=[0 EC_RTRAN zeros(1,lt)];
EC_RVERO=[EC_RVERO EC_RVEROTRAN];
end
EC3m_R=[0 0];
for i=3:lt;
summ=sum(EC_RVERO(i-2:i));
EC3m_R=[EC3m_R summ];
end;
EC3m=[EC(1) sum(EC(1:2))];
for i=3:lt;
summ=sum(EC(i-2:i));
EC3m=[EC3m summ];
end;
DIVEC3=DIV_.*EC3m;
DIV_A2=(DIV_A).^2;
DIV=LOSS~=0;
DIV_L=DIV.*EC3m_R;
DIV_LOSS=(DIV_L-LOSS)
DIV_LOSS2=(DIV_LOSS).^2
SQ_LOSS=sum(DIV_LOSS2(1,:));
EXP_E=EC-EC_RVERO;
EXP_END=[];
for i=1:lt;
summ=sum(EXP_E(1:i));
EXP_END=[EXP_END summ];
end;
DIV1=END~=0;
DIV2=DIV1-DIV_;
DIV_E=DIV2.*EXP_END;
END_DIV=DIV2.*END;
DIV_E(1)=EXP_END(1)
SQ_EN=[(EXP_END(1))^2];
pos=find(DIV2==1);
for i=4:3:pos(length(pos));
DIV_END=((DIV_E(i)-DIV_E(i-3))-(END_DIV(i)-END_DIV(i-3)))^2;
SQ_EN=[SQ_EN DIV_END];
end;
SQ_END=sum(SQ_EN(1,:));
SQ_END=SQ_END+(END(length(END))-EXP_END(length(EXP_END)))^2;

### Alan Weiss (view profile)

on 1 May 2019

To solve an optimization problem you have to force your problem into the form required by optimization solvers. Sorry, that's just the way it is. You have to decide first which variables you can optimize over. These are usually called "decision variables", and are typically called x = [x(1),x(2),...,x(n)], where n is the number of decision variables.
Then you write your objective function, the thing to minimize, as a function of x. See, for example, Optimizing Nonlinear Functions.
Then you guess a starting value of x, say x0, and call an optimization solver, such as
sol = fminsearch(@objective,x0)
If you have constraints on your decision variables, such as they must all be positive, then you need an Optimization Toolbox™ license, and use an appropriate solver.
Good luck,
Alan Weiss
MATLAB mathematical toolbox documentation

Riccardo Bartoloni

### Riccardo Bartoloni (view profile)

on 2 May 2019
Have i to do something like that? But matlab say "Not enough input arguments. Error in Myobjectivefunction (line 6)
DUMqA=[Param(5:8)];"
function [SQ_TOT] = Myobjectivefunction(Param)
t=0.33333333333333333333333333333333333333333333333333333333333333333333333333333333333333:0.33333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333:76.33333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333;
%Param=[0.5254010 11.4402640 29861.1830810 2.0010250 -0.051606000 -0.057373000 0.036517000 -0.099394000 1.647586000 81.195773000 1.422753000 -0.078617000 0.036363000 0.106596000 0.128553000];
%Param=Param(1:15);
DUMqA=[Param(5:8)];
expbeta=exp(DUMqA*DUM);
t1=t.^Param(4);
t2=t.^Param(4);
t2(length(t))=[];
t2=[0 t2];
tc=t1-t2;
lt=length(t);
cohort=zeros(1,lt);
for i=1:length(t1);
summ=[zeros(1,i-1) expbeta(i:lt).*tc(1:lt+1-i)];
cohort=[cohort;summ];
end;
cohort(1,:)=[];
HEY_A=cumsum(cohort,2);
F_A=(1-(Param(3)./(Param(3)+HEY_A)).^Param(2)).*(1-Param(1));
%%%%Retention Model%%%%%
DUMqR=[Param(12:15)];
expbeta_R=exp(DUMqR*DUM);
t1R=t.^Param(11);
t2R=t.^Param(11);
t2R(length(t))=[];
t2R=[0 t2R];
tc_R=t1R-t2R;
cohort_R=[];
%cohort_R= matrice riga 1000
for i=1:length(t1);
summ=[zeros(1,i-1) expbeta_R(i:length(t1)).*tc_R(1:length(t)+1-i)];
cohort_R=[cohort_R;summ];
end;
cohort_R(1,:)=[];
HEY_R=cumsum(cohort_R,2);
F_R=(1-(Param(10)./(Param(10)+HEY_R)).^Param(9));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
HouseCha1=housed;
%HouseCha1= riga 19
housed(229)=[];
HouseCha2=[0 housed];
ChaninPop=HouseCha1-HouseCha2;
F_AZ=[zeros(lt,1) F_A];
EXP_A=[];
EXP_AA=[];
EC_RVERO=zeros(1,lt);
EC=[];
%EC= riga 976 oppure colonna IA e riga 1712
F_RZ=[zeros(228,1) F_R];
F_AZ(:,end)=[];
F_RZ(:,end)=[];
ChaninPopR=ChaninPop(1:end-1);
EC_RVERO=[0];
for i=1:lt
EXP_AA=(F_A-F_AZ).*(ChaninPop+EC_RVERO)';
EC=sum(EXP_AA);
EC(lt)=[];
EXP_RR=(F_R-F_RZ).*(EC)';
EC_RVERO=sum(EXP_RR);
EC_RVERO=[0 EC_RVERO];
EC_RVERO(end)=[];
end
EC=[EC sum(EXP_AA(:,end))];
EC_RVERO(1)=[];
EC_RVERO=[EC_RVERO sum(EXP_RR(:,end))];
EC3m_R=[0 0];
for i=3:lt;
summ=sum(EC_RVERO(i-2:i));
EC3m_R=[EC3m_R summ];
end;
EC3m=[EC(1) sum(EC(1:2))];
for i=3:lt;
summ=sum(EC(i-2:i));
EC3m=[EC3m summ];
end;
DIVEC3=DIV_.*EC3m;
DIV_A2=(DIV_A).^2;
DIV=LOSS~=0;
DIV_L=DIV.*EC3m_R;
DIV_LOSS=(DIV_L-LOSS);
DIV_LOSS2=(DIV_LOSS).^2
SQ_LOSS=sum(DIV_LOSS2(1,:));
EXP_E=EC-EC_RVERO;
EXP_END=cumsum(EXP_E);
DIV1=END~=0;
DIV2=DIV1-DIV_;
DIV_E=DIV2.*EXP_END;
END_DIV=DIV2.*END;
DIV_E(1)=EXP_END(1)
SQ_EN=[(EXP_END(1))^2];
pos=find(DIV2==1);
for i=4:3:pos(length(pos));
DIV_END=((DIV_E(i)-DIV_E(i-3))-(END_DIV(i)-END_DIV(i-3)))^2;
SQ_EN=[SQ_EN DIV_END];
end;
SQ_END=sum(SQ_EN(1,:));
SQ_END=SQ_END+(END(length(END))-EXP_END(length(EXP_END)))^2;
end
Riccardo Bartoloni

### Riccardo Bartoloni (view profile)

on 2 May 2019
Sorry now it work thank you very much to you all!!
Riccardo Bartoloni

### Riccardo Bartoloni (view profile)

on 2 May 2019
My last problem is to find the minimum. "fminsearch" arrive at a minimum of 500000 but the real minimum is more or less 180000. How can i find this minimum? This is my code. I put the data in only one file as you said.
function [SQ_TOT] = Myobjectivefunction(Param)
t=0.33333333333333333333333333333333333333333333333333333333333333333333333333333333333333:0.33333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333:76.33333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333;
%Param=[0.5254010 11.4402640 29861.1830810 2.0010250 -0.051606000 -0.057373000 0.036517000 -0.099394000 1.647586000 81.195773000 1.422753000 -0.078617000 0.036363000 0.106596000 0.128553000];
%Param=Param(1:15);
DUMqA=[Param(5:8)];
expbeta=exp(DUMqA*DUM);
t1=t.^Param(4);
t2=t.^Param(4);
t2(length(t))=[];
t2=[0 t2];
tc=t1-t2;
lt=length(t);
cohort=zeros(1,lt);
for i=1:length(t1);
summ=[zeros(1,i-1) expbeta(i:lt).*tc(1:lt+1-i)];
cohort=[cohort;summ];
end;
cohort(1,:)=[];
HEY_A=cumsum(cohort,2);
F_A=(1-(Param(3)./(Param(3)+HEY_A)).^Param(2)).*(1-Param(1));
%%%%Retention Model%%%%%
DUMqR=[Param(12:15)];
expbeta_R=exp(DUMqR*DUM);
t1R=t.^Param(11);
t2R=t.^Param(11);
t2R(length(t))=[];
t2R=[0 t2R];
tc_R=t1R-t2R;
cohort_R=[];
%cohort_R= matrice riga 1000
for i=1:length(t1);
summ=[zeros(1,i-1) expbeta_R(i:length(t1)).*tc_R(1:length(t)+1-i)];
cohort_R=[cohort_R;summ];
end;
cohort_R(1,:)=[];
HEY_R=cumsum(cohort_R,2);
F_R=(1-(Param(10)./(Param(10)+HEY_R)).^Param(9));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
HouseCha1=housed;
%HouseCha1= riga 19
housed(229)=[];
HouseCha2=[0 housed];
ChaninPop=HouseCha1-HouseCha2;
F_AZ=[zeros(lt,1) F_A];
EXP_A=[];
EXP_AA=[];
EC_RVERO=zeros(1,lt);
EC=[];
%EC= riga 976 oppure colonna IA e riga 1712
F_RZ=[zeros(228,1) F_R];
F_AZ(:,end)=[];
F_RZ(:,end)=[];
ChaninPopR=ChaninPop(1:end-1);
EC_RVERO=[0];
for i=1:lt
EXP_AA=(F_A-F_AZ).*(ChaninPop+EC_RVERO)';
EC=sum(EXP_AA);
EC(lt)=[];
EXP_RR=(F_R-F_RZ).*(EC)';
EC_RVERO=sum(EXP_RR);
EC_RVERO=[0 EC_RVERO];
EC_RVERO(end)=[];
end
EC=[EC sum(EXP_AA(:,end))];
EC_RVERO(1)=[];
EC_RVERO=[EC_RVERO sum(EXP_RR(:,end))];
EC3m_R=[0 0];
for i=3:lt;
summ=sum(EC_RVERO(i-2:i));
EC3m_R=[EC3m_R summ];
end;
EC3m=[EC(1) sum(EC(1:2))];
for i=3:lt;
summ=sum(EC(i-2:i));
EC3m=[EC3m summ];
end;
DIVEC3=DIV_.*EC3m;
DIV_A2=(DIV_A).^2;
DIV=LOSS~=0;
DIV_L=DIV.*EC3m_R;
DIV_LOSS=(DIV_L-LOSS);
DIV_LOSS2=(DIV_LOSS).^2
SQ_LOSS=sum(DIV_LOSS2(1,:));
EXP_E=EC-EC_RVERO;
EXP_END=cumsum(EXP_E);
DIV1=END~=0;
DIV2=DIV1-DIV_;
DIV_E=DIV2.*EXP_END;
END_DIV=DIV2.*END;
DIV_E(1)=EXP_END(1)
SQ_EN=[(EXP_END(1))^2];
pos=find(DIV2==1);
for i=4:3:pos(length(pos));
DIV_END=((DIV_E(i)-DIV_E(i-3))-(END_DIV(i)-END_DIV(i-3)))^2;
SQ_EN=[SQ_EN DIV_END];
end;
SQ_END=sum(SQ_EN(1,:));
SQ_END=SQ_END+(END(length(END))-EXP_END(length(EXP_END)))^2;
end

### Alan Weiss (view profile)

on 2 May 2019

I think that you did not understand the comment about not loading the data within your objective function. You should load it once before you call the objective, and just pass in the values, like this: