# Single-Image based unsupervised joint segmentation and denoising

Nadja Gruber , , Corresponding author<EMAIL_ADDRESS>of
Mathematics, University of Innsbruck, AustriaVASCage-Research Centre on
Vascular Ageing and Stroke, Innsbruck, Austria Johannes Schwab MRC Laboratory
of Molecular Biology, Cambridge, UK Noémie Debroux Institut Pascal, Université
Clermont Auvergne, Clermont-Ferrand, France Nicolas Papadakis Institut de
Mathématiques de Bordeaux, Bordeaux, France Markus Haltmeier 22footnotemark: 2

###### Abstract

In this work, we develop an unsupervised method for the joint segmentation and
denoising of a single image. To this end, we combine the advantages of a
variational segmentation method with the power of a self-supervised, single-
image based deep learning approach. One major strength of our method lies in
the fact, that in contrast to data-driven methods, where huge amounts of
labeled samples are necessary, our model can segment an image into multiple
meaningful regions without any training database. Further, we introduce a
novel energy functional in which denoising and segmentation are coupled in a
way that both tasks benefit from each other. The limitations of existing
single-image based variational segmentation methods, which are not capable of
dealing with high noise or generic texture, are tackled by this specific
combination with self-supervised image denoising. We propose a unified
optimisation strategy and show that, especially for very noisy images
available in microscopy, our proposed joint approach outperforms its
sequential counterpart as well as alternative methods focused purely on
denoising or segmentation. Another comparison is conducted with a supervised
deep learning approach designed for the same application, highlighting the
good performance of our approach.

##  1 Introduction

Image denoising and segmentation are fundamental problems in image processing
[37, 5, 23]. In many biomedical applications, such as fluorescence microscopy
or transmission electron cryomicroscopy, one is interested in the segmentation
of objects. However, training data for this task is typically scarce and hard
to obtain due to the intrinsic complexity and high noise of such images as
well as the long time required by experts to label them. Therefore, there is a
need for unsupervised methods for tackling the two imaging tasks in a unified
way. In this work, we propose such a framework, and apply it to a subset of a
popular, public available dataset of microscopy images.

The objective of segmentation is to divide a given image into different,
meaningful regions, while denoising describes the task of removing noise from
a corrupted image. The main difficulty in noise removal is to flatten the
unwanted, high frequency corruption, while preserving essential features such
as edges. At first glance, denoising and segmentation are two different
applications. Nevertheless, both tasks share relationships, as very similar
models can be used to solve both problems [8]. As we demonstrate in this work,
denoising and segmentation can benefit a lot from each other. By identifying
edges segmentation guides the denoising process to preserve sharp structures
while smoothing the unwanted high frequency residuals. Also, by removing
unnecessary and misleading information from images, denoising helps and
improves the segmentation accuracy.

There exist at least two main kinds of approaches to tackle the two tasks
individually. The first class of methods involves the minimisation of an
energy functional within graph or variational frameworks. The second type of
approaches that recently became popular considers deep learning techniques,
especially convolutional neural networks [29]. In the following, we give a
short overview of the most important and related variational and deep learning
based methods.

###  1.1 Variational Methods

Standard imaging methods for segmentation and denoising are based on an energy
functional that captures the desired characteristics of the output image. The
energy functional typically consists of a data fitting term, and a
regularisation term that encourages properties of the output image, such as
smoothness or sparsity. The energy functional is then minimised using
optimisation techniques such as gradient descent or proximal splitting
algorithms.

##### Denoising

One of the best known variational model for image denoising is the Rudin-
Osher-Fatemi (ROF) [36] model. This model improves the region-based Mumford-
Shah [31] functional that realises a piecewise smooth approximation of an
input image. The ROF model and its extensions reduce noise by penalizing the
total variation of the image. Such methods thus promote piecewise constant
images with undesirable staircase effects in homogeneous regions and they are
unable to recover image details and patterns with higher variation. In case of
severe input noise, they provide poor denoising results as image contours are
confused with noise [16].

On the other hand, the resulting image being piecewise constant, it can be
used for segmentation, by choosing regions of the same value, or thresholding
the image. More details about the link between the ROF based denoising models
and segmentation can for example be found in [8].

##### Segmentation

In their seminal paper [11], Chan and Vese proposed to solve the MumfordShah
problem with a levelset reformulation. Let us denote by $\Omega$ a bounded
subset of $\mathbb{R}^{2}$ with Lipschitz boundary, where the given image
$f:\Omega\rightarrow[0,1]$ is defined, and by
$u\colon\Omega\rightarrow\\{0,1\\}$ the desired binary mask, separating $f$
into two different areas corresponding to the two mean intensity values
$c_{1}$ and $c_{2}$, belonging to foreground and background region,
respectively. In 2D, this framework involves a surface $\phi$ whose zero level
represents the contour of interest, and the mask is obtained as
$u(x)=H\left(\phi(x)\right),$ where $H(\cdot)$ is the Heaviside function. The
proposed energy for binary segmentation is given by

