fminsearch error while optimizing over one value in neoclassical growth model

I am trying to solve the stochatic neoclassical growth model using fminsearch. My main function has the following code:

clear all;
close all;
tic
global v0 beta delta alpha kmat k0 s prob a0 s j

%set up the fised paramters
alpha = 0.33; %capital share
beta = 0.95;%discount rate
delta = 0.1; %rate of depreciation of capital
s = 2;


tol = 0.01;%tolerance
maxits = 1000;%maximum iteration
dif = tol + 100;
its = 0;

kstar = (alpha/(1/beta-(1-delta)))^(1/(1-alpha));
cstar = kstar^(alpha)-delta*kstar;
istar = delta * kstar;
ystar = kstar ^ alpha;

kmin = 0.25*kstar;
kmax = 1.75*kstar;
kgrid = 99;
grid = (kmax - kmin)/kgrid;
kmat = kmin:grid:kmax;
kmat = kmat';
[N,n] = size(kmat);

polfun = zeros(kgrid+1,3);
v0 = zeros(N,1);

amat = [0.9 1 1.1]';
prob = (1/3) * ones (3,3);


while dif>tol && its < maxits
    
    
    for j = 1:3
        for i = 1:N
            k0 = kmat(i,1);
            a0 = amat(j,1);
            %k1 = fminsearch(@valfun_stoch,[2.3706],[kmin],[kmax]);
            %k1 = fminsearchbnd(@(x)SolveModelParallel(x),xval,[],[],opts);
            %fun = @(k) valfun_stoch;
            kset=(kmax+kmin)/2;
            k1 = fminsearch(@valfun_stoch,kset);
            v1(i,j) = -valfun_stoch(k1);
            k11(i,j) = k1;
        end
    end
    %g = abs(v1−v0);
    dif = norm(v1-v0);
    v0 = v1;
    its = its+1;
end

for i = 1:N
    con(i,1) = kmat(i,1)^alpha - k11(i,1) + (1-delta)*kmat(i,1);
    polfun(i,1) = kmat(i,1)^(alpha) - k11(i,1) + (1-delta)*kmat(i,1);
end

The valfun_stoch is as follows:

function val=valfun_stoch(k)

% This program gets the value function for a stochastic growth model with
% CRRA utility

global v0 beta delta alpha kmat k0 prob a0 s j

g = interp1(kmat,v0,k,'linear'); % smooths out previous value function

c = a0*k0^alpha - k + (1-delta)*k0; % consumption
if c <= 0
    val = -8888888888888888 - 800*abs(c); % keeps it from going negative
