🎉 Add supplement materials

appendix/README.md

@@ -0,0 +1,7 @@
+# Supplement materials
+
+[:page_facing_up: Appendix A.](appx_A) A simple example of FracGM with global optimal guarantees
+
+[:page_facing_up: Appendix B.](appx_B) Tightness of linear relaxation
+
+[:page_facing_up: Appendix C.](appx_C) Sensitivity of initial guess

appendix/appx_A/.gitignore

@@ -0,0 +1,2 @@
+# Byte-compiled / optimized / DLL files
+__pycache__/

appendix/appx_A/README.md

@@ -0,0 +1,92 @@
+### A simple example of FracGM with global optimal guarantees
+
+Here we give a simple example that one can verify the differentiability and Lipschitz continuity of $\psi(\boldsymbol{\alpha},x_{\boldsymbol{\alpha}})$. Suppose we have an optimization problem:
+
+$$
+    \min_x\ \frac{f(x)}{h(x)}=\frac{x^2}{x^2+1},
+$$
+
+where the residual $r(x)=x$ and $c=1$ in sense of Geman-McClure. The dual problem $(6)$ is as follows:
+
+$$
+    \min_x\ \mu\big(f(x)-\beta h(x)\big)=\mu\big(x^2-\beta(x^2+1)\big)=\mu x^2-\mu\beta x^2-\mu\beta,
+$$
+
+where $\beta\in\mathbb{R}$ and $\mu\in\mathbb{R}$. Let $\boldsymbol{\alpha}=(\beta,\mu)^\top\in\mathbb{R}^2$ and $x_{\boldsymbol{\alpha}}$ be the optimal solution of Problem $(6)$ given $\boldsymbol{\alpha}$, we obtain $x_{\boldsymbol{\alpha}}=0$ by setting the derivative to zero as follows:
+
+$$
+    \frac{\partial}{\partial x}(\mu x^2-\mu\beta x^2-\mu\beta)=2\mu(1-\beta)x=0\implies x^*=x_{\boldsymbol{\alpha}}=0.
+$$
+
+Thereby we can write 
+
+$$
+    \psi(\boldsymbol{\alpha},x_{\boldsymbol{\alpha}})=
+    \begin{pmatrix}
+    -f(x_{\boldsymbol{\alpha}})+\beta h(x_{\boldsymbol{\alpha}})\\
+    -1+\mu h(x_{\boldsymbol{\alpha}})
+    \end{pmatrix}=
+    \begin{pmatrix}
+    -x_{\boldsymbol{\alpha}}^2+\beta x_{\boldsymbol{\alpha}}^2+\beta\\
+    -1+\mu x_{\boldsymbol{\alpha}}^2+\mu
+    \end{pmatrix}=
+    \begin{pmatrix}
+    \beta\\
+    \mu-1
+    \end{pmatrix}.
+$$
+
+It is trivial that both component functions $\psi_1(\boldsymbol{\alpha},x_{\boldsymbol{\alpha}})=\beta$ and $\psi_2(\boldsymbol{\alpha},x_{\boldsymbol{\alpha}})=\mu-1$ are differentiable, thus $\psi(\boldsymbol{\alpha},x_{\boldsymbol{\alpha}})$ is differentiable. Given $\boldsymbol{\alpha}_1=(\beta_1,\mu_1)^\top$ and $\boldsymbol{\alpha}_2=(\beta_2,\mu_2)^\top$,
+
+$$
+    \\|\psi(\boldsymbol{\alpha}\_1,x_{\boldsymbol{\alpha}\_1})-\psi(\boldsymbol{\alpha}\_2,x_{\boldsymbol{\alpha}_2})\\|=
+    \begin{Vmatrix}
+    \beta_1-\beta_2\\\\\mu_1-1-(\mu_2-1)
+    \end{Vmatrix}\leq\begin{Vmatrix}
+    \beta_1-\beta_2\\\\\mu_1-\mu_2
+    \end{Vmatrix},
+$$
+
+then $\psi(\boldsymbol{\alpha},x_{\boldsymbol{\alpha}})$ is Lipschitz continuous. Solving such (simple) case by FracGM guarantees that the solution is global optimal.
+
+
+### Usage
+```
+cd appendix/appx_A
+python ./main.py
+```
+
+We test various random initial guess between -100000 and 100000, all solutions converge to 0.
+
+```
+=========================
+initial guess: -90760
+solution: 0.0
+=========================
+initial guess: 90371
+solution: 0.0
+=========================
+initial guess: 71362
+solution: 0.0
+=========================
+initial guess: 25425
+solution: 0.0
+=========================
+initial guess: -82805
+solution: 0.0
+=========================
+initial guess: 42735
+solution: 0.0
+=========================
+initial guess: -174
+solution: 0.0
+=========================
+initial guess: -12185
+solution: 0.0
+=========================
+initial guess: 33929
+solution: 0.0
+=========================
+initial guess: -23746
+solution: 0.0
+```

appendix/appx_A/fracgm.py