$\displaystyle\mathcal{E}(\phi,c_{1},c_{2})=\int_{\Omega}\lvert\nabla
H(\phi(x))\rvert dx$ $\displaystyle+\lambda\int_{\Omega}\lvert
f(x)-c_{1}\rvert^{2}H(\phi(x))dx$ (1)
$\displaystyle+\lambda\int_{\Omega}\lvert
f(x)-c_{2}\rvert^{2}(1-H(\phi(x))dx,$

where $\lambda>0$ is a regularization parameter to tune. Slight modifications
of this method have already been used in microscopy [34, 40], as it is well
adapted to cell segmentation, where labeled data are scarce. which computes
the intensity averages by using constant intensity information across the
region. However, the Chan-Vese model computes the intensity averages by using
constant information across the region, and thus is a global region-based
model. Therefore, it does not deal well with intensity inhomogeneities and the
presence of high noise. To mitigate this problem, many local region-based
extensions of the active contour of the piecewise constant active contour
model have been proposed [27, 45], but the methods remain sensible to the
considered hand-crafted features and the initial contour. In another line of
works, pre-filtering tools are considered to better prepare the image for
segmentation [9, 28, 43] in a sequential pipeline: denoise before segment. In
[7], a three-stage approach is proposed, consisting of smoothing, lifting, and
segmentation using thresholding.

In our work, we tackle the afore-mentioned issues by introducing a generalized
Chan-Vese segmentation functional including a robust data-fidelity term that
is jointly learned with self supervised deep denoising techniques.

###  1.2 Deep Learning Methods

We now review some of the most relevant works using neural networks for the
denoising and segmentation of images.

##### Denoising

While variational image denoising techniques focus on explicitly modeling data
noise, modern deep learning approaches directly learn how to map the noisy
image to its clean counterpart. In the literature, several types of deep
learning based denoising methods can be found. In particular, supervised
approaches require pairs of noisy images and corresponding clean ground truth
data (see [21, 44, 41]). However, the presence of such noisy and clean data
pairs is rare in reality and often artificial, and it therefore makes such
methods useless for practical applications.

To overcome the requirement of clean images, internal statistical methods
(i.e. methods, where image patches from the same image are used for the noise
reduction), have been introduced [46]. In [39], Ulyanov et al. exploit the
fact that the internal structure of CNNs, inherently resonates with the
distribution of natural images, and utilize this observation for image
restoration without the need of additional training data. For each single
image to restore, this method thus proposes to train a CNN to reconstruct the
considered image. The idea is that an early stopping of the training allows to
recover a regularised, denoised image. A different strategy is proposed in the
Noise2Noise [25] method, where noisy image pairs are mapped to one another.
The drawback of these type of methods is that it still relies on the
availability of such pairs. In practice, even the acquisition of two noisy
realisations of the same image content is often difficult [5]. To this end,
self-supervised training methods operating on one single noisy image, such as
Noise2Void [24], Noise2Self [3], and more recently Noise2Fast [26], are
promising alternatives. This self-supervision is accomplished by
excluding/masking the center (blind spot) of the receptive field of the
network. In this type of training, it is assumed that the noise is pixelwise
independent and that the true intensity of a pixel can be predicted from the
local image context, with the exception of the blind spots mentioned
previously [24]. In this work, We utilised Noise2Fast [26] in our study due to
its optimal combination of computational speed and performance. The method
itself will be explained in more detail in Sections 2, and 4.

##### Segmentation

Among existing deep learning based approaches addressing image segmentation,
the U-Net [35] first introduced for microscopy cell segmentation, is one of
the most successful network architectures. Next to it, we mention Mask-RCNN
[20], a two-stage object detection and segmentation, extending the popular
faster R-CNN architecture [17]. DeepLab [13] is a family of methods that use
atrous convolution (also known as dilated convolution) to capture multi-scale
context information. It has been shown to achieve state-of-the-art results on
several segmentation benchmarks. Still, even the best existing methods offer
plenty of scope for improvements, motivating further research in this field
[21,23,7]. Their high performance comes along with a price to pay. A common
trait to the mentioned approaches is their requirement for tremendous amounts
of labeled ground truth training data, the creation of which is time-consuming
and prone to subjective errors.

As already mentioned, for many applications such as microscopy, the available
image data is very noisy and the available ground truth training data are
scarce. It is thus of great interest to tackle both, the segmentation and the
denoising, in a unified manner. In the following, we will review variational
methods, as well as deep learning based approaches tackling both, the
segmentation and the denoising in a joint manner.

###  1.3 Joint denoising and segmentation methods

In [6], Cai et al. design a model tackling the segmentation of images with a
high level of noise or blurriness. To this end, they propose a variational
approach, coupling an extension of the piecewise constant Mumford Shah model
with an image restoration model, making it more robust in the processing of
the given corrupted image $f$. In [14], the authors propose a variational
approach for the joint reconstruction and segmentation. Therefore, they derive
a model consisting of a total variation regularised reconstruction from
undersampled data, and a Chan-Vese based segmentation. The authors show the
improvement of joint reconstruction and segmentation performance compared to
the sequential approach. In another work [33], Ramlau et. al illustrate that
the Mumford–Shah level-et method can enhance the quality of reconstructed
images and improve the accuracy of segmentation results.

In the context of microscopy data, purely deep learning based approaches
dealing with both segmentation and denoising are [32] and [5]. In [32],
Prakesh et al. demonstrated on various microscopy datasets that the use of
self-supervised denoising priors improves the segmentation results, especially
when only a few ground truth segmentation masks are available for training. In
a similar work [5], the authors propose DenoiSeg, consisting of a U-Net for
image segmentation, and the self-supervised denoising scheme Noise2Void [24],
which are combined and trained with a common loss. The authors demonstrate,
that the global optimisation outperforms the sequential counterpart, where the
image is denoised first and then segmented. This method requires labeled data.
To reach high denoising performance, a huge amount of noisy data is also
required. Moreover, the loss function is just the sum of the segmentation and
denoising losses. There exists no coupling between the two tasks in the
objective to optimise. Segmentation can therefore benefit from the noise
reduction, but the reverse is not possible.

To overcome these limitations, we a propose single-image method with a new
joint loss that allows full interaction between segmentation and noise
reduction.

###  1.4 Contributions

In this work, we propose a new model for joint image denoising and
segmentation that combines advantages from variational models and deep
learning. In contrast to the aforementioned state-of-the-art deep learning
based methods which require a large cohort of labeled and clean training
images, we obtain comparable results using only one single image. While
supervised deep learning based methods are trained with hand-made annotations
which are prone to subjective labelling errors, the proposed combination of a
variational segmentation model with a self-supervised denoising CNN does not
require any labeled data or a representative dataset leading to the
elimination of pre-training.

We combine denoising and segmentation tasks in such a way, that both of them
benefit from each other. This is a main difference between existing deep joint
approaches such as [5], where the denoising task solely aims at improving
segmentation. More specifically, we design two dedicated denoising networks
for the foreground and background regions to improve the overall image
denoising performance, and use the difference of the denoising performances in
the two regions to find the segmentation mask.

Our method can be seen as a flexible generalization of existing Chan-Vese
models. Standard Chan-Vese models as well as joint variational methods for
segmentation and denoising [6, 14] rely on the piecewise constant assumption
of hand-crafted features. Thus, these methods are struggling with intensity
inhomogeneties, textures or high noise levels. Further, methods of that kind
strongly depend on the initialisation due to the non-convexity of the
functional. In this work we propose to learn the structural patterns of
different regions in the image without any prior information, by designing a
novel energy functional, where feature information is captured in a natural
way by a denoising network. The paper is organised as follows. We start with
toy examples that illustrate the region specific denoisers as the motivation
of our method in 1.5. We then formulate the two problems that we aim to solve
and revise the necessary ingredients that make up our algorithm in Section 2.
Our proposed new algorithm is described and analysed in Section 3, and the
numerical implementation is presented in 4. Section 5 shows the application of
the method to microscopy data. Further, we apply our proposed model to natural
images, and demonstrate, that with manual guidance by the user roughly
indicating the two regions, our method can successfully be applied to more
complex problems. The paper ends with a conclusion and outlook to possible
future work.

###  1.5 Motivation

In the following, we give an intuition of how the two denoising and
segmentation tasks we aim to solve are coupled in a way that both of them have
a positive influence on each other. We present some toy examples, showing how
segmentation can successfully guide the denoising process. In a toy example,
we generate an image of size 256$\times 256$ that consists of stripe patterns
that are aligned differently in different areas; see Figure 1. This image is
further corrupted by manually adding Gaussian noise with a noise level of 50.
Here, we use two linear neural networks (2 ”experts” respectively dedicated to
the foreground and the background) consisting of a single convolutional layer
with one filter of size $15\times 15$, which is trained using a slight
modification of the Noise2Fast training strategy [26] described in Section 4.
More precisely, we restrict the training of the network to the two regions of
the image, respectively, by masking the loss function, and restricting the
training on boxes of size $30\times 30$ which are depicted in Figure 1. We
find that learned filters are adapted to the structure of the signal in the
corresponding region. As a result, the error patterns have higher values in
the region where the denoiser has not been trained. This provides the basis
for exploiting region-specific denoising for segmentation. The experimental
details for this toy example are provided in Section 5.

Figure 1: Visualisation of the idea behind the proposed joint denoising and
segmentation model. Here, we trained two networks consisting of one single
filter using the Noise2Fast [25] strategy and restricted the training to the
two boxes marked in the noisy image $f$. From the two right binary images in
the bottom row, we observe that the two denoising experts perform much better
in the region they have been trained on. The difference images (noisy image
minus Denoised by Expert 1 (resp. 2) can then be used in the segmentation
process, by exploiting the fact that regions with a small denoising error for
the first (resp. second) expert can be assigned as foreground (resp.
background).  Figure 2: Given noisy RGB input image (corrupted with Gaussian
noise, noise level = 0.75), denoised image using Noise2Fast on the whole
image, region-specific experts, and ground truth image. We clearly observe
sharper edges, and better recovered color information in the “two-
experts”-example.

The positive effect of segmentation on the denoising process is even more
evident in the natural image shown in Figure 2. We used the network
architecture proposed by the authors in [26], resulting in a non-linear neural
network. First, the denoising network was trained and subsequently applied on
the whole image. The second image shows the result obtained with two
separately trained neural networks. This strategy yields a better visual
result, which is further confirmed by PSNR values of 19.69 and 19.12,
respectively. The noisy image has been generated by scaling the given clean
RGB input image to $[0,1]$, and adding randomly distributed noise scaled with
the maximum pixel value, ensuring that the noise is proportional to the image
intensity. Here, we used a manually generated mask of the zebra and background
region, and during training, computed the MSE restricted to the two different
regions, respectively.

In the next section, we will fix the notation, formalize the problem, and
describe the main ingredients that are used for our proposed method.

##  2 Problem Description

We now present in detail the background that is relevant for our algorithm.
First we set our notations, and describe our proposed energy functional for
our unified denoising and segmentation framework.

In the following, we denote by $\Omega\subset\mathbb{R}^{2}$ a bounded set
with Lipschitz boundary, and by $\mathbb{F}$ a space of functions
$f\colon\Omega\rightarrow\mathbb{R}^{d}$, with $d=1$ for grayscale images, and
$d=3$ in the RGB case. We consider a given (noisy) image $f\in\mathbb{F}$,
which we want to jointly denoise, and split up into $C$ different regions.

######  Problem 1 (Image Denoising).

The goal of image denoising is to recover a clean image $g$ from a noisy
observation $f$ which follows an image degradation model $f=g+n$, where $n$ is
the signal degrading noise which we want to remove.

Note that although other degradation types are possible, we assume an additive
model here and specifically we will consider noise with an expected value of
zero.

######  Problem 2 (Image Segmentation).

Image segmentation refers to the process of automatically dividing an image
into meaningful regions. Based on specific pre-defined characteristics of a
given image $f\in\mathbb{F},$ one is interested in splitting the image domain
into two (in the case of binary segmentation) regions $\Sigma,$ and
$\Omega\setminus\Sigma$. In the case of multiclass segmentation, the objective
is to build a partition $\Omega=\bigcup_{i=1}^{C}\Sigma_{i}$ of the image
domain into $C$ disjoint regions (classes), where each of the regions
$\Sigma_{1},\dots,\Sigma_{C-1}$ represents a specific structure of objects in
$f$ and
$\Omega\setminus(\Sigma_{1}\uplus\Sigma_{2}\uplus\dots\uplus\Sigma_{C-1})$
represents the background.

In this work, we address these two problems simultaneously by designing an
energy functional in a way that both tasks benefit from each other. Next, we
discuss the two main components from the literature that form the basis of our
approach.

###  2.1 Convex Chan-Vese Formulation

In [12], Chan et al propose to relax the binary Chan-Vese segmentation problem
(1) and let the desired solution $u(x)$ take values in $[0,1]$. The resulting
convex energy is

$\displaystyle\min_{0\leq u\leq 1}\int_{\Omega}\lvert\nabla
u\rvert+\lambda\int_{\Omega}\left((c_{1}-f(x))^{2}-(c_{2}-f(x))^{2}\right)u(x)dx.$
(2)

The authors showed that, for any fixed constants $c_{1},c_{2}\in\mathbb{R},$ a
global minimiser for the non-convex problem can be found by carrying out the
minimisation in (2), and setting $\Sigma=\\{x:u(x)>\tau\\}$ for a.e.
$\tau\in[0,1]$.

Though the model is convex it still suffers from difficulties in segmenting
images where the piecewise constant assumption is not a relevant prior for the
different regions in the image, or if the image is corrupted by severe noise.
These issues are the main problems to solve in the current paper.

###  2.2 Self-supervised single-image based denoising

For a given noisy image $f=g+n$, self-supervised denoising methods are based
on some variant of the self supervised loss

$\displaystyle\mathcal{L}_{f}(\theta)=\int_{\Omega}(\Phi_{\theta}(f)(x)-f(x))^{2}.$
(3)

Clearly, such a strategy cannot work without restricting the class of
functions $\Phi_{\theta}$, since $\Phi_{\theta}=Id$ would be the minimiser and
does not yield a denoised image. One strategy to overcome this problem is the
method introduced in [39] where a generative model $\Phi_{\theta}$ is trained
minimising (3). In this framework, the convolutional structure and early
stopping prevents $\Phi_{\theta}$ to learn the fine image features (noise) to
obtain a denoised image. Another strategy is linear filtering with a
restriction on the filter [3, 24]. For example, a filter which is zero in its
central position, and therefore not taking into account the information of
this pixel but only the surrounding areas can be used to denoise an image
minimising (3). Another type of method tries to minimise a slightly different
functional. Motivated by the previous example, the authors of [26] introduce
$N$ random binary masks $\mathcal{H}_{k}$ that delete information in the
image. Training is then done using the loss function

$\displaystyle\mathcal{L}_{f}(\theta)=\frac{1}{N}\sum_{k=1}^{N}\int(\Phi_{\theta}(\mathcal{H}_{k}\cdot
f)(x)-f(x))^{2}\cdot(1-\mathcal{H}_{k}).$ (4)

This training strategy also prevents the network from learning the identity
operator. Although not directly minimising (3), we use a variant of this
method named Noise2Fast [26]. This variant uses regular masks
$\mathcal{H}_{1},\mathcal{H}_{2},\mathcal{H}_{3},\mathcal{H}_{4}$, which
consist of horizontal and vertical stripes on the even and odd image indices.

##  3 Proposed Joint Denoising and Segmentation

We now introduce our joint model, inspired by the observations described in
Section 1.5. To control binary segmentation, we propose to train two denoising
neural networks, each focusing on performing well in one of the regions to be
segmented (cf. Figure 1 and 2). We denote these “experts” by
$\Phi_{\theta^{F}}$ for the foreground, and $\Phi_{\theta^{B}}$ for the
background. These experts are neural networks with parameters $\theta^{F}$ and
$\theta^{B}$, that are trained with a modified denoising strategy. Let us
mention that the model is presented in the case of two regions, but the
extension to multi-class is straightforward, following for instance the
framework in [19, 1, 30]

In Section 3.1, we present the proposed joint energy function designed for the
combined denoising and segmentation process. This energy generalizes the
convex Chan-Vese functional (2) with a data-fidelity term defined from the
self-supervised denoising method . The optimisation scheme is performed in an
alternating way, as presented in Section 3.2. We finally provide theoretical
convergence results for our algorithm in Section 3.3.

###  3.1 Joint energy functional

In the following, we denote by $BV(\Omega)$ the space of all integrable
functions $u:\Omega\rightarrow\mathbb{R}$ with bounded total variation
$|u|_{\text{TV}}$, and consider the admissible set

$\displaystyle\mathbb{A}\coloneqq\\{u\in BV(\Omega)\mid 0\leq u\leq 1\\}.$

Further, let $i_{\mathbb{A}}:BV(\Omega)\rightarrow[0,\infty]$ denote the
associated indicator function, which is 0 inside $\mathbb{A}$, and $\infty$
elsewhere. The parameters of the two denoising experts, $\Phi_{\theta^{F}}$
and $\Phi_{\theta^{B}}$ are denoted by
$\boldsymbol{\theta}=(\theta^{F},\theta^{B})\in\mathbb{R}^{L\times L},$ and
are respectively dedicated to the foreground and the background. These two
experts are neural networks trained using the strategy proposed in [26]. We
consider the joint model

$\begin{split}\mathcal{E}_{f,\lambda}(u,\boldsymbol{\theta})=i_{\mathbb{A}}(u)+\lambda\lvert
u\rvert_{\text{TV}}&+\int_{\Omega}\left(f(x)-\Phi_{\theta^{F}}(f)(x)\right)^{2}u(x)dx\\\
&+\int_{\Omega}\left(f(x)-\Phi_{\theta^{B}}(f)(x)\right)^{2}(1-u(x))dx\,.\end{split}$
(5)

Note that for fixed network parameters $\boldsymbol{\theta}$, the proposed
energy is convex in $u$. Moreover, we can threshold the result and still have
a global optima (see Theorem 3.2). Further we point out that in the case where
the noise2fast training strategy is used, the energy functional for the
denoising step is not exactly functional (5).

Figure 3 illustrates the idea behind the proposed segmentation model. For
grayscale images, one can initialise the algorithm by thresholding image
values. In more complex cases, a user can be asked to provide representative
boxes for the background and foreground regions. Then, alternately, the
denoising experts are trained on subsets of the two different segmented
regions and the segmentations are updated. In practice, the data-fidelity term
in (5) is updated given the denoising performance of the two experts
$\Phi_{\theta^{F}}$ and $\Phi_{\theta^{B}}$. For fixed network parameters
$\boldsymbol{\theta}$, the energy (5) is minimised. Repeating this procedure
until a convergence criteria is met, we obtain the segmentation mask $u$, as
well as the denoised image $g\approx
u\odot\Phi_{\theta^{F}}+(1-u)\odot\Phi_{\theta^{B}}$.

Figure 3: The first image shows the given grayscale input image $f$, and user
defined boxes representing rough foreground and background regions. The third
image highlights pixels where the foreground expert denoiser performs better
than the background one, while the last image is the segmentation result
obtained by minimising the proposed energy (5).

######  Example 3.

Here, we give examples for neural networks that act as denoisers and relate to
existing approaches.

  * •

Constant Background: In case where the background is assumed constant, one
could simply assume that $\Phi_{\theta^{B}}(f)=\theta^{B}\mathbbm{1}$, which
corresponds to estimate a scalar value $\theta^{B}$ being the mean value of
the given image inside the corresponding region as in the original Chan and
Vese model.

  * •

Linear filter: In this case, the network is linear with respect to the network
parameters $\theta^{B}$, more precisely,
$\Phi_{\theta^{B}}(f)=\omega_{\theta^{B}}\ast f$, leading to a bi-convex
energy functional (5). In our toy example in Figure 1, we have applied such a
linear network consisting of one single filter of kernel size $15\times 15$.

  * •

Filtering of data fidelity term: When one of the regions is assumed to be
constant and high noise levels are present, mean filtering improves the
results. The data fidelity terms of energy (5) can then be replaced by
$\int_{\Omega}\left[K_{\sigma}\ast\left(f-\phi_{\theta^{F}}(f)\right)\right]^{2}u$
and
$\int_{\Omega}\left[K_{\sigma}\ast\left(f-\phi_{\theta^{B}}(f)\right)\right]^{2}(1-u)$,
respectively, where $K_{\sigma}$ is a mean filter with kernel size $\sigma$. A
similar approach has been done in [27], where a more robust version of the
Chan-Vese model [11] has been proposed by introducing Gaussian convolution in
the data fitting terms, in order to make the method robust to non homogeneous
regions.

  * •

Generic CNN: Any typical state of the art denoising neural network (Deep image
prior [39], Noise2Void [24]) can be used in our framework. Note, that in this
case the bi-convexity of energy (5) is not ensured anymore.

In the next paragraph, we discuss in more detail the joint alternating
optimisation procedure we propose to minimise energy (5).

###  3.2 Joint optimisation

We propose to iteratively optimise problem (5) with an alternating procedure
[15]. In case the denoising step does not exactly minimise energy (5), we
actually alternate between minimising two slightly different functionals. For
the sake of readability, this is not indicated in the notation. We start with
the initialiseation of the segmentation mask $u.$ This is either achieved by
thresholding for grayscale images, or as shown in Figure 3, manually choosing
boxes representing the different regions to segment in the image. Then, based
on the initial guess, the denoising expert(s) $\Phi_{\theta^{F}},$ and
$\Phi_{\theta^{B}}$ are trained on the given initial masks. To this end, we
use the ADAM optimiser [22] until convergence. As a next step, for fixed
network parameters $\boldsymbol{\theta}$, we update the segmentation mask $u$.
For fixed $\boldsymbol{\theta},$ the energy functional (5) is convex, and all
the necessary assumptions for the application of the primal dual algorithm
[10] are fulfilled. A more detailed description on the considered discrete
schemes is provided in Section 4 (see Algorithm 2). These alternate steps are
repeated as long as the decrease of energy (5) is greater than $p=15$ percent,
which we empirically found to give a good compromise between computation speed
and quality of the results.

The overall joint optimisation scheme is presented in Algorithm 1. A sketch of
the alternating procedure is provided in Figure 4.

Algorithm 1 Alternating optimisation scheme.

Initialise $u^{0}\leftarrow\boldsymbol{1}_{\\{f>\epsilon\\}}$ and
$\boldsymbol{\theta}^{0}=\boldsymbol{\theta}_{0}$

while
$\mathcal{E}^{k}_{f,\lambda}(u^{k},\boldsymbol{\theta}^{k})/\mathcal{E}^{k-1}_{f,\lambda}(u^{k-1},\boldsymbol{\theta}^{k-1})\geq
p\cdot\mathcal{E}^{k-1}_{f,\lambda}(u^{k-1},\boldsymbol{\theta}^{k-1})/\mathcal{E}^{k-2}_{f,\lambda}(u^{k-2},\boldsymbol{\theta}^{k-2})$
do

$\boldsymbol{\boldsymbol{\theta}}^{k+1}\leftarrow\operatorname*{argmin}_{\boldsymbol{\theta}}\mathcal{E}^{k}_{f,\lambda}(u^{k+1},\boldsymbol{\theta})$
{with a few ADAM iterations for $\theta^{F}$ and Chan and Vese update for the
background if $\Phi_{\theta^{B}}(f)=\theta^{B}\mathbbm{1}$}

$u^{k+1}\leftarrow\operatorname*{argmin}_{u}\mathcal{E}^{k}_{f,\lambda}(u,\boldsymbol{\theta}^{k})$
{with Algorithm 2)}

end while

In the following paragraph, we will discuss the convergence property of
Algorithm 1.

Figure 4: Alternating optimisation scheme. As a first step, regions are
provided for the training of the two denoising experts using the strategy.
These regions can be obtained by thresholding image values or by manually
choosing boxes. The differences between the given noisy image $f$ and network
outputs $\Phi_{\theta^{F}}(f)$ and $\Phi_{\theta^{B}}(f)$, are used in the
subsequent segmentation step, minimising
$\mathcal{E}_{\lambda,f}(\cdot,\boldsymbol{\theta})$ with Algorithm 2.

###  3.3 Theoretical Results

In this section, we discuss some theoretical results of the proposed energy
functional and the presented alternating algorithm. Note that these results
hold if the denoiser is trained by minimising (3).

######  Remark 4 (Monotonicity of alternating minimisation).

The proposed energy functional (5) is continuous and bounded from below.
Therefore, for each $k\geq 0$, the following relations hold

$\displaystyle\mathcal{E}_{f,\lambda}(u^{(k)},\boldsymbol{\theta}^{(k+1)})$
$\displaystyle\leq\mathcal{E}_{f,\lambda}(u^{(k-1)},\boldsymbol{\theta}^{(k)})$
$\displaystyle\mathcal{E}_{f,\lambda}(u^{(k+1)},\boldsymbol{\theta}^{(k)})$
$\displaystyle\leq\mathcal{E}_{f,\lambda}(u^{(k)},\boldsymbol{\theta}^{(k-1)}).$

Hence, the generated sequence
$\\{\mathcal{E}_{f,\lambda}(u^{(k)},\boldsymbol{\theta}^{(k)})\\}_{k\in\mathbb{N}}$
converges monotonically.

######  Theorem 3.1 (Convergence of Algorithm 1).

Assume that the level set
$S^{0}=\\{(u,\boldsymbol{\theta}):\mathcal{E}_{f,\lambda}(u,\boldsymbol{\theta})\leq\mathcal{E}_{f,\lambda}(u^{0},\boldsymbol{\theta}^{0})\\}$
of $\mathcal{E}_{f,\lambda}$ defined in (5) is compact and that
$\mathcal{E}_{f,\lambda}$ is continuous on $S^{0}$. Then, the sequence
$\\{(u^{k},\boldsymbol{\theta}^{k}\\}$ generated by Algorithm 1 is defined and
bounded. Moreover, every cluster point of
$\\{(u^{k},\boldsymbol{\theta}^{k})\\}$ is a stationary point of
$\mathcal{E}_{f,\lambda}$.

###### Proof.

This is a direct application of Theorem 4.1 in [38], using that (i) we only
alternate between two variables $u$ and $\boldsymbol{\theta}$, (ii) the
coupling between $u$ and $\boldsymbol{\theta}$ in $\mathcal{E}_{f,\lambda}$ is
smooth. ∎

######  Remark 5.

The energy (5), which is convex for fixed network parameters
$\boldsymbol{\theta}=(\theta^{F},\theta^{B})$ is a relaxation of the fully
non-convex problem

$\displaystyle\mathcal{E}(\Sigma,\boldsymbol{\theta})=\text{Per}(\Sigma,\Omega)+\int_{\Sigma}(f-\Phi_{\theta^{F}}(f))^{2}dx+\int_{\Omega\setminus\Sigma}(f-\Phi_{\theta^{B}}(f))^{2}dx,$
(6)

where $\Sigma\subset\mathbb{R}^{2}$, and $\Omega\setminus\Sigma$ are the two
regions of the given image $f(x),$ and Per$(\Sigma,\Omega)$ is the perimeter
of the interface seperating these two regions.

######  Theorem 3.2 (Thresholding).

For any fixed $\boldsymbol{\theta}$, a global minimiser for the non-convex
problem
$\min_{\Sigma,\boldsymbol{\theta}}\mathcal{E}(\cdot,\boldsymbol{\theta})$ in
(6) can be found by carrying out the minimisation of
$\min_{u}\mathcal{E}_{f,\lambda}(\cdot,\boldsymbol{\theta})$, and then setting
$\Sigma=\\{x:u(x)\geq\tau\\}$ for a.e. $\tau\in[0,1]$.

###### Proof.

The proof is similar to the one in [12](Theorem 2). The only difference is in
the data fidelity term, where instead of the fixed constants $c_{1}$, and
$c_{2}$, we look at fixed network outputs $\Phi_{\theta^{F}}(f)$, and
$\Phi_{\theta^{B}}(f)$. As the problem is one-homogeneous in $u$, thanks to
the co-area formula, we show that
$\mathcal{E}_{f,\lambda}(u,\boldsymbol{\theta})=\int_{\epsilon}\mathcal{E}_{f,\lambda}(\mathbbm{1}_{u>\tau},\boldsymbol{\theta})$,
so we can thus conclude that if $u$ is a minimiser of the energy (5) for fixed
$\boldsymbol{\theta}$, then for a.e. $\tau\in[0,1]$ the set $\Sigma(\tau)$ has
to be a minimiser of (6). ∎

##  4 Numerical Implementation

In the following, we describe the numerical implementation of the proposed
method.

###  4.1 Segmentation Step

We can rewrite our segmentation sub-problem in the form

$\displaystyle\min_{u\in\mathbb{X}}\mathcal{F}(K(u))+\mathcal{G}(u),$ (7)

where $K(u)\coloneqq\nabla u$, $\mathcal{F}(v)\coloneqq\|v\|_{1,2}$ and
$\mathcal{G}(u)\coloneqq
i_{\mathbb{A}}(u)+\int_{\Omega}(f-\Phi_{\theta_{F}}(f))^{2}u+\int_{\Omega}(f-\Phi_{\theta_{B}}(f))^{2})(1-u)$.
It holds that $K:\mathbb{X}\rightarrow\mathbb{Y}$ is a linear mapping between
Hilbert spaces $\mathbb{X},\mathbb{Y}$ and
$\mathcal{F}:\mathbb{Y}\rightarrow[0,\infty]$ and
$\mathcal{G}:\mathbb{X}\rightarrow[0,\infty]$ are convex and lower semi-
continuous functionals, i.e. all the necessary assumptions for the application
of the primal dual algorithm framework proposed in [10] are fulfilled.

####  4.1.1 Discretisation

In the following, we fix the notation which we use throughout this Section. We
work with discrete images in $\mathbb{H}\coloneqq\mathbb{R}^{N_{1}\times
N_{2}}$, denoting a finite dimensional Hilbert space equipped with an inner
product $\langle u,v\rangle=\sum_{i}u[i]v[i]$ for $u,v\in\mathbb{H}$ with
$i=(i_{1},i_{2})\in\\{1,\dots,N_{1}\\}\times\\{1,\dots,N_{2}\\}.$ The discrete
gradient
$\nabla=(\nabla_{1},\nabla_{2}):\mathbb{H}\rightarrow\mathbb{H}\times\mathbb{H}$
is defined by forward differences with Neumann boundary conditions,

$\displaystyle(\nabla_{1}u)[i]$
$\displaystyle\coloneqq\begin{cases}(u[i_{1}+1,i_{2}]-u[i_{1},i_{2}])/h&\text{if
}i_{1}<N_{1}\\\ 0&\text{if }i_{1}=N_{1}\end{cases}$
$\displaystyle(\nabla_{2}u)[i]$
$\displaystyle\coloneqq\begin{cases}(u[i_{1},i_{2}+1]-u[i_{1},i_{2}])/h&\text{if
}i_{2}<N_{2}\\\ 0&\text{if }i_{2}=N_{2}\,.\end{cases}$

Its adjoint is given by
$\nabla^{*}(v_{1},v_{2})=\nabla^{*}_{1}v_{1}+\nabla_{2}^{*}v_{2}=:-\operatorname{div}(v_{1},v_{2})$
where $\operatorname{div}\colon\mathbb{H}\times\mathbb{H}\to\mathbb{H}$ is the
discrete divergence operator and for
$(v_{1},v_{2})\in\mathbb{H}\times\mathbb{H}$ we have

$\displaystyle(\nabla^{*}_{1}v_{1})[i]$
$\displaystyle=\begin{cases}-(v_{1}[i_{1},i_{2}]-v_{1}[i_{1}-1,i_{2}])/h&\text{if
}1<i_{1}<N_{1}\\\ -v_{1}[1,i_{2}]&\text{if }i_{1}=1\\\
\phantom{-}v_{1}[N_{1}-1,i_{2}]&\text{if }i_{1}=N_{1}\end{cases}$
$\displaystyle(\nabla^{*}_{2}v_{2})[i]$
$\displaystyle=\begin{cases}-(v_{2}[i_{1},i_{2}]-v_{2}[i_{1},i_{2}-1])/h&\text{if
}1<i_{2}<N_{2}\\\ -v_{2}[i_{1},1]&\text{if }i_{2}=1\\\
\phantom{-}v_{2}[i_{1},N_{2}-1]&\text{if }i_{2}=N_{2}\,.\end{cases}$

The discrete, isotropic TV semi-norm of an image $u\in\mathbb{H}$ is defined
as

$\lVert\nabla
u\rVert_{2,1}\coloneqq\sum_{i}\sqrt{(\nabla_{1}u[i])^{2}+(\nabla_{2}u[i])^{2}}\,.$

The discrete versions of the admissible set and the corresponding indicator
function, are $\mathbb{A}=\\{u\in\mathbb{H}|0\leq u\leq 1\\}$, and
$i_{\mathbb{A}}.$ The discretisation of the data fidelity term of energy (5)
is written as $\sum_{i}{D}(u[i],\boldsymbol{\theta})$, where

$\displaystyle{D}(u,\boldsymbol{\theta})$ $\displaystyle\coloneqq
d(u,\theta^{F})+d(1-u,\theta^{B})$ (8) $\displaystyle d(u,\theta^{F})$
$\displaystyle\coloneqq u\cdot\left(\Phi_{\theta^{F}}(f)-f\right)^{2}$
$\displaystyle d(1-u,\theta^{B})$
$\displaystyle\coloneqq(1-u)\cdot\left(\Phi_{\theta^{B}}(f)-f\right)^{2}.$

Using these notations, the discrete version of energy (5) reads

$\displaystyle\mathcal{E}_{f,\lambda}^{CV}(u,\boldsymbol{\theta})=i_{\mathbb{A}}(u)+\lambda\lVert\nabla
u\rVert_{1,2}+\sum_{i}{D}(u[i],\boldsymbol{\theta})\,.$ (9)

The optimisation problem (9) is in general a non-convex and challenging
problem to be solved. We will use alternating minimisation, where we employ
for the update step of the segmentation mask $u$ the Chambolle-Pock algorithm
[10], while for updating the network parameters
$\boldsymbol{\theta}=(\theta^{F},\theta^{B})$, we apply ADAM optimisation
[22].

####  4.1.2 Segmentation algorithm

We here detail the minimisation of the functional (5) with respect to $u$ for
fixed $\boldsymbol{\theta}$, which corresponds to solve problem (7) with

$\displaystyle\mathbb{X}$ $\displaystyle=\mathbb{H}$ $\displaystyle\mathbb{Y}$
$\displaystyle=\mathbb{H}^{2}$ $\displaystyle\mathcal{F}$
$\displaystyle=\lambda\lVert v\rVert_{2,1}$ $\displaystyle K$
$\displaystyle=\nabla$ $\displaystyle\mathcal{G}$
$\displaystyle=i_{\mathbb{A}}+\sum_{i}{D}(u[i],\boldsymbol{\theta}).$

As the operator $K$ is linear, and the functionals $\mathcal{F}$ and
$\mathcal{G}$ are convex and lower semi-continuous, all requirements for the
application of the primal dual algorithm proposed in [10] are fulfilled.

To implement this algorithm, it is required to compute the Fenchel conjugate
$\mathcal{F}^{*}$ of $\mathcal{F}$, as well as the proximal mappings of
$\mathcal{F}^{*}$ and $\mathcal{G}$. We start with the derivation of the
Fenchel conjugate of $\mathcal{F}$. For $\lVert\cdot\rVert_{2,1}$ it
corresponds to the indicator function of the unit ball of the dual norm,
resulting in $i_{2,\infty}=\lVert\cdot\rVert^{\ast}_{2,\infty}$. Hence we have
$\mathcal{F}^{*}(v)=i_{2,\infty}(v/\lambda)$, the indicator function of
$\\{v|\lVert\boldsymbol{v}\rVert_{2,\infty}\leq 1\\}\subset(\mathbb{H})^{2}$.
As a next step, we compute the proximal operators of $\mathcal{F}^{\ast}$ and
$\mathcal{G}$. Recall that the proximal operator of the indicator function
$i_{C}$ of some set $C$ is given by the orthogonal projection on $C$. The
projection $P_{2,\infty}\colon\mathbb{H}\rightarrow\mathbb{H}^{2}$ onto the
unit ball in the $(2,\infty)-$norm is thus obtained by

$\displaystyle(P_{2,\infty}(\boldsymbol{v}))[i,k]=\frac{v[i,k]}{\max\\{1,(v_{1}[i,k]^{2}+v_{2}[i,k]^{2})^{1/2}\\}}.$

Thus, the proximal operator of $\mathcal{F}^{\ast}$ results in

$\displaystyle\operatorname{prox}_{\mathcal{F}^{\ast}}(\boldsymbol{v})=P_{2,\infty,\lambda}(\boldsymbol{v})\coloneqq
P_{2,\infty}(\boldsymbol{v}/\lambda).$

Further, by introducing
$\tilde{f}=\left(f-\Phi_{\theta^{F}}(f)\right)^{2}-\left(f-\Phi_{\theta^{B}}(f)\right)^{2}$,
one can finally show that

$\operatorname{prox}_{\tau\mathcal{G}}(u_{0}[i])=P_{\mathbb{A}}\left(u_{0}[i]-\tau\tilde{f}[i]\right).$

The overall primal dual Algorithm 2 is summarised below.

Algorithm 2 Segmentation algorithm based on the minimisation of the energy
functional (5) with respect to $u$ for a fixed $\theta$.

Input: noisy input image $f\in\mathbb{H}$

initialisation: $v^{0}\in\mathbb{H}$, $u^{0},\bar{u}^{0}\in\mathbb{H}$

while $\lVert u^{n+1}-u^{n}\rVert>\epsilon$ do

$v^{n+1}\leftarrow
P_{2,\infty,\lambda}(v^{n}+\sigma\boldsymbol{\nabla}\bar{u}^{n})$

$u^{n+1}\leftarrow
P_{\mathbb{A}}(u^{n}-\tau\boldsymbol{\nabla}^{\intercal}v^{n+1}-\tau\tilde{f})$

$\bar{u}^{n+1}\leftarrow u^{n+1}+\eta(u^{n+1}-u^{n})$.

end while

return $u^{n+1}$

###  4.2 Acceleration with a mask prior

In our experiments, we observed, that one of the denoisers (the one that is
trained on the more complex region), tends to improve on the other region.
Once this progress has started, it is quite difficult to stop the segmentation
mask from expanding and converging to an undesired minimum being the constant
segmentation result. Inspired by the work in [2, 42], we now propose to
overcome the problem of finding an appropriate stopping criteria, by adding a
fidelity term, ensuring that the updated segmentation mask $u^{k}$ does not
deviate too far from its initial guess. Assume that we have a reference mask
$u_{R}^{0}$, then, in the continuous setting, we consider the successive
problems:

$\displaystyle\mathcal{E}^{k}_{f,\lambda}(u,\theta)\coloneqq
i_{\mathbb{A}}(u)$ $\displaystyle+\lambda\lvert
u\rvert_{\text{TV}}+\int_{\Omega}\left(f-\Phi_{\theta^{F}}(f)\right)^{2}u$
(10)
$\displaystyle+\int_{\Omega}\left(f-\Phi_{\theta^{B}}(f)\right)^{2}(1-u)+\frac{\mu}{2}||u-u_{R}^{k}||^{2}\,.$

We can therefore optimise iteratively problem (10) with the alternate
procedure presented in Algorithm 3. Note, that in this case, as the global
energy is changed at each iteration, we do not have convergence guarantee
anymore for the alternating procedure.

Algorithm 3 Alternating optimisation scheme with acceleration.

Initialise $u^{0}\leftarrow f$ and $\theta^{0}=\theta_{0}$ and $u_{R}^{0}$

for $k=1,\dots,N$ do

$\boldsymbol{\boldsymbol{\theta}}^{k+1}\leftarrow\operatorname*{argmin}_{\boldsymbol{\theta}}\mathcal{E}^{k}_{f,\lambda}(u^{k+1},\boldsymbol{\theta})$
{with a few ADAM iterations for $\theta^{F}$ and Chan and Vese update for the
background if $\Phi_{\theta^{B}}(f)=\theta^{B}\mathbbm{1}$}

$u^{k+1}\leftarrow\operatorname*{argmin}_{u}\mathcal{E}^{k}_{f,\lambda}(\boldsymbol{u},\theta^{k})$
{with Algorithm 4)}

