3

I am trying to perform an empirical likelihood estimation in a regression setting using JuliaOpt (Convex or JuMP) and ran into difficulties using either API.

The problem:

Empirical likelihood for regression can be written as a maximization problem in the form:

$$max_{\theta, \pi_t} \sum\log(\pi_t)\\s.th\ \pi_t > 0, \sum \pi_t = 1, \sum \pi_t(y_i - h(x_t;\theta)) = 0$$

where $y$ is the depended variable, $x$ the explanatory variable, $h$ is in this case a linear function and $\theta$ is a vector of the regression parameters.

Here is my code:

function mele(  )
    x = rand(10,2)
    y = x*[2.1,12.7]

    π = Convex.Variable(10)
    b = Convex.Variable(2)
    p = maximize( sum(log(π)) )
    p.constraints += [ sum(π) == 1 ]
    p.constraints += [ π >= 0 ]
    p.constraints += [ sum( π'*(y - x*b) ) == 0  ]
    solve!(p, SCSSolver(verbose = 0))
    println( b.value )

    #returns
    > WARNING: Expression not DCP compliant. 
    > Trying to solve non-DCP compliant problems 
    > can lead to unexpected behavior.
    > WARNING: Problem not DCP compliant: constraint 3 is not DCP

    #If ceteris paribus (the \pi' is gone)
    p.constraints += [ sum( y - x*b ) == 0  ]

    #then julia returns
    > [6.48,7.57]

    #which is unfortunately wrong.
    #hence since I can not multiply the third constraint with
    #\pi, Convex fails yo solve this problem.

    #hence I tried JuMP using linear and non linear solvers:

    #m = Model( solver = NLoptSolver( algorithm = :LD_LBFGS ) )
    m = Model( solver = SCSSolver( verbose = 0 ) )
    @defVar( m, b[1:2] )
    @defVar( m, π[1:10] >= 0 )
    @setObjective( m, Max, sum{ log(π[i]),i=1:10 } )
    @addConstraint( m, sum{ π[i],i=1:10 } == 0 )
    @addConstraint( m, sum{π[i]*(y[i]- sum{x[i,j]*b[j], j=1:2}), i=1:10} == 0 )
    solve(m)
    println( getValue( b ) )
end
mele()


    #where the SCSSolver variant returns:
    > ERROR: LoadError: log is not defined for type Variable. 
    > Are you trying to build a nonlinear problem? 
    > Make sure you use @addNLConstraint/@setNLObjective.

    #hence:
    @setNLObjective( m, Max, sum{ log(π[i]),i=1:10 } )
    @addNLConstraint( m, sum{ π[i],i=1:10 } == 0 )
    @addNLConstraint( m, sum{π[i]*(y[i]- sum{x[i,j]*b[j], j=1:2}), i=1:10} == 0 )

    #and now:
    > ERROR: LoadError: MethodError: `loadnonlinearproblem!` has 
    > no method matching loadnonlinearproblem!

    #and hence the variant with the NLoptSolver:
    > ERROR: LoadError: ArgumentError: invalid NLopt arguments

From my understanding Convex has the the problem that $\pi*x*\theta$ has two unknown paramters in the multiplication which is not DCP compliant (correct me if I am wrong). JuMP on the other hand can not hande the log function in the problem statement.

Any ideas how to solve this?

kjetil b halvorsen
  • 63,378
  • 26
  • 142
  • 467
Vincent
  • 253
  • 1
  • 2
  • 7
  • 1
    Maybe you should post over at https://groups.google.com/forum/#!forum/julia-opt . – Mark L. Stone Oct 17 '15 at 18:45
  • $\pi x \theta$ is not only not DCP compliant, it's not (jointly) convex, so you're not going to be able to use a convex solver for this formulation. – Danica Oct 17 '15 at 19:43

0 Answers0