.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/basics/demo_optimizing_physics_parameter.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_basics_demo_optimizing_physics_parameter.py: Solving blind inverse problems / estimating physics parameters ============================================================== This demo shows you how to use :class:`deepinv.physics.Physics` together with automatic differentiation to optimize your operator. Consider the forward model .. math:: y = \noise{\forw{x, \theta}} where :math:`N` is the noise model, :math:`\forw{\cdot, \theta}` is the forward operator, and the goal is to learn the parameter :math:`\theta` (e.g., the filter in :class:`deepinv.physics.Blur`). In a typical blind inverse problem, given a measurement :math:`y`, we would like to recover both the underlying image :math:`x` and the operator parameter :math:`\theta`, resulting in a highly ill-posed inverse problem. In this example, we only focus on a much more simpler problem: given the measurement :math:`y` and the ground truth :math:`x`, find the parameter :math:`\theta`. This can be reformulated as the following optimization problem: .. math:: \min_{\theta} \frac{1}{2} \|\forw{x, \theta} - y \|^2 This problem can be addressed by first-order optimization if we can compute the gradient of the above function with respect to :math:`\theta`. The dependence between the operator :math:`A` and the parameter :math:`\theta` can be complicated. DeepInverse provides a wide range of physics operators, implemented as differentiable classes. We can leverage the automatic differentiation engine provided in Pytorch to compute the gradient of the above loss function w.r.t. the physics parameters :math:`\theta`. The purpose of this demo is to show how to use the physics classes in DeepInverse to estimate the physics parameters, together with the automatic differentiation. We show 3 different ways to do this: manually implementing the projected gradient descent algorithm, using a Pytorch optimizer and optimizing the physics as a usual neural network. .. GENERATED FROM PYTHON SOURCE LINES 34-36 Import required packages .. GENERATED FROM PYTHON SOURCE LINES 36-44 .. code-block:: Python import deepinv as dinv import torch from tqdm import tqdm import matplotlib.pyplot as plt device = dinv.utils.get_freer_gpu() if torch.cuda.is_available() else "cpu" dtype = torch.float32 .. GENERATED FROM PYTHON SOURCE LINES 45-50 Define the physics ------------------ In this first example, we use the convolution operator, defined in the :class:`deepinv.physics.Blur` class. We also generate a random convolution kernel of motion blur .. GENERATED FROM PYTHON SOURCE LINES 50-57 .. code-block:: Python generator = dinv.physics.generator.MotionBlurGenerator( psf_size=(25, 25), rng=torch.Generator(device), device=device ) true_kernel = generator.step(1, seed=123)["filter"] physics = dinv.physics.Blur(noise_model=dinv.physics.GaussianNoise(0.02), device=device) .. GENERATED FROM PYTHON SOURCE LINES 58-68 .. code-block:: Python x = dinv.utils.load_url_image( dinv.utils.demo.get_image_url("celeba_example.jpg"), img_size=256, resize_mode="resize", ).to(device) y = physics(x, filter=true_kernel) dinv.utils.plot([x, y, true_kernel], titles=["Sharp", "Blurry", "True kernel"]) .. image-sg:: /auto_examples/basics/images/sphx_glr_demo_optimizing_physics_parameter_001.png :alt: Sharp, Blurry, True kernel :srcset: /auto_examples/basics/images/sphx_glr_demo_optimizing_physics_parameter_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 69-79 Define an optimization algorithm -------------------------------- The convolution kernel lives in the simplex, ie the kernel must have positive entries summing to 1. We can use a simple optimization algorithm - Projected Gradient Descent - to enforce this constraint. The following function allows one to compute the orthogonal projection onto the simplex, by a sorting algorithm (Reference: `Large-scale Multiclass Support Vector Machine Training via Euclidean Projection onto the Simplex -- Mathieu Blondel, Akinori Fujino, and Naonori Ueda `_) .. GENERATED FROM PYTHON SOURCE LINES 79-102 .. code-block:: Python @torch.no_grad() def projection_simplex_sort(v: torch.Tensor) -> torch.Tensor: r""" Projects a tensor onto the simplex using a sorting algorithm. """ shape = v.shape B = shape[0] v = v.view(B, -1) n_features = v.size(1) u = torch.sort(v, descending=True, dim=-1).values cssv = torch.cumsum(u, dim=-1) - 1.0 ind = torch.arange(n_features, device=v.device)[None, :].expand(B, -1) + 1.0 cond = u - cssv / ind > 0 rho = ind[cond][-1] theta = cssv[cond][-1] / rho w = torch.maximum(v - theta, torch.zeros_like(v)) return w.reshape(shape) # We also define a data fidelity term data_fidelity = dinv.optim.L2() .. GENERATED FROM PYTHON SOURCE LINES 103-106 Run the algorithm Initialize a constant kernel .. GENERATED FROM PYTHON SOURCE LINES 107-137 .. code-block:: Python kernel_init = torch.zeros_like(true_kernel) kernel_init[..., 5:-5, 5:-5] = 1.0 kernel_init = projection_simplex_sort(kernel_init) n_iter = 1000 stepsize = 0.7 kernel_hat = kernel_init losses = [] for i in tqdm(range(n_iter)): # compute the gradient with torch.enable_grad(): kernel_hat.requires_grad_(True) physics.update(filter=kernel_hat) loss = data_fidelity(y=y, x=x, physics=physics) / y.numel() loss.backward() grad = kernel_hat.grad # gradient step and projection step with torch.no_grad(): kernel_hat = kernel_hat - stepsize * grad kernel_hat = projection_simplex_sort(kernel_hat) losses.append(loss.item()) dinv.utils.plot( [true_kernel, kernel_init, kernel_hat], titles=["True kernel", "Init. kernel", "Estimated kernel"], suptitle="Result with Projected Gradient Descent", ) .. image-sg:: /auto_examples/basics/images/sphx_glr_demo_optimizing_physics_parameter_002.png :alt: Result with Projected Gradient Descent, True kernel, Init. kernel, Estimated kernel :srcset: /auto_examples/basics/images/sphx_glr_demo_optimizing_physics_parameter_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none 0%| | 0/1000 [00:00` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: demo_optimizing_physics_parameter.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: demo_optimizing_physics_parameter.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_