Question about optimize() function in Mata

[Nhan Nguyen]

Dear Stata users,

I am playing around with data and I just wonder why I get inconsistent results for one optimization problem when I try different approaches and how to use the KKT condition in optimize(). Let say I have the optimization problem as follow. I have the Var-Covariance matrix as

m Var_Covar = (2,1.3,0.35,0\1.3,1.5,0.8,0\0.35,0.8,0.6,0\0,0,0,1.2)

and I would want to find the combination of stocks (by adjusting how many of them that I should hold in a portfolio) that minimizes the global variance. Therefore, the optimization problem can be expressed as min{wVw'} Subject to rowsum(w) = 1 and w_i >= 0 where w is the (1x4) vector of weights. This is the convex optimization problem because the objective function is convex and continuous. Here is how I use optimize() to tackle the problem.

1. If I relax the inequality constraints and proceed with optimize() then I have two ways of coding which should return the same (or at least similar to some extent) results.
- The first way is using the optimize_init_constraints() to specify the equality constraint (the sum of the proportions of each stock in the portfolio is one).

mata
w1 = (0.2,0.2,0.2,0.2)'
void globalminvar2(minvar, params, sig, y, g, H)
    { 
        y = -params*sig*params'
        if(minvar==1){
            g = (-2:*params*sig)
        }
    }
T = optimize_init()
optimize_init_evaluator(T,&globalminvar2())
optimize_init_evaluatortype(T,"d1")
optimize_init_params(T,(w1'))
optimize_init_argument(T,1,Var_Covar)
optimize_init_constraints(T,(1,1,1,1,1))
p = optimize(T)
p
end

- The second way is adding the Lagrange multiplier to the objective function which means I no longer have to declare my constraint optimize_init_constraints(T,(1,1,1,1,1)). Below is my coding for the problem.

mata
void globalminvar3(minvar, params, sig, y, g, H)
{
    w = params[1..rows(sig)]
    lambda = params[cols(params)]
    y = -w*sig*w'-(rowsum(w)-1)*lambda
    if(minvar==1){
         g[1] = (-2*w*sig[.,1]:-lambda)
         g[2] = (-2*w*sig[.,2]:-lambda)
         g[3] = (-2*w*sig[.,3]:-lambda)
         g[4] = (-2*w*sig[.,4]:-lambda)
         g[5] = -rowsum(w)+1
    }
}
U = optimize_init()
optimize_init_evaluator(U,&globalminvar3())
optimize_init_evaluatortype(U,"d1")
optimize_init_params(U,(w1',0.2))
optimize_init_argument(U,1,Var_Covar)
p = optimize(U)
p
end

The result from using optimize_init_constraints(T,(1,1,1,1,1))is correct because rowsum(p) is equal to 1. However, the result I got from the second way is not binding. You can see by checking rowsum(p[1..4]). At this point, I would want to know which part that I understand wrong. My conjecture is that the problem that results in inconsistent outcome is the initial guess on the parameters because every time I feed in a different parameter for lambda, the result changes. This optimization problem is quite basic and I can do it in other statistical softwares but I mostly use Stata so I truly want to know how to do it.

2. Now I want the weight of each stock to be non-negative (no-short selling), which means now I have to use the KKT conditions. Let theta_i be the Lagrange multiplier for each of the inequality condition, we have the complementary slackness conditions as theta_i > 0, theta_i = 0 if w_i > 0. Below is the implementation of the above idea.

mata
Return = (2.5,2,1,1.5)'
Var_Covar = (2,1.3,0.35,0\1.3,1.5,0.8,0\0.35,0.8,0.6,0\0,0,0,1.2)
Corr_mat = Var_Covar:/(sqrt(diagonal(Var_Covar))#sqrt(diagonal(Var_Covar)'))
w1 = (0.2,0.2,0.2,0.2)'
void globalminvar1(minvar, params, sig, y, g, H)
{
    w = params[1..rows(sig)]
    lambda = params[rows(sig)+1]
    theta = params[rows(sig)+2..cols(params)]
    y = -w*sig*w'+(rowsum(w)-1)*lambda+w*theta'
    if(minvar==1){
        theta[selectindex(w:>0)]=J(rows(selectindex(w:>0)),cols(selectindex(w:>0)),0)
        g = (-2:*w*sig:+lambda:+theta,rowsum(w)-1,theta)
    }
}
S = optimize_init()
optimize_init_evaluator(S,&globalminvar1())
optimize_init_evaluatortype(S,"d1")
optimize_init_params(S,(w1',1,0,0,0,0))
optimize_init_argument(S,1,Var_Covar)
p = optimize(S)
p
end

However, mata could not find the optimum value of this solution. Can anyone help me with this. Thank you in advance.