$\boldsymbol{u}_{R}^{k+1}=\boldsymbol{u}^{k+1}$ (update reference mask)

end for

To solve the segmentation problem, we reformulate the optimisation of problem
(10) for fixed $\boldsymbol{\theta}$ as
$\min_{u}\mathcal{F}(K(u))+\mathcal{G}^{k}(u)$, with

$\mathcal{G}^{k}(u)=i_{\mathbb{A}}(u)+\frac{\mu}{2}||u-u_{R}^{k}||^{2}+\int_{\Omega}\left(f(x)-\Phi_{\theta^{F}}(f)(x)\right)^{2}u(x)+\int_{\Omega}\left(f(x)-\Phi_{\theta^{B}}(f)(x)\right)^{2}(1-u(x)).$

Recalling that
$\tilde{f}=\left(f(x)-\Phi_{\theta^{F}}(f)(x)\right)^{2}-\left(f-\Phi_{\theta^{B}}(f)(x)\right)^{2}$,
we can show that:

$\operatorname{prox}_{\tau\mathcal{G}^{k}}(u^{0}[i])=P_{\mathbb{A}}\left(\frac{u^{0}[i]+\tau\mu
u_{R}^{k}[i]-\tau\tilde{f}[i]}{1+\tau\mu}\right).$