@@ -0,0 +1,63 @@
+import numpy as np
+from scipy.optimize import minimize_scalar
+
+
+class FracGM:
+    """
+    FracGM solver for the optimization problem:
+        min x^2 / (x^2 + 1),
+    where the optimal solution is x* = 0.
+    """
+
+    def __init__(self, max_iteration=1000, tolerance=1e-7, verbose=True):
+        self.max_iter = max_iteration
+        self.tol = tolerance
+        self.verbose = verbose
+
+    @staticmethod
+    def f(x):
+        return x**2
+
+    @staticmethod
+    def h(x):
+        return x**2 + 1
+
+    def solve_beta_mu(self, x):
+        beta = self.f(x) / self.h(x)
+        mu = 1 / self.h(x)
+
+        return beta, mu
+
+    def solve_x(self, beta, mu):
+        def dual_problem(x):
+            return mu * (self.f(x) - beta * self.h(x))
+
+        # Scipy minimizer for scalar function of one variable
+        res = minimize_scalar(dual_problem)
+
+        return res.x
+
+    def compute_psi_norm(self, beta, mu, x):
+        psi_1 = -self.f(x) + beta * self.h(x)
+        psi_2 = -1 + mu * self.h(x)
+
+        return np.linalg.norm([psi_1, psi_2])
+
+    def solve(self, initial):
+        if self.verbose:
+            print("initial guess:", initial)
+
+        x = initial
+
+        for i in range(self.max_iter):
+            beta, mu = self.solve_beta_mu(x)
+            x = self.solve_x(beta, mu)
+
+            norm = self.compute_psi_norm(beta, mu, x)
+            if norm < self.tol:
+                break
+
+        if self.verbose:
+            print("solution:", x)
+
+        return x

appendix/appx_A/main.py

@@ -0,0 +1,9 @@
+import random
+from fracgm import FracGM
+
+# random initial guess
+initial_list = random.sample(range(-100000, 100000), 10)
+
+for initial_guess in initial_list:
+    print("=" * 25)
+    FracGM().solve(initial=initial_guess)

appendix/appx_B/README.md

@@ -0,0 +1,30 @@
+### Tightness of linear relaxation
+
+For FracGM-based rotation and registration solvers, the relaxation is tight if solutions of the relaxed program in $\mathbb{R}^{3\times 3}$ and $\mathbb{R}^{4\times 4}$ are also within $\text{SO}(3)$ and $\text{SE}(3)$ respectively.
+In practice, we can evaluate the tightness of the relaxation by checking the orthogonality and determinant of a solution given by the relaxed program.
+Taking our rotation estimation tasks on synthetic dataset for example, we quantify the tightness by an orthogonal loss $\mathcal{L}_o$ and a determinant loss $\mathcal{L}_d$ defined as follows:
+
+$$
+\begin{aligned}
+\mathcal{L}_o&=\|(\mathbf{R}^\prime)^\top\mathbf{R}^\prime - \mathbf{I}_3\|_F\\
+\mathcal{L}_d&=|\text{det}(\mathbf{R}^\prime)-1|
+\end{aligned}
+$$
+
+where $\|\cdot\|_F$ is the Frobenius norm, and we verify the orthogonality and determinant of the rotation part of solutions $\mathbf{R}^\prime$ before SVD projection:
+
+![](./docs/rotation_loss.png)
+
+Similarly, we verify the rotation part of solutions for registration tasks on synthetic dataset as follows:
+
+![](./docs/registration_loss.png)
+
+Above evaluations reveal that FracGM-based rotation and registration solvers are empirically tight mostly in our experiments.
+We conjecture that the solution tends to approach the original solution space with sufficient number of measurements, in which a similar discussion can be found in [[34](#ref1), Section 2.1].
+Though, unfortunately we are unable to clarify the conditions to theoretically guarantee tightness and global optimality for both FracGM-based rotation and registration solvers.
+Besides the issue of tightness, the global optimality of the relaxed problem is also subject to the conditions specified in Proposition 3, which is difficult to verify as it depends on the input data and requires discussions on various situations.
+We appreciate the reviewer for raising this important issue, which we also consider to be one of the difficulties but crucial research problems in our future work.
+
+
+<a id="ref1">[34]</a> 
+C. Olsson and A. P. Eriksson, “Solving quadratically constrained geometrical problems using lagrangian duality,” in ICPR. IEEE Computer Society, 2008, pp. 1–5.

appendix/appx_C/README.md

@@ -0,0 +1,11 @@
+### Sensitivity of initial guess
+
+If a FracGM method satisfies Proposition 3, then it is a global solver, and it should not be sensitive to initial guesses.
+We do not think the FracGM-based rotation and registration solvers are insensitive to initial guesses for any kinds of input data, due to the fact that we are unable to verify Proposition 3 at the moment. Nevertheless, empirical studies on experiments show that our FracGM-based solvers are mostly insensitive to initial guesses.
+
+![](./docs/perturbation.png)
+
+We compare our FracGM-based rotation solver with the other Geman-McClure solver, GNC-GM [[17]](#ref1), on the synthetic dataset. We construct various initial guesses by adding some degree of perturbation to the ground truth rotation matrix. With an outlier rate of 50\%, both FracGM and GNC-GM seems insensitive to initial guesses. However, in extreme cases where the outlier rate is 90\%, GNC-GM produces different solutions due to different initial guesses, while FracGM is still insensitive to initial guesses.
+
+<a id="ref1">[17]</a> 
+H. Yang, P. Antonante, V. Tzoumas, and L. Carlone, “Graduated nonconvexity for robust spatial perception: From non-minimal solvers to global outlier rejection,” IEEE Robotics Autom. Lett., vol. 5, no. 2, pp. 1127–1134, 2020.