The function value and the fval don't match after fminsearch

Hi everyone,

I seem to be having an issue with fminsearch. I'm minimizing over a function and it seems to be working, but the value of the function (with the optimal values plugged in) and the number it gives me for the fval don't match. I can't seem to find a problem, but maybe someone out here can take a look. I really appreciate it. I've pasted the two sets of code below. Also notice that the value of the minimized function is one of the variables I solve for in the second file.

This first one sets up the function

%% Solving for the Steady State Values for WP3 % We solve a system of 12 equations using 2 additional to % to test the residituals. % We choose two parameters to start and use this function to % solve for the resulting residuals.

    function[resid] = Base_Steady(choice, phi_b, alpha, phi_k, tau,...
                                          pi_ss, z_ss, x_a, theta_m,...
                                          delta, eta, beta, psi,...
                                          theta_y, phi_m1, phi_m2, g_ss)
%%Solving 7 variables by hand

    r = pi_ss*z_ss/beta;
    w = (theta_y-1)/theta_y*z_ss;
    lam3 = eta/w;
    lam1 = r*lam3;
    r_k = phi_k+1;
    r_u = phi_m1*r^phi_m2;
    r_l = (1-psi)*r_k+psi*r_u;
%%Bounding the Choice Variables
    a_s = 1/(1+exp(choice(1)));     % a_s in (0,1)
    l = 2/(1+exp(choice(2)));       % l in (0,2)
%%Base Steady State Equations
    l_s = l/a_s;                                                %((A.2.2 ))
    s_b = phi_b*a_s^alpha;                                      %((A.2.1 ))
    r_d = a_s*(1-tau)*r_l+(1-a_s*(1-tau))*r_u-(r_u-1)*s_b*(1-tau);%(A.2.5))
    d = l_s/(1-tau);                                            %((A.2.3 ))
    l_d = -l_s*log(1-(l/l_s));                                  %((A.2.23))
    a_h = l/l_d;                                                %((A.2.7 ))
    k = (1-a_h)*l_d;                                            %((A.2.8 ))
    lam2 = (lam3*(a_h*r_l+(1-a_h)*r_k)-lam1*a_h)/(1-a_h);       %((A.2.17))
    c = 1/lam2;                                                 %((A.2.19))
    n = c-k;                                                    %((A.2.10))
    m_a = (x_a^(1/theta_m)*n^((theta_m-1)/theta_m)...
           +(1-x_a)^(1/theta_m)*d^((theta_m-1)/theta_m))...
           ^(theta_m/(theta_m-1));                              %((A.2.12))
    lam4 = delta/m_a;                                           %((A.2.16))
    resid = abs(d-(1-x_a)*m_a*(lam4/(lam1-r_d*lam3))^theta_m)...%((A.2.14))
            +abs(n-x_a*m_a*(lam4/(lam1-lam2))^theta_m);         %((A.2.15))

The second set here runs the fminsearch and solves the system at the optimal values.

%% Solving for the Steady State Values of WP3 % Here, we call on Base_Steady.m in order to minimize the residual. % Then we relist all the steady state equations to solve for the % steady state values of the model.

%% Setting the Parameter Values % Only those needed to solve the steady state equations are listed.

phi_b = 0.01; alpha = 1.00; tau = 0.002; phi_k = 0.0333; pi_ss = 1.005; z_ss = 1.005; x_a = 0.70; theta_m = 0.50; delta = 2.00; eta = 2.50; beta = 0.995; psi = 0.75; theta_y = 6.00; phi_m1 = 0.999375; phi_m2 = 1.0; g_ss = 1.005;

%% Using fminsearch to solve Base_Steady.m

options = optimset('Display','iter','TolX',1e-12,'TolFun',1e-12,'MaxFunEvals',1e10,'MaxIter',1e10); [choice_star, fval, exitflag] = fminsearch(@(choice)Base_Steady(choice,... phi_b, alpha, tau, phi_k, pi_ss,... z_ss, x_a, theta_m, delta, eta, beta,... psi, theta_y, phi_m1, phi_m2, g_ss),... [0,0], options);

%% Retransformation of the SOlution Variables

a_s = 1/(1+exp(choice_star(1))) l = 2/(1+exp(choice_star(2)))

%% Solving the Rest of the Model

    r = pi_ss*z_ss/beta
    w = (theta_y-1)/theta_y*z_ss
    lam3 = eta/w
    lam1 = r*lam3
    r_k = phi_k+1
    r_u = phi_m1*r^phi_m2
    r_l = (1-psi)*r_k+psi*r_u
    l_s = l/a_s                                                
    s_b = phi_b*a_s^alpha                                      
    r_d = a_s*(1-tau)*r_l+(1-a_s*(1-tau))*r_u-(r_u-1)*s_b*(1-tau)%(A.2.5))
    d = l_s/(1-tau)                                            %((A.2.3 ))
    l_d = -l_s*log(1-(l/l_s))                                  %((A.2.23))
    a_h = l/l_d                                                %((A.2.7 ))
    k = (1-a_h)*l_d                                            %((A.2.8 ))
    lam2 = (lam3*(a_h*r_l+(1-a_h)*r_k)-lam1*a_h)/(1-a_h)       %((A.2.17))
    c = 1/lam2                                                 %((A.2.19))
    n = c-k                                                    %((A.2.10))
    m_a = (x_a^(1/theta_m)*n^((theta_m-1)/theta_m)...
           +(1-x_a)^(1/theta_m)*d^((theta_m-1)/theta_m))...
           ^(theta_m/(theta_m-1))                              %((A.2.12))
    lam4 = delta/m_a                                           %((A.2.16))
    resid = abs(d-(1-x_a)*m_a*(lam4/(lam1-r_d*lam3))^theta_m)...%((A.2.14))
            +abs(n-x_a*m_a*(lam4/(lam1-lam2))^theta_m)         %((A.2.15))

    n_u = (1-a_s)*l_s+tau*d
    m = pi_ss*z_ss*(d+n-l)
    h = (m-r_d*d+r_l*l+r_k*k)/z_ss
    y = z_ss*h
    f = (1-w/z_ss)*y
    q = f/(r-1)

[Matt J]

I can't reproduce the problem. On my machine, I get the same values precisely:

    K>> func=@(choice)Base_Steady(choice,...
                 phi_b, alpha, tau, phi_k, pi_ss,... 
                 z_ss, x_a, theta_m, delta, eta, beta,...
                 psi, theta_y, phi_m1, phi_m2, g_ss);

    K>> isequal(func(choice_star),fval)

    ans =

         1