Observing that $\mathcal{G}^{k}$ is $\mu$-strongly convex in $u$, we consider
the accelerated primal dual algorithm of [10] to solve problem (10).

Algorithm 4 Segmentation algorithm based on the minimisation of the energy
functional (10) with respect to $u$ for a fixed $\theta$.

Input: noisy input image $f\in\mathbb{H}$

Parameters: $\lambda,\sigma,\tau,\theta$

Initialisation: $v^{0}\in\mathbb{H}$, $u^{0},\bar{u}^{0}\in\mathbb{H}$

while $\lVert u^{n+1}-u^{n}\rVert>\epsilon$ do

$v^{n+1}\leftarrow
P_{2,\infty,\lambda}(v^{n}+\sigma\boldsymbol{\nabla}\bar{u}^{n})$

$u^{n+1}\leftarrow
P_{\mathbb{A}}\left((u^{n}-\tau\boldsymbol{\nabla}^{\intercal}v^{n+1}+\tau\mu
u_{R}^{k}-\tau\tilde{f})/(1+\tau\mu)\right)$

$\eta=\frac{1}{1+2\mu\tau},\tau=\tau\eta,\sigma=\frac{\sigma}{\eta}$

$\bar{u}^{n+1}\leftarrow u^{n+1}+\eta(u^{n+1}-u^{n})$.

