MATLAB: Failure of do loop after first step when applying globalsearch to a nonlinear least squares problem

do loopfminconglobalsearchMATLAB

Dear All,
I am trying to generate bootstrap standard errors for a nonlinear least-squares optimization.
The optimization works for the first run, but my do loop for generating repeated runs of the optimization stops after the first run.
I am using global search from the Global Optimization toolbox.
The error message is:
Index in position 1 exceeds array bounds (must not exceed 1).
Looking for help to solve this do loop problem.
Details are below.
Best,
Srinivasan
My dependent variable y is an 18,026 by 1 vector
My independent variable matrix X is 18,026 rows by 14 columns
My parameter vector is x2 and is an 1 by 15 vector.
I use globalsearch and nonlinear squares on the following function to generate an estimate of x2:
function F2 = myfun2(x2,y,X)
F2 = sum((((y - (1 + x2(1)*(x2(2)/0.0025662))*X(:,1) - (x2(2)/x2(3))*(sqrt(1-x2(1)^2))*X(:,2) - x2(4)*X(:,3) - x2(5)*X(:,4) - x2(6)*X(:,5) - x2(7)*X(:,6) - x2(8)*X(:,7) - x2(9)*X(:,8) - x2(10)*X(:,9) - x2(11)*X(:,10) - x2(12)*X(:,11) - x2(13)*X(:,12) - x2(14)*X(:,13) - x2(15)*X(:,14)))).^2);
end
gsproblem = createOptimProblem('fmincon','objective',@(x2) myfun2(x2, y, X),'x0',x0,'lb',lb,'ub',ub);
gs = GlobalSearch('Display','iter');
gs = GlobalSearch(gs,'XTolerance',1e-3,'StartPointsToRun','all', 'NumTrialPoints', 1000, 'NumStageOnePoints', 200);
[x2,F2,solutions] = run(gs,gsproblem);
I get nice estimates of x2.
So far so good.
Next I generate yhat (predicted y) based on the model estimate of x2 and the formula.
yhat = (1 + x2(1)*(x2(2)/0.0025662))*X(:,1) + (x2(2)/x2(3))*(sqrt(1-x2(1)^2))*X(:,2) + x2(4)*X(:,3) + x2(5)*X(:,4) + x2(6)*X(:,5) + x2(7)*X(:,6) + x2(8)*X(:,7) + x2(9)*X(:,8) + x2(10)*X(:,9) + x2(11)*X(:,10) + x2(12)*X(:,11) + x2(13)*X(:,12) + x2(14)*X(:,13) + x2(15)*X(:,14);
I then compute the residual vector as
resid = y - yhat;
My bootstrapping set-up is as follows:
I define the number of bootstrap samples:
nboot = 10;
Then create 10 column vectors of bootstrap indices:
[~, bootIndices] = bootstrp(nboot, [], resid);
I then define the 18026 by 10 matrix of 10 residual column vectors:
bootresid = resid(bootIndices);
I then define the 18026 by 10 matrix of 10 pseudo-y vectors based on the bootstrapped residual vectors.
yboot = (1 + x2(1)*(x2(2)/0.0025662))*X(:,1) + (x2(2)/x2(3))*(sqrt(1-x2(1)^2))*X(:,2) + x2(4)*X(:,3) + x2(5)*X(:,4) + x2(6)*X(:,5) + x2(7)*X(:,6) + x2(8)*X(:,7) + x2(9)*X(:,8) + x2(10)*X(:,9) + x2(11)*X(:,10) + x2(12)*X(:,11) + x2(13)*X(:,12) + x2(14)*X(:,13) + x2(15)*X(:,14) + bootresid;
I initialize my parameter matrix of 10 rows and 15 columns to be a matrix of zeros.
xb = zeros([10 15]);
I modify my original function to accommodate yboot (instead of y) and xb instead of x2
function Fboot = mybootfun(xb,yboot,X)
Fboot = sum((((yboot - (1 + xb(1)*(xb(2)/0.0025662))*X(:,1) - (xb(2)/xb(3))*(sqrt(1-xb(1)^2))*X(:,2) - xb(4)*X(:,3) - xb(5)*X(:,4) - xb(6)*X(:,5) - xb(7)*X(:,6) - xb(8)*X(:,7) - xb(9)*X(:,8) - xb(10)*X(:,9) - xb(11)*X(:,10) - xb(12)*X(:,11) - xb(13)*X(:,12) - xb(14)*X(:,13) - xb(15)*X(:,14)))).^2);
end
Then I run:
for i = 1:10
problem_b = createOptimProblem('fmincon','objective',@(xb) mybootfun(xb(i,:), yboot(:,i), X),'x0',x0,'lb',lb,'ub',ub);
[xb(i,:),Fboot,solution] = run(gs,problem_b);
end
Results are below. The first run works.
GlobalSearch stopped because it analyzed all the trial points.
All 15 local solver runs converged with a positive local solver exit flag.
After this, the error message that appears when the do loop increments to 2 is as follows.
Index in position 1 exceeds array bounds (must not exceed 1).
Error in @(xb)mybootfun(xb(i,:),yboot(:,i),X)
Error in fmincon (line 546) initVals.f = feval(funfcn{3},X,varargin{:});
Error in globalsearchnlp
Error in GlobalSearch/run (line 340) globalsearchnlp(FUN,X0,A,B,Aeq,Beq,LB,UB,NONLCON,options,localOptions);
Caused by: Failure in initial objective function evaluation. FMINCON cannot continue. Failure in initial call to fmincon with user-supplied problem structure.

Best Answer

Dear All,
There was a simple solution.
I replaced:
for i = 1:100
problem_b = createOptimProblem('fmincon','objective',@(xb) mybootfun(xb(i,:),yboot(:,i), X),'x0',x0,'lb',lb,'ub',ub);
[xb(i,:),Fboot,solution] = run(gs,problem_b);
end
with
for i = 1:100
problem_b = createOptimProblem('fmincon','objective',@(xb) mybootfun(xb,yboot(:,i), X),'x0',x0,'lb',lb,'ub',ub);
[xb(i,:),Fboot,solution] = run(gs,problem_b);
end
and that did the trick.
Thanks so much! Srini