How to avoid creating many binary switching variables in GEKKO

I am solving for 14 variables by minimizing on the order of thousands of equations with IMODE = 3 in GEKKO.

Each equation is the squared error between the true response and the prediction of a P-spline model (i.e., penalized B-spline):

eq[i] = m.Minimize((y_true[i] - spline(coeffs, knots, vars)[i]) ** 2).

The spline models are constituted of their coefficients and knots (which are previously calculated) along with the 14 variables to optimize.

When building the P-spline models for GEKKO, I need to check in-between which knots the value of a variable lies. I tried using both m.if2 and m.if3 to achieve this; however, both of these logical functions create many additional binary switching variables, especially for splines with many pieces. In the end, I end up with tens of thousands of binary switching variables. These outnumber the equations, resulting in the number of degrees of freedom being less than 0.

Question: How can I avoid using m.if2 or m.if3 to build my splines?

Note: I am aware that GEKKO has the prebuilt object m.bspline; however, it appears that it can only do 2D B-splines with two independent variables, while my splines can have over ten independent variables.

The m.if2() function uses MPCCs instead of binary variables. Additional slack variables and equations are created when the m.if2() function is used. The m.if2() function can be numerically unreliable. The m.if3() function is typically more reliable, but it requires the additional binary variables and an MINLP solver (APOPT) to find the solution. Negative degrees of freedom is a warning and may still solve because inequality constraints may not be active at the solution.

One additional thing to consider is the spline(coeffs, knots, vars)[i] is not re-evaluated during the solution. The Gekko model is compiled into byte-code and there are no call-backs into the original functions. Because the spline is evaluated once during the model initialization, could the position be determined either before or after the solve?

A final suggestion is to use a machine learned model such as a Neural Network, Gaussian Processes, Support Vector Regressor, or similar form that can be imported into Gekko with ML functions after it has been trained with sklearn or tensorflow / keras.