end while

return $u^{n+1}$

As we have discussed the numerical implementation of the segmentation step, we
now present the discrete setting and implementation of the denoising step.

###  4.3 Denoising step using Noise2Fast strategy

We here detail the denoising of a discretized 2D image
$f\in\mathbb{R}^{m\times n}$ composed of a clean signal
$g\in\mathbb{R}^{n\times n}$ and noise $n\in\mathbb{R}^{m\times n}$, i.e.

$\displaystyle f=g+n.$

For completeness, we introduce $u_{B}^{k},$ which for $k=0$ corresponds to the
Initialisation of the background region. These masks can either be obtained by
thresholding the image, or can be given in form of user-provided boxes. For
the next update steps, i.e. $k=1,\dots,N$ it holds that $u_{B}^{k}=1-u^{k}.$
Using these notations, for fixed $\boldsymbol{u}^{k}=(u^{k},u_{B}^{k})$, in
the $k$-th denoising step of our alternating procedure, the energy functional
(5), reduces to

$\displaystyle\min_{\boldsymbol{\theta}}\sum_{i}{D}(\boldsymbol{u}^{k}[i],\boldsymbol{\theta})=\min_{\boldsymbol{\theta}}\sum_{i}\left(\Phi_{\theta^{F}}(f)[i]-f[i]\right)^{2}\cdot
u^{k}[i]+\cdot\left(\Phi_{\theta^{B}}(f)[i]-f[i]\right)^{2}\cdot
u_{B}^{k}[i],$ (11)

