Generation of point sets by convex optimization for interpolation in reproducing kernel Hilbert spaces

Abstract

We propose algorithms to take point sets for kernel-based interpolation of functions in reproducing kernel Hilbert spaces (RKHSs) by convex optimization. We consider the case of kernels with the Mercer expansion and propose an algorithm by deriving a second-order cone programming (SOCP) problem that yields n points at one sitting for a given integer n. In addition, by modifying the SOCP problem slightly, we propose another sequential algorithm that adds an arbitrary number of new points in each step. Numerical experiments show that in several cases the proposed algorithms compete with the P-greedy algorithm, which is known to provide nearly optimal points.

Appendix. Remarks on implementation of SOCP (4.6) by MOSEK

1.1 MOSEK

As a software to solve SOCP problems, we use the MOSEK Optimization Toolbox for MATLAB 8.1.0.56, which is provided by MOSEK ApS in Denmark (https://www.mosek.com/, last accessed on 19 October 2018).

1.1.1 Expression of rotated second-order cones

Recall that $p = \lceil \log _{2} \ell \rceil$. In SOCP (4.6) in which ℓ is even, we need to consider the rotated second-order cone constraints

$$ \begin{array}{@{}rcl@{}} && z_{ij}^{2} \leq t_{ij} w_{i} \qquad (i=1,\ldots,m, j=1,\ldots, \ell), \end{array} $$

(5.6)

$$ \begin{array}{@{}rcl@{}} && {u_{i}^{2}} \leq u_{2i} u_{2i+1} \qquad (i=1,\ldots, 2^{p-1}-1), \end{array} $$

(5.7)

$$ \begin{array}{@{}rcl@{}} && {u_{i}^{2}} \leq \tilde{g}_{2i-2^{p}+1} \tilde{g}_{2i-2^{p}+2} \qquad (i = 2^{p-1}, \ldots, 2^{p-1} + \ell/2 - 1). \end{array} $$

(5.8)

Furthermore, in the counterpart in which ℓ is odd, we need to consider the constraints

$$ \begin{array}{@{}rcl@{}} && z_{ij}^{2} \leq t_{ij} w_{i} \qquad (i=1,\ldots,m, j=1,\ldots, \ell), \end{array} $$

(5.9)

$$ \begin{array}{@{}rcl@{}} && {u_{i}^{2}} \leq u_{2i} u_{2i+1} \qquad (i=1,\ldots, 2^{p-1}-1), \end{array} $$

(5.10)

$$ \begin{array}{@{}rcl@{}} && {u_{i}^{2}} \leq \tilde{g}_{2i-2^{p}+1} \tilde{g}_{2i-2^{p}+2} \qquad (i = 2^{p-1}, \ldots, 2^{p-1} + (\ell-3)/2), \end{array} $$

(5.11)

$$ \begin{array}{@{}rcl@{}} && {u_{i}^{2}} \leq \tilde{g}_{2i-2^{p}+1} u_{1} \qquad (i = 2^{p-1} + (\ell-1)/2). \end{array} $$

(5.12)

However, MOSEK permits only the expression ${\eta _{1}^{2}} + {\cdots } + {\eta _{N}^{2}} \leq 2 \xi \zeta$ for a rotated second-order cone, where we just need the case that N = 1. Therefore, we employ the following variable transformations:

$$ \begin{array}{@{}rcl@{}} && w_{i} = 2 \hat{w}_{i} \quad (i = 1,\ldots, m), \end{array} $$

(5.13)

$$ \begin{array}{@{}rcl@{}} && \tilde{g}_{i} = \sqrt{2} \hat{g}_{i} \quad (i = 1, \ldots, \ell-1), \end{array} $$

(5.14)

$$ \begin{array}{@{}rcl@{}} && \tilde{g}_{\ell} = \left\{\begin{array}{ll} \sqrt{2} \hat{g}_{\ell} & (\ell \text{ is even}), \\ 2 \hat{g}_{\ell} & (\ell \text{ is odd}). \end{array}\right. \end{array} $$

(5.15)

We do not change u_i because of the following reason. If we change the constraints ${u_{i}^{2}} \leq u_{2i} u_{2i+1}$ to ${u_{i}^{2}} \leq 2 u_{2i} u_{2i+1}$, they just change the optimal value by $2^{(p-1)2^{p-1}}$ times.⁵ Consequently, we deal with the SOCP problem

(5.16)

in the case that ℓ is even, and

(5.17)

in the case that ℓ is odd.

1.1.2 Restriction about cone constraints

MOSEK does not permit that a variable belongs to plural cones. However, we need to consider the cone constraints in (5.16) and (5.17) in which each u_i appears 2 times and each $\hat {w}_{i}$ appears ℓ times. Therefore, we prepare the copies of these variables:

$$ \begin{array}{@{}rcl@{}} u_{i} \ \leftrightarrow \ u_{i1}, u_{i2}, \qquad \hat{w}_{i} \ \leftrightarrow \ \hat{w}_{i1}, \hat{w}_{i2}, \ldots, \hat{w}_{i \ell}, \end{array} $$

and set the linear constraints between these variables:

$$ \begin{array}{@{}rcl@{}} && u_{i1} = u_{i2}, \end{array} $$

(5.18)

$$ \begin{array}{@{}rcl@{}} && \hat{w}_{i1} = \hat{w}_{i2} = {\cdots} = \hat{w}_{i \ell}. \end{array} $$

(5.19)

Then, we introduce a vector

$$ \begin{array}{@{}rcl@{}} v = (&& u_{11}, u_{12}, \ u_{21}, u_{22}, \ \ldots, \ u_{(2^{p}-1), 1}, u_{(2^{p}-1), 2}; \\ && \hat{g}_{1}, \hat{g}_{2}, \ldots, \hat{g}_{\ell}; \\ && z_{11}, \ldots, z_{m1}, \ z_{12}, \ldots, z_{m2}, \ \ldots, \ z_{1\ell}, \ldots, z_{m\ell}; \\ && t_{11}, \ldots, t_{m1}, \ t_{12}, \ldots, t_{m2}, \ \ldots, \ t_{1\ell}, \ldots, t_{m\ell}; \\ && \hat{w}_{11}, \ldots, \hat{w}_{m1}, \ \hat{w}_{12}, \ldots, \hat{w}_{m2}, \ \ldots, \ \hat{w}_{1\ell}, \ldots, \hat{w}_{m\ell})^{T}. \end{array} $$

(5.20)

1.2 Linear constraints

To describe the linear constraints in (5.16) and (5.17) in terms of the vector v, we provide a matrix A and vectors b_l and b_u that express the constraints as b_l ≤ Av ≤ b_u.

We begin with the case that ℓ is even. First, we express the constraints in (5.18) and (5.19) by

$$ \begin{array}{@{}rcl@{}} \underbrace{ \left[\begin{array}{llllllll} 1 & -1 & 0 & 0 & {\cdots} & 0 & 0 & 0_{\ell + 3m\ell}^{T} \\ 0 & 0 & 1 & -1 & {\cdots} & 0 & 0 & 0_{\ell + 3m\ell}^{T} \\ {\vdots} & & & & {\ddots} & & {\vdots} & {\vdots} \\ 0 & 0 & 0 & 0 & {\cdots} & 1 & -1 & 0_{\ell + 3m\ell}^{T} \end{array}\right]}_{A_{1}} v = 0_{2^{p} - 1}, \end{array} $$

(5.21)

and

$$ \begin{array}{@{}rcl@{}} \underbrace{ \left[\begin{array}{lllllll} O_{m, 2(2^{p}-1) + \ell + 2m\ell} & I_{m} & -I_{m} & O_{m,m} & {\cdots} & O_{m,m} & \\ O_{m, 2(2^{p}-1) + \ell + 2m\ell} & O_{m,m} & I_{m} & -I_{m} & & O_{m,m} & \\ {\vdots} & & & {\ddots} & & & \\ O_{m, 2(2^{p}-1) + \ell + 2m\ell} & O_{m,m} & {\cdots} & O_{m,m} & I_{m} & -I_{m} & \end{array}\right]}_{A_{2}} v = 0_{m (\ell-1)}, \end{array} $$

(5.22)

respectively. Next, we consider the constraint

$$ \begin{array}{@{}rcl@{}} (a_{1}, \ldots, a_{m}) Z = \left[\begin{array}{llll} \sqrt{2} \hat{g}_{1} & 0 & {\cdots} & 0 \\ \ast & \sqrt{2} \hat{g}_{2} & {\ddots} & {\vdots} \\ \ast & \ast & {\ddots} & 0 \\ \ast & \ast & \ast & \sqrt{2} \hat{g}_{\ell} \end{array}\right]. \end{array} $$

By using v and the expression

$$ \begin{array}{@{}rcl@{}} (a_{1}, \ldots, a_{m}) = \left[\begin{array}{lll} \boldsymbol{\alpha}_{1} & {\cdots} & \boldsymbol{\alpha}_{\ell} \end{array}\right]^{T} \qquad (\boldsymbol{\alpha}_{j} \in \mathbf{R}^{m}), \end{array} $$

we express the linear constraint as follows:

(5.23)

Next, we express the constraint $\hat {w}_{1\ell } + {\cdots } + \hat {w}_{m\ell } = n/2$ by

$$ \begin{array}{@{}rcl@{}} \underbrace{ \left[\begin{array}{ll} 0_{2(2^{p}-1) + \ell + 3m(\ell - 1)}^{T} & {1_{m}^{T}} \end{array}\right]}_{A_{4}} v = \frac{n}{2}. \end{array} $$

(5.24)

Next, we express the inequality constraints

$$ \begin{array}{@{}rcl@{}} \sum\limits_{i=1}^{m} t_{ij} \leq \sqrt{2} \hat{g}_{j} \quad (j = 1,{\ldots} , \ell) \end{array} $$

by

(5.25)

Finally, we consider the constraints

$$ \begin{array}{@{}rcl@{}} {u_{i}^{2}} \leq {u_{1}^{2}} \qquad (i = 2^{p-1} + \ell/2, \ldots, 2^{p}-1), \end{array} $$

(5.26)

which are equivalent to the linear constraints u_i ≤ u₁ owing to the non-negative constraints u_i ≥ 0. These constraints are expressed by

(5.27)

Consequently, by setting

$$ \begin{array}{@{}rcl@{}} A & =& \left[\begin{array}{llllll} {A_{1}^{T}} & {A_{2}^{T}} & {A_{3}^{T}} & {A_{4}^{T}} & {A_{5}^{T}} & {A_{6}^{T}} \end{array}\right]^{T}, \end{array} $$

(5.28)

$$ \begin{array}{@{}rcl@{}} b_{\mathrm{l}} & =& \left( 0_{2^{p} - 1}^{T}, \ 0_{m (\ell-1)}^{T}, \ 0_{\ell (\ell+1)/2}^{T}, \ n/2, \ -\infty \cdot 1_{\ell}^{T} , \ -\infty \cdot 1_{2^{p-1}-\ell/2} \right)^{T}, \end{array} $$

(5.29)

$$ \begin{array}{@{}rcl@{}} b_{\mathrm{u}} & =& \left( 0_{2^{p} - 1}^{T}, \ 0_{m (\ell-1)}^{T}, \ 0_{\ell (\ell+1)/2}^{T}, \ n/2, \ 0_{\ell}^{T} , \ 0_{2^{p-1}-\ell/2}^{T} \right)^{T}, \end{array} $$

(5.30)

we can express the linear constraints as b_l ≤ Av ≤ b_u.

In the case that ℓ is odd, we need to modify A₃, A₅, A₆, b_l, and b_u. We replace the component “$-~\sqrt {2}$” in the last rows of A₃ and A₅ by “− 2”. Furthermore, we modify A₆, b_l, and b_u as

respectively.

1.3 Constraints for the ranges of the variables

Besides the constraints u_i ≥ 0 and $0 \leq \hat {w}_{j} \leq 1/2$ in (5.16) and (5.17), we can add $\hat {g}_{i} \geq 0$ and t_ij ≥ 0 because any solution that does not satisfy these constraints is useless. Therefore, by using the vectors given by

$$ \begin{array}{@{}rcl@{}} \tilde{b}_{\mathrm{l}} & =& \left( 0_{2(2^{p}-1) + \ell}^{T}, \ -\infty \cdot 1_{m\ell}^{T}, \ 0_{2m\ell}^{T} \right)^{T}, \end{array} $$

(3.4)

$$ \begin{array}{@{}rcl@{}} \tilde{b}_{\mathrm{u}} & =& \left( \infty \cdot 1_{2(2^{p}-1) + \ell + 2m\ell}^{T}, \ (1/2) \cdot 1_{m\ell}^{T} \right)^{T}, \end{array} $$

(3.5)

we can express the constraints for the ranges of the variables as $\tilde {b}_{\mathrm {l}} \leq v \leq \tilde {b}_{\mathrm {u}}$.

1.4 Cone constraints

To express the cone constraints in (5.16) and (5.17) by MOSEK, we need to specify the components of the vector v that correspond to the constraints. For the constraints

$$ \begin{array}{@{}rcl@{}} {u_{i}^{2}} \leq 2 u_{2i} u_{2i+1} \qquad (i = 1, \ldots, 2^{p-1} - 1), \end{array} $$

we take the components v_2i, v_2(2i)− 1, and v_{2(2i+ 1)− 1} for u_i, u_2i, and u_2i+ 1, respectively. For the constraints

$$ \begin{array}{@{}rcl@{}} {u_{i}^{2}} \leq 2 \hat{g}_{2i - 2^{p} + 1} \hat{g}_{2i - 2^{p} + 2} \qquad \left\{\begin{array}{ll} (i = 2^{p-1} , \ldots, 2^{p-1} + \ell/2 + 1) & (\ell \text{ is even}), \\ (i = 2^{p-1} , \ldots, 2^{p-1} + (\ell-3)/2) & (\ell \text{ is odd}), \end{array}\right. \end{array} $$

we take the components v_2i, $v_{2^{p} + 2i - 1}$, and $v_{2^{p} + 2i}$ for u_i, $\hat {g}_{2i - 2^{p} + 1}$, and $\hat {g}_{2i - 2^{p} + 2}$, respectively. In the case that ℓ is odd, we take the components $v_{2^{p} + \ell - 1}$, $v_{2(2^{p} -1) + \ell }$, and v₁ for the constraint $u_{2^{p-1} + (\ell -1)/2}^{2} \leq 2 \hat {g}_{\ell } u_{1}$. Finally, for the constraints

$$ \begin{array}{@{}rcl@{}} z_{ij}^{2} \leq 2 t_{ij} \hat{w}_{i} \qquad (i=1,\ldots, m, \ j=1,\ldots, \ell), \end{array} $$

we take the components $v_{2(2^{p}-1) + \ell + ij}$, $v_{2(2^{p}-1) + \ell + ij + m\ell }$, and $v_{2(2^{p}-1) + \ell + ij + 2m\ell }$.