from typing import TYPE_CHECKING
from pylops.utils.typing import NDArray
from pyproximal.ProxOperator import ProxOperator, _check_tau
if TYPE_CHECKING:
from pylops.linearoperator import LinearOperator
[docs]
class Orthogonal(ProxOperator):
r"""Proximal operator of any function of the product between an orthogonal
matrix and a vector (plus summation of a vector).
Proximal operator of any quadratic function
:math:`g(\mathbf{x})=f(\mathbf{Qx} + \mathbf{b})` where :math:`\mathbf{Q}`
is an orthogonal operator and :math:`\mathbf{b}` is a vector in the data
space of :math:`\mathbf{Q}`.
Parameters
----------
f : :obj:`pyproximal.ProxOperator`
Proximal operator
Q : :obj:`pylops.LinearOperator`
Orthogonal operator
partial : :obj:`bool`, optional
Partial (``True``) of full (``False``) orthogonality
b : :obj:`numpy.ndarray`, optional
Vector
alpha : :obj:`float`, optional
Positive coefficient for partial orthogonality. It will be ignored if
``partial=False``
Notes
-----
The Proximal operator of any function of the form
:math:`g(\mathbf{x}) = f(\mathbf{Qx} + \mathbf{b})` with the operator
:math:`\mathbf{Q}` satisfying the following condition
:math:`\mathbf{Q}\mathbf{Q}^T = \alpha \mathbf{I}` (partial orthogonality)
is [1]_:
.. math::
\prox_{\tau g}(x) = \frac{1}{\alpha} ((\alpha \mathbf{I} -
\mathbf{Q}^H \mathbf{Q}) \mathbf{x} + \mathbf{Q}^H
(\prox_{\alpha \tau f}(\mathbf{Qx} + \mathbf{b}) - \mathbf{b}))
A special case arises when :math:`\mathbf{Q}\mathbf{Q}^T =
\mathbf{Q}^T\mathbf{Q} = \mathbf{I}`
(full orthogonality), and the proximal operator reduces to:
.. math::
\prox_{\tau g}(x) = \mathbf{Q}^H (\prox_{\tau f}(\mathbf{Qx} +
\mathbf{b}) - \mathbf{b}))
.. [1] Daniel O'Connor, D., and Vandenberghe, L., "Primal-Dual
Decomposition by Operator Splitting and Applications to Image
Deblurring", SIAM J. Imaging Sciences, vol. 7, pp. 1724–1754. 2014.
"""
def __init__(
self,
f: ProxOperator,
Q: "LinearOperator",
partial: bool = False,
b: NDArray | None = None,
alpha: float = 1.0,
) -> None:
super().__init__(None, False)
self.f = f
self.Q = Q
self.partial = partial
self.alpha = alpha
self.b = b if b is not None else 0.0
def __call__(self, x: NDArray) -> bool | float | int:
y = self.Q.matvec(x)
y += self.b
f = self.f(y)
return f
@_check_tau
def prox(self, x: NDArray, tau: float) -> NDArray:
y = self.Q.matvec(x)
if self.partial:
z = (1.0 / self.alpha) * (
self.alpha * x
- self.Q.rmatvec(y)
+ self.Q.rmatvec(self.f.prox(y + self.b, self.alpha * tau) - self.b)
)
else:
y = y + self.b
z = self.f.prox(y, tau)
z = z - self.b
z = self.Q.rmatvec(z)
return z
PyProximal