where $\Phi_{\theta^{F}}$ and $\Phi_{\theta^{B}}$ are (deep) experts
respectively dedicated to the denoising of the foreground and background.

We build our denoisers on top of the the Noise2Fast method introduced by
Lequyer et al. in [26]. In this paper, the authors propose a fast single image
blind denoiser, using a special downsampling strategy. More precisely, their
method consists in splitting a given image into smaller parts by using a
checkerboard downsampling strategy. From a single image, four images are thus
generated, by removing one half of all pixels, and shifting the remaining
pixels to fill in the gaps left behind. Then, a network is trained to learn
the mappings between the resulting downsampled image pairs. Due to the
internal redundancy in form of recurrent patches present in images, and the
high degree of self-similarity, the neural network will also be able to
denoise the whole image instead of the downsampled ones [4, 46, 18]. For a
more detailed description of the Noise2Fast training strategy, such as the
network architecture, we refer the reader to [26].

In our approach, we use a different loss function as the one described in the
work of Lequyer et al [26]. Instead of considering the whole image domain for
training, we restrict the optimisation process for the foreground
$\Phi_{\theta^{F}}$ (resp. background $\Phi_{\theta^{B}}$) expert to the
current segmentation masks $u^{k}$ (resp. $1-u^{k}$) obtained by Algorithm 2.

In a first step, as in [26] the downsampled training images are generated in
the following way

$\displaystyle
f_{\text{even}}(i,j)=f\left(i,2j+(i\,\text{mod}\,2)\right)\in\mathbb{R}^{m\times\frac{n}{2}}$
$\displaystyle
f_{\text{odd}}(i,j)=f\left(i,2j+(i\,\text{mod}\,2)+1)\right)\in\mathbb{R}^{m\times\frac{n}{2}}$
$\displaystyle
f^{\prime}_{\text{even}}(i,j)=f\left(2i+(i\,\text{mod}\,2),j\right)\in\mathbb{R}^{\frac{m}{2}\times
n}$ $\displaystyle
f^{\prime}_{\text{odd}}(i,j)=f\left(2i+(i\,\text{mod}\,2)+1),j\right)\in\mathbb{R}^{\frac{m}{2}\times
n},$

and we repeat this downsampling procedure for the segmentation masks $u^{k}$
and $u_{B}^{k}$, for $k=0,\dots,N$ as well. We denote as

$\displaystyle\mathcal{J}^{k}=\\{(f_{\text{even}},f_{\text{odd}},u_{\text{odd}}^{k},u_{B\text{,odd}}^{k}),(f_{\text{odd}},f_{\text{even}},u_{\text{even}}^{k},u_{B\text{,even}}^{k}),$
$\displaystyle(f^{\prime}_{\text{even}},f^{\prime}_{\text{odd}},u_{\text{odd}}^{k^{\prime}},u_{B\text{,odd}}^{k^{\prime}}),(f^{\prime}_{\text{odd}},f^{\prime}_{\text{even}},u_{\text{even}}^{k^{\prime}},u_{B\text{,even}}^{k^{\prime}})\\}$

the set of training data for $k=0,\dots,N$, with $N$ being the number of
iterations of the alternating minimisation.

We then train the two denoising networks, $\Phi_{\theta^{F}}$ and
$\Phi_{\theta^{B}}$, restricted to the given regions, $u^{k}$, and
$u^{k}_{B}$, i.e. for
$(\tilde{f},\tilde{g},\tilde{u},\tilde{u}_{B})\in\mathcal{J}^{k}$ we minimise

$\displaystyle\mathcal{L}_{\boldsymbol{u}}^{k}(\boldsymbol{\theta})=\sum_{i}\left(\Phi_{\theta^{F}}(\tilde{f})[i]-\tilde{g}[i])\right)^{2}\cdot\tilde{u}[i]+\left(\Phi_{\theta^{B}}(\tilde{f})[i]-\tilde{g}([i])\right)^{2}\cdot\tilde{u}_{B}[i].$
(12)

Thus the self-supervised denoisers learn to reconstruct even (resp. odd) lines
and columns $\tilde{f}$ of the image thanks to the odd (resp. even) ones
$\tilde{g}$. As mentioned above, caused by the self-similarity redundancy, by
minimising (12), $\mathcal{L}_{u}^{k}$, we also solve problem (11).

In the next paragraph, we demonstrate the possible applications of three
different variants of the proposed joint denoising and segmentation method.

##  5 Experiments and Results

The code which was used to obtain the results presented in this work is
provided on GitHub (https://github.com/Nadja1611/Single-Image-based-
unsupervised-joint-segmentation-and-denoising.git). As a first application, we
test our method on the microscopy cell nuclei dataset from the DSB2018
dataset111https://www.kaggle.com/c/data-science-bowl-2018 stemming from the
Kaggle 2018 Data Science Bowl challenge. The data consists of a diverse
collection of cell nuclei imaged by various fluorescence microscopes. The
patches are of size $128\times 128$, and come with manually generated
segmentation ground truths. More precisely, we use the noise free data and
manually add gaussian noise with three different noise levels, namely 10, 30,
and 50. In our experiments, we considered the same subset of images than the
one used in [5], where the authors demonstrated that the segmentation of noisy
data can be improved by addressing denoising and segmentation in a cooperative
(but not fully joint) manner.

In the following experiments, for the evaluation of the segmentation
performance we use the Dice metric, and for capturing the denoising
performance in the experiments, we choose peak signal to noise ratio (PSNR)
and structural similarity metric (SSIM).

We stop our alternating Algorithm 1, as soon as the decrease of energy (5) is
less than 15 percent of the previous decrease rate. We tried out a few
different strategies, and this one turned out to be the most promising one. We
indeed observed that the a criteria based on the change in the energy decay is
robust to different scales for the regularisation parameter $\lambda$, and it
also adapts to different type of images.

We compare the segmentation performance of our joint denoising and
segmentation approach with the convex Chan-Vese model from [12] applied either
on the noisy images directly, or on the previously denoised data within a
sequential approach. For both the proposed joint approach and the sequential
one, we use the same denoised image as starting point for fair comparisons.
Further, we test our method against the partially joint denoising and
segmentation framework in [5].

### Segmentation with the constant background assumption

We start with the evaluation of our method on a subset of the DSB2018 cell
nuclei data which were manually corrupted (noise levels 10, 30 and 50). To
this end, we train a foreground denoiser, $\Phi_{\theta^{F}}$, and we assume
the background to be constant, i.e. $\Phi_{\theta^{B}}=\theta^{B}\mathbbm{1}.$
For this particular type of images, this assumption is useful, while for
images with more structural patterns, this may not be a reasonable choice, and
two denoising experts might be necessary.

To apply our joint approach, we first denoise the given image using the
Noise2Fast strategy in the way as described in Section 4, and use the
thresholded denoised image (with the threshold $\epsilon$ set to $0.5$) as
initialisation. For noise level 10, we applied the segmentation Algorithm 1
with the constant background assumption, while a Noise2Fast expert was
considered for higher noise levels. We recall that the overall process for
solving the joint segmentation and denoising is presented in Algorithm 1.
Depending on the type of image, for the alternate process between two and six
iterations are required to meet the convergence criteria.

For each method, we conducted the experiments with ten different values of the
regularisation parameter $\lambda$ evenly distributed in the interval $[0,1]$,
and then selected for each image the result with the highest Dice value.

As a further comparison, we applied the convex Chan-Vese model from [12]
directly on the noisy images. The obtained results are depicted in Figures 5
to 7, while the segmentation evaluation metrics are summarised in Table 1. We
observe that for all three noise levels, the sequential and Chan-Vese method
from [12] struggle with intensity inhomogenities of the cells. These examples
highlight the strength of the proposed unified approach, which is capable of
segmenting cells with intensities close to the mean value of the background.
Notice that the proposed approach does not perform well on the last example
due to the presence of intensity inhomogeneities, ascribed to a spatially
varying field, the bias field, in the upper left corner of the image. Please
note that in this case, evaluating the denoising performance might not be
appropriate, as we are assuming a constant background and not applying
denoising to the background.

