# Non-linear optimization programming, with step function in constraint

I want to optimize a non-linear function $$f(x)$$, $$f: \mathcal{R}^{n} \to \mathcal{R}$$ (being a log-likelihood over $$m$$ observations, i.e. $$i$$ being the observation index) under constraints numerically, where one of the constraints contains a non-linear function $$g(x)$$ that is an expectation over a step function. The optimization problem can be stated as

\begin{align} &\text{max}_{x \in \mathcal{R}^{n}} \; f(x) \\ &\text{s.t.} \; g(x) = 0 \\ & \; x \in [a,b], \; a,b \in \mathcal{R}^{n}\\ & \text{where} \; g(x) = \frac{1}{m} \sum_{i=1}^{m} s(x)_{i} - c, \; c \in [0, 1], \; g: \mathcal{R}^{n} \to [-1,1] \\ & \text{with} \; s(x)_{i} = \begin{cases} 0.1 & 0 \leq j(x)_{i} < 0.2 \\ \vdots & \vdots \\ 0.9 & 0.8 \leq j(x)_{i} \leq 1. \end{cases} \end{align}

$$j(x)_{i}$$ is the $$i-th$$ component of the vector-valued, monotonic function $$j: \mathcal{R}^{n} \to [0,1]^{m}$$. While $$s(x)_{i}$$ is the $$i-th$$ component of the vector-valued function $$s: \mathcal{R}^{n} \to [0.1, \dots, 0.9]^{m}$$. I'm not familiar enough with optimization to know how to proceed i.e. using the appropriate optimization technique, and if it is solveable in the first place. I tried to solve it with SLSQP but got no where, after 5 iterations in scipy I got the error: Singular matrix C in LSQ subproblem. I found that as soon as indicator variables are involved in non-linear programming (optimization) one needs Mixed-Integer Nonlinear Programs which are quite computationally intensive. Some information (reference, links) that help with tackling this problem would be really cool.

• Do you need a provably optimal solution or would you settle for a "good" solution? Commented Feb 4 at 16:08
• I would settle for any solution, so a "good" solution is definitely enough :D. Maybe also some background, so essentially i want to train parameters of a binary classifier however to the constraint that the class probabilites are mapped to expert based probabilities (the step function) and the average expert based probabilites should match with some long run probabilitiy of the observed data (the constraint that $g$ is zero). Commented Feb 6 at 9:09

One possibility would be to try a penalty approach. For a suitably large value of $$\lambda > 0$$ you can try maximizing $$h(x) = f(x) -\lambda g(x)^2$$ over the hyperrectangle $$[a, b].$$ If the maximizer has a nonzero value for $$g(x),$$ ratchet up $$\lambda$$ and try again. Since $$h$$ is nondifferentiable, I would be inclined to try a method like Nelder-Mead.