else
    val = (1/(1 - s))*(c^(1 - s) - 1) + beta*(g*prob(j,:)');
end
val = -val; %

I keep getting the following error:

Unable to perform assignment because the size of the left side is 1-by-1 and the size of the right side is 3-by-1.

Error in fminsearch (line 335)
fv(:,end) = fxe;

Error in main2 (line 63)
k1 = fminsearch(@valfun_stoch,kset);

I am not sure where I am going wrong. I am only trying to optimize over capital k, so I am not sure why I should get 3-1 output.

Any help would be greatly appreciated.

[Torsten]

You overwrite v0 here:

v0 = v1;

So you have the same problem as in your previous question.

Furthermore, the line

val = (1/(1 - s))*(c^(1 - s) - 1) + beta*(g*prob(j,:)');

assigns a (3x1) vector to "val", but "val" must be a scalar (1x1).

[Prerna Mishra]

I expanded the problem to a 3 dimensional one as follows:

tic

global disc alpha beta R f_t M sigma W_t theta chi c kappa avec kvec dvec rho_a sigma_a 
global i j l blog_mean blog_var k0 d0 v0 a0 b_sorted b d pi v1 val kd1
 
na = 7; %number of points on shock grid
nk = 11; %number points on capital grid
nd = 3; %number of points on debt grid

%fixed parameters%
disc = 0.96; %discount rate
alpha = 0.364; %elasticity of capital
beta = 0.115; %elasticity of labor
R = 3.131696599;%market size 
f_t = -0.500185113; %fixed cost of production
M = log(1056);%Mass of firms
sigma = -1.336449027; %degree of product substitutability
W_t = log(29448.08908); %Wage data from the WageDataCompiles

%parameters we would calibrate using SMM
theta = 0.3;%cost of bankruptcy
chi = 0.5;%level of debt at which bankruptcy happens
c = 0.5; %parameters
kappa = 0.5; %managerial rent

%The next portion creates grids points for the TFP process, the capital
%grid and the debt grid

%generating the shocks and the transition probability matrix
multiple = 2.575; %how many std dev's should the prod grid span?
emean = 0; 
rho_a = 0.5; %to be estimated using SMM, persistence of the TFP process
sigma_a = 0.5;  %to be estimated using SMM, variace of the TFP process
[avec, pi] = tauchen(0.0, sigma_a, rho_a, multiple, na); %na = 11
avec = exp(avec)';

%creating grid for capital and debt
kmin = log(.5);
kmax = log(5000000);
kvec = linspace(kmin,kmax,nk)';
kvec = exp(kvec);
dmin = log(.25);
dmax = log(25000);
dvec = linspace(dmin,dmax,nd)';
dvec = exp(dvec);

%Parameters for the firm level shock, we will calibrate this using
%structural estimation
blog_mean = 0.01;
blog_var = 0.3;
blog_sigma = sqrt(blog_var);
[b_sorted,b] = RandomLogNormal(blog_mean,blog_sigma);

%Setting the parameters for iterations
tol = 0.1;%tolerance
maxits = 5000;%maximum iteration
dif = tol + 100;
its = 0;

%polfun = zeros(kgrid+1,3);
v0 = zeros(nd,nk,na);
v1 = zeros(nd,nk,na);

while dif>tol && its<maxits
    for i = 1:na
      for j = 1:nk
         for l = 1:nd
            k0 = kvec(j,1);
            %d0 = dvec(l,1);
            a0 = avec(i,1);
            kset = (kmin+kmax)/2;
            dset = (dmin+dmax)/2;
            %f0 = [kset,dset];
            %k1,d1 = fminbnd(@RmatVal,kmin,kmax);
            kd = [kset,dset];
            kd = fminsearch(@RmatVal,kd);
            kopt = kd1(1);
            dopt = kd1(2);
            v1(l,j,i) = -RmatVal(kd);
            k11(l,j,i) = kopt;
            d11(l,j,i) = dopt;
    
         end
      end
    end
dif = norm(v1-v0);
v0 = v1;
its = its + 1;
end

toc

In each iteration of i, I will get, 7*11 values of k and 7*11 values of d as output of fminsearch. I am struggling as to how to appropriately store it in the k11 and d11 matrices. I am not sure how to exactly extract values from the output of fminsearch.

For reference, a small most of my RmatVal looks as follows:

function val = RmatVal(kd)

global disc alpha beta R f_t M sigma W_t theta chi c kappa avec kvec dvec rho_a sigma_a 
global i j l blog_mean blog_var k0 d0 v0 a0 b_sorted b pi v1 na nk nd
%g = interp1(kvec,v0,k,'linear');
%g = interp2(kvec,dvec,v0,k,d,'linear');
val = First_term+second_term+third_term; %which have k and d in them

val = val + disc * (val*pi(i,:)'); %pi is an na*na matrix of transition probability.
val = -val; 

[Torsten]

Maybe you could state the optimization problem in a mathematical way, not with MATLAB code. I think advice in your case can only be given if the underlying problem is understood.

Answer 1

[Gokul Nath S J]

Dear Prerna Mishra,

I have tried running your code on my machine. However, I faced the same issue. There are a few possible workarounds that I can suggest.

1. The error is triggered in line 40 while calling “fminsearch”. While looking at “valfun_stoch”, it is evident that you have two conditions for the final variable 'val'. The conditions are ‘c’ would be either less than or equal to zeros and if it’s greater than zero. The error is triggered once the variable ‘c’ is greater than zero.
2. Line number 10 in function file “valfun_stoch” might be one of the reasons why the code is throwing out error. In some scenariosthe value returned by the function is simply 1 x 1 array while in some scenarios it is a 3 x 1 array. Hence this mismatch is reflected in the main function. Note that the ‘prob’ variable in "valfun_stoch" function (line 10) is of dimension 3 x 3.
3. Also, confirm again if your global declaration is correct. In certain cases, the value of the arrays is changing (variable ‘g’ & ‘v0’ if 'val' is forced to 1 x 1 by changing prob(j,:) to a sample value of prob(2, 2)) since some of the variables determining its size is globally declared.