In Table 1, the results obtained by the supervised joint denoising and
segmentation method DenoiSeg [5], are summarised. Here, we ran the provided
code using 10 training images. More precisely, we used the DSB2018 dataset
with noise level zero and added Gaussian noise in the same way as before to
all of the 4320 images, among which 3750 were used for training, 670 for
validation and the same 50 as we used for our experiments for testing. It has
to be mentioned, that for validation all 570 annotated validation images are
used in [5], resulting in a total number of 580 annotated images in total
during the training process. As displayed in Table 1, this method performs the
best. To have a fairer comparison in term of training data, we decided to
adapt their method by using 10 images in total (DenoiSeg(10 images in Table
1), 7 for training and 3 for validation. In this setting, all available data
are still used for the training of the denoiser, whereas for the segmentation
network, the ground truth masks for all but the ten training images are zeroed
out. With this smaller level of supervision, our approach outperforms the
method of [5].

Figure 5: Visual comparison of the segmentation results of data with noise
level 10. From left to right, this figure shows: the noisy input, the results
obtained with the proposed joint approach, the sequential approach, the chan-
Vese baseline and the ground truth segmentation masks. For all compared
methods, the $\lambda$ maximising the Dice score has been selected. Figure 6:
Visual comparison of the segmentation results of data with noise level 30.

For higher noise levels, it is required to filter the background fidelity
term. This prevents from considering higher values of the regularisation
parameter $\lambda$, that may lead to an over-segmentation of the background
and an overall decrease of the segmentation performance. For noise level 30
and 50, as mentioned in Section 3.1, we therefore minimise

$\begin{split}\mathcal{E}_{f,\lambda}(u,\boldsymbol{\theta})=i_{\mathbb{A}}(u)+\lambda\lvert
u\rvert_{\text{TV}}&+\int_{\Omega}\left[K_{\sigma}\ast\left(f-\phi_{\theta^{F}}(f)\right)\right]^{2}u(x)\\\
&+\int_{\Omega}\left[K_{\sigma}\ast\left(f-\phi_{\theta^{B}}(f)\right)\right]^{2}(1-u(x))dx\,\end{split}$

with $K_{\sigma}$ being a mean filter with $\sigma=3$.

The next paragraph shows experimental results which were obtained applying our
idea of training denoising experts for both regions.

Figure 7: Visual comparison of the segmentation results of data with noise level 50. noise level | $n10$ | $n30$ | $n50$
---|---|---|---
baseline | 0.820 | 0.773 | 0.582
sequential | 0.799 | 0.777 | 0.735
proposed | 0.851 | 0.825 | 0.786
DenoiSeg [5] | 0.864 | 0.848 | 0.818
DenoiSeg (10 images) | 0.843 | 0.820 | 0.750
Table 1: Dice values obtained on 50 images of the DSB2018 dataset for the
compared methods, and three different noise levels. Here, baseline is the
convex Chan-Vese [12] method directly applied to the noisy data, while for the
sequential method, we first denoise the image using Noise2Fast [26]. Our
unsupervised method almost reaches the performance of the fully supervised
approach [5].

### Segmentation using two denoisers

In the toy example in Figure 1 from Section 1.5, we trained two denoising
experts (in this case we used a linear network consisting of one filter of
size $15\times 15$) initialised by the yellow and purple boxes of size
$30\times 30$. We iterated between the denoising and segmentation three times,
until the energy decrease was less then 10 percent. For segmentation, we set
the regularisation parameter $\lambda$ to 0.02. After the first segmentation
step, the loss functions of the denoisers were restricted to $u$ and $1-u$
respectively.

Figure 8 is a typical example showing the strength of the proposed algorithm
compared to intensity-based approaches. In this experiment, we preprocessed
the given image of size 256$\times$256 in a way, that both regions have the
same mean value, and added gaussian noise as described before, with a noise
level of 10. As a consequence, the classical Chan-Vese algorithm totally fails
on this example. This model can nevertheless perform well with an adapted
hand-crafted transformation of the image to segment. As illustrated in the two
last image of Figure 8, when fed with the map of normalized image gradient
instead of the original image intensities, the Chan-Vese model is able to
segment the two part of the image.

On the other hand, our approach is able to automatically learn a relevant
transformation of image data and provides excellent segmentation without any
previous trick. The reason for that is again, that the weights learnt by the
two denoising experts strongly depend on the true underlying signal, which, in
contrast to the mean intensity, is different in the two regions. Here, both
denoising experts were initialised by boxes of size 50$\times$50 centered in
the regions. We used a regularisation parameter $\lambda$ of 0.06, and set the
learning rate to 0.001. Using the same stopping criterion as in the cell
example, these results were obtained after $3$ iterations of the alternating
procedure involving denoising and segmentation steps.

Figure 8: Segmentation of a noisy Brodatz image consisting of two different
textures. The first three images show the noisy input $f$, the minimiser of
energy (5), and the result obtained by directly applying the active contour
algorithm [11]. The fourth image shows the normalized gradient of $f$, and the
last one is the result obtained when applying the classical Chan-Vese
algorithm on the normalized gradient map.

In Figure 9, we display the clean image considered in the experiment of Figure
8, as well as different denoised images with their corresponding quantitative
metrics. More precisely, the second image in the figure is obtained by
applying the Noise2Fast strategy to the whole image, while the third image is
the result of the proposed joint optimisation procedure, where the image is
composed using the segmentation mask $u$ and the denoised images from the two
denoising experts. Especially in the left region, we can observe a better
denoising performance of the proposed method, which is more evident by
examining the PSNR (20.36 vs 19.815) and SSIM (0.753 vs 0.696) values.

Figure 9: Comparison of denoising performance with different Noise2Fast
strategies. On the middle image, Noise2Fast is applied to the whole image. On
the right image, we present the final denoised image obtained from the two
separate denoisers learned with the proposed framework.

### Segmentation with a reference mask using Algorithm 3

In Figure 10, we show another example of image segmentation for three
different noise levels using Algorithm 4. The main difficulty of this image
lies in the intensities which are shared by the object to be segmented and the
background. Therefore, we chose a representative box for initialising the
squirrel, which includes both, dark and bright areas, in order to enable the
foreground denoising expert to better generalize on the foreground region
consisting of dark and bright areas. Naturally, as the squirrel and background
do not differ much in terms of their structural properties, the foreground
denoiser, $\Phi_{\theta^{F}}$ also performs well on the background, causing
the segmentation mask $u$ to grow. In order to control this behaviour, we
applied our second strategy that includes a recursive reference mask as
described in Algorithm 3, thus preventing the segmentation mask obtained at
iteration $k+1$ from deviating too much from the previous one at iteration
$k$. More precisely, the parameters that we used for noise level 10 were
$\mu=0.0001,\lambda=0.005$, for noise level 30 we set $\mu=0.005$,
$\lambda=0.005$, while for a noise level 50 $\mu=0.00015$, $\lambda=0.005$.

Figure 10: Segmentation results obtained on the noisy images showing a
squirrel corrupted with three different noise levels. The first column shows
clean input image, and initialisation for foreground and background regions,
while in the second column the noisy versions of the given image are depicted.
The remaining ones present the segmentation results obtained using the
proposed strategy with the segmentation Algorithm 4, the segmentation masks
using the Chan-Vese algorithm provided by skimage with checkerboard
initialisation and box initialisation, respectively. The last row shows the
denoised images which are obtained by composing the obtained segmentation mask
and expert denoiser outputs.

In the following we discuss some possible extensions, and current limitations
of the proposed joint denoising and segmentation approach.

##  6 Extensions and limitations

First, our proposed unified framework can be extended to the (multichannel)
multiclass segmentation case, as we discuss in the following paragraph.

###  6.1 Vector-valued multi-class model

In order to segment a noise-corrupted vector-valued image represented as
$\boldsymbol{f}=(f_{1},\dots,f_{L})$ into $C$ different regions, we can
consider $C$ dedicated neural networks acting as denoising experts for each
region. In this case, the objective is to estimate $C$ segmentation masks
$\\{u_{i}\\}_{i=1}^{C}$ satisfying the simplex constraint, i.e.
$\sum_{k=1}^{C}u_{i}=1$, as well as the set of network parameters
$\boldsymbol{\theta}^{\text{MC}}=(\theta^{\text{MC}}_{1},\dots,\theta^{\text{MC}}_{C})$.
With these notations, the energy (5) can be extended to segment noise-
corrupted, vector-valued images $\boldsymbol{f}$ as

$\displaystyle\mathcal{E}_{f,\lambda}(\boldsymbol{u},\boldsymbol{\theta})\coloneqq
i_{\mathbb{A}}(\boldsymbol{u})$
$\displaystyle+\lambda\lvert\boldsymbol{u}\rvert_{\text{TV}}+\sum_{i=1}^{C}\sum_{j=1}^{L}\int_{\Omega}\left(f_{j}-\Phi_{\theta_{i}^{\text{MC}}}(f_{j})\right)^{2}u_{i}\,.$
(13)

As before, it may not be necessary to train $C$ different denoising networks,
as some regions may be assumed to be constant and in this case the “expert”
for region $i$ can be replaced by the mean value of the image inside region
$i$.

###  6.2 Limitations

A limitation of the current work lies in the training strategy in the case
where two denoisers are applied. In our experiments, we observed that once the
denoising experts have been trained on the initial boxes and the subsequent
segmentation step has been realised, it may occur that one of the two classes
in the image contains an important part of the other class. As a result,
during the next denoising step, one of the network is trained on parts of both
regions present in the image. With the influence of the total variation
regularisation, $u$ may converge to an undesired constant mask. With the
recursive integration of a reference mask, we already proposed in section 4.2
a strategy to overcome this drawback. One interesting alternative would be to
include an additional constraint enforcing the denoisers to perform better in
their initial regions than in the other initial ones.

Next, in some of our experiments we have observed that the Noise2Fast denoiser
is not well suited for the segmentation of certain images such as the zebra
image in Figure 3. The reason for that is that the filters of the learned
experts operate locally, and are not good in capturing global information of
the regions. As a consequence, in the case of the zebra, in regions where the
stripes are thicker, and the intensity values are closer to the one in the
background, the background expert outperforms the one of the foreground,
resulting in an undesired result as obtained by the piecewise constant Chan-
Vese model. To overcome this limitation, we modified the checkerboard strategy
of the Noise2Fast method, and instead of squeezing the image to half of its
width/height, we divided its size by a factor of four. In addition to include
different denoisers, such as for instance the deep image prior [39], an
interesting perspective would be to define new data fitting terms focusing on
structural similarities within the different classes.

##  7 Conclusion

In this work, we have proposed a novel energy functional for the joint
denoising and segmentation of images. Our framework combines the advantages of
well-established variational models with modern self-supervised deep learning
strategies. A major strength of the method lies in the fact that it can handle
single images without the need of ground truth segmentation masks or noisy-
clean training pairs. Further, the energy functional is designed in a such a
way that both tasks benefit from each other, which has also been confirmed by
experiments.

## References

  * [1] Egil Bae and Ekaterina Merkurjev.  Convex variational methods on graphs for multiclass segmentation of high-dimensional data and point clouds.  Journal of Mathematical Imaging and Vision, 58:468–493, 2017.
  * [2] Antonio Baeza, Vicent Caselles, Pau Gargallo, and Nicolas Papadakis.  A narrow band method for the convex formulation of discrete multilabel problems.  Multiscale Modeling & Simulation, 8(5):2048–2078, 2010.
  * [3] Joshua Batson and Loic Royer.  Noise2Self: Blind denoising by self-supervision.  In International Conference on Machine Learning, pages 524–533. PMLR, 2019.
  * [4] Antoni Buades, Bartomeu Coll, and Jean-Michel Morel.  Self-similarity-based image denoising.  Communications of the ACM, 54(5):109–117, 2011.
  * [5] Tim-Oliver Buchholz, Mangal Prakash, Deborah Schmidt, Alexander Krull, and Florian Jug.  DenoiSeg: joint denoising and segmentation.  In Computer Vision–ECCV 2020 Workshops: Glasgow, UK, August 23–28, 2020, Proceedings, Part I, pages 324–337. Springer, 2021.
  * [6] Xiaohao Cai.  Variational image segmentation model coupled with image restoration achievements.  Pattern Recognition, 48(6):2029–2042, 2015.
  * [7] Xiaohao Cai, Raymond Chan, Mila Nikolova, and Tieyong Zeng.  A three-stage approach for segmenting degraded color images: Smoothing, lifting and thresholding (slat).  Journal of Scientific Computing, 72:1313–1332, 2017.
  * [8] Xiaohao Cai, Raymond Chan, Carola-Bibiane Schonlieb, Gabriele Steidl, and Tieyong Zeng.  Linkage between piecewise constant Mumford–Shah model and Rudin–Osher–Fatemi model and its virtue in image segmentation.  SIAM Journal on Scientific Computing, 41(6):B1310–B1340, 2019.
  * [9] Xiaohao Cai, Raymond Chan, and Tieyong Zeng.  A two-stage image segmentation method using a convex variant of the Mumford–Shah model and thresholding.  SIAM Journal on Imaging Sciences, 6(1):368–390, 2013.
  * [10] Antonin Chambolle and Thomas Pock.  A first-order primal-dual algorithm for convex problems with applications to imaging.  Journal of mathematical imaging and vision, 40:120–145, 2011.
  * [11] Tony Chan and Luminita Vese.  An active contour model without edges.  In Scale-Space Theories in Computer Vision: Second International Conference, Scale-Space’99 Corfu, Greece, September 26–27, 1999 Proceedings 2, pages 141–151. Springer, 1999.
  * [12] Tony F Chan, Selim Esedoglu, and Mila Nikolova.  Algorithms for finding global minimizers of image segmentation and denoising models.  SIAM journal on applied mathematics, 66(5):1632–1648, 2006.
  * [13] Liang-Chieh Chen, George Papandreou, Iasonas Kokkinos, Kevin Murphy, and Alan L Yuille.  Deeplab: Semantic image segmentation with deep convolutional nets, atrous convolution, and fully connected crfs.  IEEE transactions on pattern analysis and machine intelligence, 40(4):834–848, 2017.
  * [14] Veronica Corona, Martin Benning, Matthias J Ehrhardt, Lynn F Gladden, Richard Mair, Andi Reci, Andrew J Sederman, Stefanie Reichelt, and Carola-Bibiane Schönlieb.  Enhancing joint reconstruction and segmentation with non-convex bregman iteration.  Inverse Problems, 35(5):055001, 2019.
  * [15] Imre Csiszár.  Information geometry and alternating minimization procedures.  Statistics and Decisions, Dedewicz, 1:205–237, 1984.
  * [16] Linwei Fan, Fan Zhang, Hui Fan, and Caiming Zhang.  Brief review of image denoising techniques.  Visual Computing for Industry, Biomedicine, and Art, 2:1–12, 2019.
  * [17] Ross Girshick.  Fast R-CNN.  In Proceedings of the IEEE international conference on computer vision, pages 1440–1448, 2015.
  * [18] Daniel Glasner, Shai Bagon, and Michal Irani.  Super-resolution from a single image.  In 2009 IEEE 12th international conference on computer vision, pages 349–356. IEEE, 2009.
  * [19] Nadja Gruber, Johannes Schwab, Sebastien Court, Elke Gizewski, and Markus Haltmeier.  A joint variational multichannel multiphase segmentation framework.  arXiv preprint arXiv:2202.04680, 2022.
  * [20] Kaiming He, Georgia Gkioxari, Piotr Dollár, and Ross Girshick.  Mask R-CNN.  In Proceedings of the IEEE international conference on computer vision, pages 2961–2969, 2017.
  * [21] Viren Jain and Sebastian Seung.  Natural image denoising with convolutional networks.  Advances in neural information processing systems, 21, 2008.
  * [22] Diederik P Kingma and Jimmy Ba.  Adam: A method for stochastic optimization.  arXiv preprint arXiv:1412.6980, 2014.
  * [23] Sreedhar Kollem, Katta Rama Linga Reddy, and Duggirala Srinivasa Rao.  A review of image denoising and segmentation methods based on medical images.  International Journal of Machine Learning and Computing, 9(3):288–295, 2019.
  * [24] Alexander Krull, Tim-Oliver Buchholz, and Florian Jug.  Noise2Void-learning denoising from single noisy images.  In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pages 2129–2137, 2019.
  * [25] Jaakko Lehtinen, Jacob Munkberg, Jon Hasselgren, Samuli Laine, Tero Karras, Miika Aittala, and Timo Aila.  Noise2Noise: Learning image restoration without clean data.  arXiv preprint arXiv:1803.04189, 2018.
  * [26] Jason Lequyer, Reuben Philip, Amit Sharma, and Laurence Pelletier.  Noise2Fast: Fast self-supervised single image blind denoising.  arXiv preprint arXiv:2108.10209, 2021.
  * [27] Chunming Li, Chiu-Yen Kao, John C Gore, and Zhaohua Ding.  Implicit active contours driven by local binary fitting energy.  In 2007 IEEE conference on computer vision and pattern Recognition, pages 1–7. IEEE, 2007.
  * [28] Xu Li, Xiaoping Yang, and Tieyong Zeng.  A three-stage variational image segmentation framework incorporating intensity inhomogeneity information.  SIAM Journal on Imaging Sciences, 13(3):1692–1715, 2020.
  * [29] Geert Litjens, Thijs Kooi, Babak Ehteshami Bejnordi, Arnaud Arindra Adiyoso Setio, Francesco Ciompi, Mohsen Ghafoorian, Jeroen Awm Van Der Laak, Bram Van Ginneken, and Clara I Sánchez.  A survey on deep learning in medical image analysis.  Medical image analysis, 42:60–88, 2017.
  * [30] Niklas Mevenkamp and Benjamin Berkels.  Variational multi-phase segmentation using high-dimensional local features.  In 2016 IEEE Winter Conference on Applications of Computer Vision (WACV), pages 1–9. IEEE, 2016.
  * [31] David Bryant Mumford and Jayant Shah.  Optimal approximations by piecewise smooth functions and associated variational problems.  Communications on pure and applied mathematics, 1989.
  * [32] Mangal Prakash, Tim-Oliver Buchholz, Manan Lalit, Pavel Tomancak, Florian Jug, and Alexander Krull.  Leveraging self-supervised denoising for image segmentation.  In 2020 IEEE 17th international symposium on biomedical imaging (ISBI), pages 428–432. IEEE, 2020.
  * [33] Ronny Ramlau and Wolfgang Ring.  A Mumford-Shah level-set approach for the inversion and segmentation of X-ray tomography data.  Journal of Computational Physics, 221(2):539–557, 2007.
  * [34] R Rashmi, Keerthana Prasad, and Chethana Babu K Udupa.  Multi-channel Chan-Vese model for unsupervised segmentation of nuclei from breast histopathological images.  Computers in Biology and Medicine, 136:104651, 2021.
  * [35] Olaf Ronneberger, Philipp Fischer, and Thomas Brox.  U-net: Convolutional networks for biomedical image segmentation.  In Medical Image Computing and Computer-Assisted Intervention–MICCAI 2015: 18th International Conference, Munich, Germany, October 5-9, 2015, Proceedings, Part III 18, pages 234–241. Springer, 2015.
  * [36] Leonid I Rudin, Stanley Osher, and Emad Fatemi.  Nonlinear total variation based noise removal algorithms.  Physica D: nonlinear phenomena, 60(1-4):259–268, 1992.
  * [37] Otmar Scherzer, Markus Grasmair, Harald Grossauer, Markus Haltmeier, and Frank Lenzen.  Variational methods in imaging, volume 167.  Springer, 2009.
  * [38] Paul Tseng.  Convergence of a block coordinate descent method for nondifferentiable minimization.  Journal of optimization theory and applications, 109(3):475, 2001.
  * [39] Dmitry Ulyanov, Andrea Vedaldi, and Victor Lempitsky.  Deep image prior.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 9446–9454, 2018.
  * [40] Yalin Wang, Ilhwan Jo, Stephen Wong, Shing-tung Yau, and Tony F Chan.  Segmentation and tracking of 3D neuron microscopy images using a pde based method and connected component labeling algorithm.  In 2006 IEEE/NLM Life Science Systems and Applications Workshop, pages 1–2. IEEE, 2006.
  * [41] Martin Weigert, Uwe Schmidt, Tobias Boothe, Andreas Müller, Alexandr Dibrov, Akanksha Jain, Benjamin Wilhelm, Deborah Schmidt, Coleman Broaddus, Siân Culley, et al.  Content-aware image restoration: pushing the limits of fluorescence microscopy.  Nature methods, 15(12):1090–1097, 2018.
  * [42] Romain Yildizoglu, Jean-François Aujol, and Nicolas Papadakis.  Active contours without level sets.  In 2012 19th IEEE International Conference on Image Processing, pages 2549–2552. IEEE, 2012.
  * [43] Tianming Zhan, Jun Zhang, Liang Xiao, Yunjie Chen, and Zhihui Wei.  An improved variational level set method for MR image segmentation and bias field correction.  Magnetic Resonance Imaging, 31(3):439–447, 2013.
  * [44] Kai Zhang, Wangmeng Zuo, Yunjin Chen, Deyu Meng, and Lei Zhang.  Beyond a gaussian denoiser: Residual learning of deep cnn for image denoising.  IEEE transactions on image processing, 26(7):3142–3155, 2017.
  * [45] Kaihua Zhang, Lei Zhang, Huihui Song, and Wengang Zhou.  Active contours with selective local or global segmentation: a new formulation and level set method.  Image and Vision computing, 28(4):668–676, 2010.
  * [46] Maria Zontak and Michal Irani.  Internal statistics of a single natural image.  In CVPR 2011, pages 977–984. IEEE, 2011.