Optimization of known function with respect to two unknown function arguments

Question

I have a data set, composed of points $(x_i, y_i)$ for $i=1,N$. I also have a known function $F$, which maps these points $x_i$ to $y_i$ as such $F(x_i, a(x_i),b(x_i)) = y_i$, where $a(x_i)$ and $b(x_i)$ are unknown functions of $x_i$. What I'd like to do is run an optimization (ideally using Python packages) in order to determine the functions $a$ and $b$ such that the error is minimized and also such that the functions $a$ and $b$ are smooth.

One, possibly naive way of doing this is to bin my data with respect to the $x_i$ variable, and then solve the optimization problem $F(x_i, a,b) = y_i$ for the single values $a,b$ for all points within the bin. Then, using this discretized solution for $a,b$, perform a numerical fit which results in a smooth curve.

However, I'm wondering if this is the best way to go about it. Is there a method that I could use which would determine the $a,b$ values simultaneously while trying to create a smooth curve from them? Or even if there were, would the effect be the same as if I perform the naive binned-optimization as I've described above?

Answer 1

Your problem is ill posed in the sense that if a solution exists, it is not unique. To see this, let ($a_i^*$, $b_i^*$), denote a solution of the equation $F(x_i, a_i,b_i)=y_i$, $i\in\{1,\ldots,N\}$ (or the minimizers of an appropriate cost function if equality cannot be achieved). Then, any pair of functions $\{a^*(x)$, $b^*(x)\}$ such that $a^*(x_i)=a_i^*$, $b^*(x_i)=b_i^*$, for all $i\in\{1,\ldots,N\}$, is a solution to your problem. Clearly, the set of all possible solutions is uncountable, even with the additional requirement that they are 'smooth'.

In order to avoid this ambiguity in the solution, one approach is to impose additional conditions. For example, you could restrict the solution functions to belong to a finite dimensional (and 'smooth') space (e.g., the space of polynomials of order $n$). In that case, the (difficult) functional optimization problem is converted to the (much simpler) problem of determination of a finite set of variables.

Answer 2

You could discretize the unknown functions $a$ and $b$, then formulate this as a least squares optimization problem with smoothing regularization. Specifically:

$$\min_{\mathbf{a}, \mathbf{b}} \frac{1}{2}\sum_{i=1}^N\left( f_i - F(x_i, \mathbf{m}_i^T \mathbf{a}, \mathbf{m}_i^T \mathbf{b}) \right)^2 + \frac{\alpha}{2} \mathbf{a}^T R \mathbf{a} + \frac{\alpha}{2} \mathbf{b}^T R \mathbf{b}$$

where

• $\mathbf{a}$ and $\mathbf{b}$ are coefficients of $a$ and $b$ with respect to some discrete set of basis functions $\{\phi_j\}_{j=1}^m$ (say, a piecewise linear basis of "hat" functions), $$a(x) = \sum_{j=1}^m \mathbf{a}_j \phi_j (x)$$
• $f_i$ is the $i$'th piece of data you are trying to fit,
• $\mathbf{m}_i^T := \begin{bmatrix}\phi(x_1) & \phi(x_2) & \dots & \phi(x_N)\end{bmatrix}$ is the "measure at $x_i$ vector" such that $$\mathbf{m}_i^T \mathbf{a} = a(x_i)$$
• R is a discrete form of the Laplacian with Neumann boundary conditions, plus some small multiple $\tau$ of the identity: $$R := \Delta_N + \tau I.$$

This could then be solved with the Gauss-Newton method, or L-BFGS, or something like that, (it is straightforward enough to figure out the derivative of the objective function by hand). In Python you could do the discretization with FEniCS, and the optimization with Scipy, though it will involve knowing a little about finite elements and a little about optimization.