# Probabilistic Graphical Models Meet Deep Learning for Semantic Segmentation

Martina Pastorino, Member, IEEE, Gabriele Moser, Fellow, IEEE, Sebastiano B. Serpico, Life Fellow, IEEE, and Josiane Zerubia, Life Fellow, IEEE

Semantic segmentation, also known as spatially dense image classification, plays a crucial role in image analysis, bridging the fields of image processing and machine learning. It has wide applications, ranging from land cover mapping in Earth observation to medical diagnostics using biomedical images, fault detection in industrial imagery, etc. This paper focuses on the mathematical connections between two pivotal families of methodological approaches: probabilistic graphical models and deep learning, and explores the potential of their integration for semantic segmentation tasks. After providing a comprehensive overview of state-of-the-art techniques from both families, the paper highlights recent developments that combine these approaches, either through theoretical equivalence or direct integration. Examples of results are provided for renowned benchmark datasets in computer vision and remote sensing, and the paper concludes with a discussion of promising future research directions.

## INTRODUCTION

Deep learning (DL) is currently the dominant approach to image semantic segmentation. Thanks to the non-parametric formulation of deep learning methods and the intrinsically multiscale processing stages that characterize convolutional neural networks, deep learning architectures can be effectively employed for several image processing tasks, often including multimodal data fusion. At the same time, probabilistic graphical models have sparked major interest in the past few years, because of the ever-growing need for structured predictions. Through the underlying graph topology over which they are defined, they offer effective and rigorous models for spatial and multiresolution information.

Specifically, a variety of DL architectures have been successfully proposed to address the task of semantic segmentation [1] and, among these, methods based on fully convolutional networks (FCNs) [2], a kind of convolutional neural networks (CNNs), and on vision transformers [1] play a primary role. In particular, FCNs share with conventional CNNs the use of convolutional and pooling layers, but are designed to produce classification maps for inputs of arbitrary size. For this purpose, they recover the spatial information lost along the downsampling operations involved in the convolutional-pooling structure, through the introduction of upsampling blocks built through unpooling and deconvolution layers [1].

However, despite their widely acknowledged experimental success, DL architectures generally exhibit two major drawbacks, as far as semantic segmentation problems are concerned. First, from a methodological perspective, as compared to other approaches to this task, whose formulations derive from well-defined optimality criteria (rooted, for example, in the theories of Bayesian decision and estimation [3], Vapnik-Chervonenkis statistical learning [4], or probably approximately correct learning [5]), the definition, design, and tuning of DL architectures, albeit intuitive, may sometimes appear empirical, when looked upon from an analytical viewpoint.

Secondly, from an operational perspective, training from scratch supervised DL models usually require large datasets, with exhaustively labeled training images to obtain accurate segmentations. Operatively, “exhaustively labeled” means that a collection of training images, in which class membership is a priori known for every individual pixel, is assumed available. These spatially dense training sets are costly and time-consuming to produce, and rarely feasible in real-world applications, for example in remote sensing or medical image processing [6], [7]. Individual examples of exhaustive label annotations are well-known in benchmark datasets collected for international contests. On the contrary, the training data that are usually available in real applications are spatially sparse and do not

M. Pastorino, G. Moser, and S. B. Serpico are with the Department of Electrical, Electronic, Telecommunications Engineering and Naval Architecture, University of Genoa, 16145 Genoa, Italy.

J. Zerubia is with the AYANA team, INRIA (National Institut for Research in Digital Science and Technology), Université Côte d'Azur, 06902 Sophia Antipolis, France.
represent spatial borders between different semantic entities – a condition that strongly affects the results of DL techniques and is a major challenge in their use for semantic segmentation. Indeed, recent progress in large-scale pretraining and foundation models has demonstrated that accurate segmentation can be achieved with much less task-specific supervision, often through self-supervised learning followed by few-shot fine-tuning [8]. Nevertheless, in the case of purely supervised training, large labeled datasets are typically required.

At the same time, stochastic models such as probabilistic graphical models (PGMs) are powerful tools for pattern recognition and image processing tasks, with special emphasis on structured prediction problems, such as image restoration [9], image denoising [10], [11], semantic segmentation, and object detection [6], [12]–[14]. These models, whose rigorous mathematical formalization is at the crossroad among graph theory, Markov modeling, and Bayesian inference, encompass Bayesian networks [15] as well as random fields [16]–[19] and characterize a structured output based on graph representations of the dependencies among random variables in a multidimensional space. For instance, Markov and conditional random fields (MRFs and CRFs) [20] are well-known and powerful families of models for spatial information, encoded according to an underlying graph topology, which can be planar [9], multiresolution [18], [21], or possibly both [6], [22]. These models also exhibit interesting analogies to statistical mechanical models (especially the well-known canonical ensemble) through the role of Boltzmann distributions (Gibbs measures) and of energy, potential, and partition functions, via the Hammersley Clifford theorem [23]. This theorem allows for not only a mathematically sound but also computationally tractable means for statistical image analysis [17].

In this framework, we address the potential of the combination of DL and PGM methodologies for semantic segmentation purposes, according to the recent literature and with special focus on two different rationales: (i) approaches allowing for a rigorous interpretation of DL models as PGMs or vice versa, through suitable equivalence theorems [24]–[29]; and (ii) methods joining PGMs and DL into common integrated frameworks [6], [26], [30]–[33]. The mathematical formulations involved in (i) and (ii) will be thoroughly discussed, *per se* and in relation to the aforementioned drawbacks of DL networks for semantic segmentation. Indeed, the opportunity to link certain DL and PGM models intrinsically allow relating those DL networks to Bayesian formulations on suitable graphs, thus supporting their interpretability from a probabilistic viewpoint. Moreover, leveraging the structured prediction capabilities of PGMs allows mitigating the training data requirements of DL networks, possibly making their use successful even in the case of highly spatially sparse training maps [6], [24], [25], [34], as shown in Figure 1. We recall that a previous survey paper on the integrations of CRFs into deep networks for semantic segmentation in computer vision was published in [26] in 2018. Differently from [26], here, the focus is on the mathematical properties that allow establishing formal equivalences between PGMs and deep models and integrating their analytical formulations. The impact of these recent developments at the crossroad between DL and PGM is also exemplified in application domains with scarce training data, such as remote sensing.

## KEY IDEAS ON THE COMBINATION OF DL AND PGM

On the one hand, deep neural and probabilistic graphical models stem from significantly different rationales, and have historically originated within distinct methodological areas of signal and image processing. Yet, in the structure and formulation of current neural networks and PGMs for image processing, it is also possible to notice important common key ideas and intuitions. On a superficial level, an obvious point of contact is the use of graphs in both methodological approaches, through explicit graph terminology in PGMs and through network terminology in DL. More in depth, for example, the notions of receptive field in a CNN and of neighborhood system in an MRF or CRF do share a partly similar intuition, related in both cases to the need for capturing local spatial context in image analysis. Similarly, the multiscale structure determined by the pooling and convolution operators in CNNs and FCNs exhibit similarities with the structure of hierarchical PGMs. In both cases, this structure meets the need to characterize multiscale and multiresolution information associated with the input imagery.

Among the families of methods at the crossroad between the DL and PGM approaches, a well-known example, which also played a pivotal role in the evolution of DL [37], is represented by restricted Boltzmann machines (RBMs). They are energy-based approaches, built corresponding to probabilistic graphical models on an undirected graph and containing a layer of observable variables and a single layer of latent variables [37]. RBMs can be considered as building blocks of deep probabilistic models – such as deep belief networks (DBNs) – and are defined on a bipartite graph, with no connections permitted between any variables in the observed layer or between
<figure><img src="image_1.png" /><figcaption>Fig. 1. Examples of semantic segmentation results with deep learning methods and PGMs in the case of scarce input training set. Input data come from the ISPRS 2D Semantic Labeling Challenge Vaihingen dataset, a well-known benchmark for semantic segmentation of aerial images. The examples refer to: stand-alone deep FCN [2]; FCN with a PGM-based post-processing [35]; FCN integrated with hierarchical PGMs [6], [36]; and PGM learned through a deep network [24], [25], [27].</figcaption></figure>

any units in the latent layer [37]. They are a special form of Boltzmann machines [38], in which the conditional probability distribution and the likelihood are both tractable because they suitably factorize. A deep Boltzmann machine (DBM) is an undirected graphical model with several layers of latent/hidden variables and, like RBMs, no connection links units of the same layer. The main difference between DBNs and DBMs lies in the connectivity: in DBNs, the connections between the top two layers are undirected, while the connections in all other layers are directed, with arrows pointing toward the layer closest to the data. DBMs can learn complex and abstract internal representations of the input. Unlike deep CNNs, they pursue the inference and training procedures in both directions, bottom-up and top-down, which helps unveiling the representations of the input structures [37].

Indeed, in the framework of semantic segmentation, Boltzmann machines, while pioneering in the field of probabilistic models and neural networks, have been largely overcome by other models, such as multilayer perceptrons (MLPs) and CNNs, due to several limitations, such as their slow speed. In particular, training DBMs is hindered by the intractability of exact maximum likelihood learning, which forces reliance on approximate methods such as mean-field inference to estimate data-dependent expectations and approximate the expected sufficient statistics by using Gibbs sampling [37]. These procedures are computationally demanding and scale poorly with network size, making joint optimization impractical for large, high-dimensional datasets. Moreover, unlike CNNs, they do not exploit inductive biases such as local connectivity and weight sharing.
Given the aforementioned drawbacks of these models, more recent solutions to embed structured predictions within a deep learning framework have involved the integration of PGMs into neural networks, or their direct formulation through deep models. Several recent approaches have been proposed to relate the formulations of CRFs and those of CNNs or recurrent neural networks (RNNs) for semantic segmentation, for instance through a CNN-based training [28], [36] or by embedding the CRF into an RNN [24]. For example, the approach in [24] proves how mean-field inference of the CRF in [35] can be modeled as an RNN and incorporated into the neural network itself. While this methodology aims to compute Gaussian CRF potentials, other techniques have also been proposed with fully non-parametric unary and pairwise potentials [25].

Then, concerning the combination of DL and PGM techniques, on one hand, a common trivial strategy is to use the latter as a post-processing stage after the former [34], to refine the classification results obtained by the neural networks. On the other hand, methods have also been proposed to address an analytically rigorous integration, often by leveraging the possible common properties of the network and graph topologies in specific families of DL and PGM models, respectively [6], [39], [40].

## REVIEW OF KEY IDEAS OF PGMS FOR SEMANTIC SEGMENTATION

Semantic segmentation is usually formulated within a probabilistic framework, in which the input image and the output segmentation map are modeled as realizations of two-dimensional random processes $\mathcal{X} = \{x_s\}_{s \in S}$ and $\mathcal{Y} = \{y_s\}_{s \in S}$, respectively, defined over the same finite pixel lattice $S \subset \mathbb{Z}^2$. The process $\mathcal{X}$ is typically vector-valued and maps each pixel $s$ to the multidimensional random vector $x_s \in \mathbb{R}^d$ of the related observations (features). The process $\mathcal{Y}$ is discrete and maps each pixel $s$ to a class label $y_s \in \Omega$ in a finite set $\Omega = \{\omega_1, \omega_2, \dots, \omega_K\}$ of semantic classes endowed with training samples [4]. The general objective is to compute an optimal estimate $\hat{y}_s$ of the label $y_s$ such that the probability of error $\hat{y}_s \neq y_s$ is minimized.

In this framework, the PGM approach aims at capturing the dependence structure associated with the random process $\mathcal{Y}$ or with its distribution conditioned on $\mathcal{X}$ [9], [17], [18], [23]. A graph $G$ (which can be planar, multiscale, bipartite, irregular, undirected, directed and acyclic, etc.) is defined, in which the edges express the desired dependencies, while probability distributions are associated with the nodes, usually involving Markovianity properties and naturally leading to Bayesian inference schemes. In this sense, PGMs generalize one-dimensional Markov chains to more general graph topologies.

From a modeling perspective, on one hand, PGMs defined over planar graphs are aimed at characterizing spatial information over the pixel lattice. The graphs are typically expressed in terms of neighborhood systems, and depending on the specific definition of Markovianity, the models can capture a variety of dependence structures, including non-causal (MRF) and causal (Markov mesh) behaviors, involving only $\mathcal{Y}$ [17], [18] or explicitly conditioning on $\mathcal{X}$ (CRF) [23]. On the other hand, PGMs defined over multiscale graphs (e.g., hierarchical PGMs on quadtrees) are aimed at modeling the multiscale information associated with multiresolution imagery, and typically involve Markovianity across the scales.

From an optimization perspective, inference on planar non-causal PGMs is usually addressed by leveraging stochastic (e.g., simulated annealing) or graph-theoretic (e.g., graph-cuts or loopy belief propagation) algorithms. Vice versa, inference on causal (planar or multiscale) PGMs can typically be formulated through sequential algorithms, which are known to be remarkably time-efficient [9], [17], [18].

## DEFINING THE EQUIVALENCE BETWEEN A DL AND A PGM MODELS FOR SEMANTIC SEGMENTATION

As mentioned above, while PGM and DL methods usually arise from different methodological rationales, common ideas and formulations can be remarked across the two fields. In this section, we discuss recent advanced approaches in which these points of contact are leveraged on a methodological basis to identify models for which a DL-PGM connection can be established in a particularly strong fashion, i.e., even formulating equivalence theorems that link a DL architecture and a PGM [24], [25], [27].

### CRFNet

In this framework, first, let us focus on the recently developed CRFNet approach in [25], which defines an equivalence relation between a CRF model and an FCN. For this purpose, it is convenient to briefly recall the basic ideas and nomenclature of CRF models. Let a neighborhood system be defined on the pixel lattice $S \subset \mathbb{Z}^2$, i.e.,
<figure><img src="image_2.png" /><figcaption>Fig. 2. Architecture of an FCN with an additional convolutional layer with K filters and softmax activations (K is the number of semantic classes of interest).</figcaption></figure>

for each pixel $s \in S$, a subset of neighbors is well defined [17]. We write $r \sim s$ to indicate that pixel $r \in S$ is a neighbor of $s$. The neighborhood relation is supposed to be:

* symmetric: $r \sim s \iff s \sim r \forall r, s \in S$;

* anti-reflexive: $s \not\sim s \forall s \in S$, i.e., no pixel is neighbor to itself. This means that, in the language of topology, the neighborhoods used in CRF modeling are always meant as punctured neighborhoods (or deleted neighborhoods).

Typical examples of neighborhood systems that are popularly used for semantic segmentation include the first- and second-order neighborhoods shown in Figure 3. Apart from obvious exceptions on the image borders, $r \sim s$ indicates that $s$ is adjacent to $r$ and that $s$ is among the eight pixels surrounding $r$ in the cases of the first- and second-order neighborhoods, respectively.

The notion of neighborhood directly reflects into a graphical structure. Specifically, let the graph $G$ be such that every pixel in $S$ is a node, and there is an edge between every pair of neighboring pixels. As a well-known result due to the Hammersley-Clifford theorem, if $\mathcal{Y}$ is a CRF with respect to this graph, then the global posterior distribution can be written as a Boltzmann distribution $P(\mathcal{Y}|\mathcal{X}) \propto \exp[-\mathcal{U}(\mathcal{Y}|\mathcal{X})]$, where the energy function $\mathcal{U}(\cdot)$ is expressed as:

$$ \mathcal{U}(\mathcal{Y}|\mathcal{X}) = \sum_{s \in S} V_s(y_s|\mathcal{X}) + \sum_{\substack{s \in S \\ r \sim s}} V_{sr}(y_s, y_r|\mathcal{X}) + \sum_{\substack{s \in S, r \sim s \\ q \sim r, q \sim s}} V_{srq}(y_s, y_r, y_q|\mathcal{X}) + \dots, \quad (1) $$

where $V_s(\cdot)$ (unary potential) is associated with the statistics of the label $y_s$ of each pixel $s$, given the observations, $V_{sr}(\cdot)$ (pairwise potential) models the spatial relation in a pair of neighboring pixels $s$ and $r$, $V_{srq}(\cdot)$ does the same for a triplet of mutual neighbors $s, r$, and $q$, and in general, the sum is extended over all subsets of mutual neighbors ("cliques").

Focusing on the CRFNet approach, let us consider an FCN in which a convolutional layer with a softmax activation is included in the output (see Figure 2). This layer is made of K convolutional filters, as many as the classes to be discriminated, and is followed by the softmax operator. The softmax layer also has K outputs, and as usual in the case of neural models for classification and segmentation, the kth softmax output returns, on each pixel $s \in S$, an estimate of the posterior probability that $s$ belongs to class $\omega_k$ ($k = 1, 2, \dots, K$). We indicate this estimated posterior as $P_{sk}^{(\text{fcn})}(\mathcal{X})$, and we denote as $h_s(\omega_k)$ and $f_s(\omega_k|\mathcal{X})$ the values, on pixel $s \in S$, of the impulse response of the kth filter and of the feature map in input to that filter (see Figure 3). These explicit notations emphasize the dependence on the class, the pixel location, and the input observations. Under mild assumptions on the support of the filters, it is possible to prove the following equivalence theorem [25]:
<figure><img src="image_3.png" /><figcaption>Fig. 3. Example of integration of spatial-contextual information: K feature maps (grey), K kernel filters (green; 3 × 3 filters are used in this example), first and second order neighborhoods (light blue; central pixel in dark blue), which can be encoded by 3 × 3 kernels.</figcaption></figure>

**Theorem 1.** The pixelwise probability distribution predicted on the output of the aforementioned FCN is equal to the local posterior distribution predicted by the CRF model defined by the following unary and pairwise potentials ($r, s \in S; r \sim s$):

$$V_s(y_s|\mathcal{X}) = -h_0(y_s)f_s(y_s|\mathcal{X}), \quad V_{sr}(y_s, y_r|\mathcal{X}) = -h_{s-r}(y_s)f_r(y_r|\mathcal{X})\delta(y_s, y_r) \quad (2)$$

and by zero higher-order potentials, in the particular case of neighboring pixels sharing the same label, i.e.:

$$P_{sk}^{(\text{fcn})}(\mathcal{X}) = P^{(\text{crf})}(y_s = \omega_k | y_r = \omega_k \forall r \sim s, \mathcal{X}) \quad \forall s \in S \quad \forall k \in \{1, 2, \dots, K\}. \quad (3)$$

<figure><img src="image_4.png" /></figure>

In (2), $\delta(\cdot)$ is the Kronecker impulse. In (3), the superscript “(crf)” emphasizes that the distribution on the right-hand side is modeled by the CRF, while $P_{sk}^{(\text{fcn})}(\mathcal{X})$ on the left-hand side is the posterior probability that the FCN assigns class $\omega_k$ at pixel $s$ given the input features $\mathcal{X}$. The neighborhood constraint in the conditioning on the right-hand side refers solely to the labels $\mathcal{Y}$ and does not involve the observation field $\mathcal{X}$. The key ideas of the proof are recalled in the Appendix.

The interpretation of the theorem is twofold. First, it allows viewing the network output as equivalent to a realization of a suitable PGM, thus providing a rigorous link between the two formulations in a probabilistic sense. In this perspective, it is worth recalling that the condition that neighboring pixels share the same label is usually satisfied within homogeneous image regions, and violated only on their borders. Then, the equivalence result for the CRFNet approach can be actively used to learn the potentials in (2) through the training of an FCN directly, thus providing a powerful computational tool to learn a CRF model directly from input image data in a fully non-parametric fashion. It is worth emphasizing that the connection between PGMs and DL established in Theorem 1 is not restricted to a particular architecture. This link arises from the fact that the convolutional operator in the final layer of the FCN has the same functional form as the pairwise energy term of a CRF. Consequently, the feature mapping function $f_s$ absorbs the non-trivial aspects of the FCN construction, while the convolution ensures that the resulting network can always be reinterpreted as the neural (non-parametric) realization of a CRF. This observation underscores the flexibility of the framework: once an FCN with a convolutional output layer is chosen, it is possible to identify a mathematical equivalence to a CRF, without requiring an architecture specifically tailored to a given graphical model. Inference in CRFNet is implicit in the forward pass of the FCN: once the network is trained, prediction reduces to a single feed-forward evaluation with no need for further (e.g., iterative) inference procedures. The equivalence theorem guarantees that this output corresponds to the local posterior distribution of a CRF. Details can be found in [25].
## CRF-RNN

In the CRF-RNN approach in [24], the dense (i.e., fully connected) CRF model in [35] is utilized and analytically linked to an RNN. The CRF model in [35] has up to pairwise nonzero potentials. Specifically, its pairwise potentials are defined as sums of Gaussian kernels, i.e. ($s, r \in S, s \sim r$):

$$V_{sr}(y_s, y_r|\mathcal{X}) = \mu(y_s, y_r) \sum_{n=1}^{N} w_n \mathcal{G}_n(\mathbf{x}_s, \mathbf{x}_r), \quad (4)$$

where the function $\mu(y_s, y_r)$ represents a measure of compatibility of the labels $y_s$ and $y_r$, $\mathcal{G}_n(\cdot)$ is a Gaussian kernel, and $w_n$ are linear combination weights. In [24], [35], the well-known Potts model is used to define the compatibility function, i.e., $\mu(y_s, y_r) = 1 - \delta(y_s, y_r)$. Indeed, the following property relates this CRF to a suitable neural model.

**Theorem 2.** The mean-field approximation of the posterior distribution $P(\mathcal{Y}|\mathcal{X})$ by a product $Q(\mathcal{Y}|\mathcal{X}) = \prod_{s \in S} Q_s(y_s|\mathcal{X})$ of independent marginals, can be expressed as an RNN. ■

More in detail, the marginal distributions in the mean-field approximation can be obtained iteratively [24]. The values $Q_{sk} = Q_s(\omega_k|\mathcal{X})$ of the marginals on all pixels and for all classes are initialized based on the unary potentials as $Q_{sk} \propto \exp[-V_s(\omega_k|\mathcal{X})]$ (with the obvious sum-to-one constraint $\sum_k Q_{sk} = 1$) and iteratively updated as follows ($s \in S; k = 1, 2, \dots, K$):

$$Q_{sk}^{\text{new}} \propto \exp \left[ -V_s(\omega_k|\mathcal{X}) - \sum_{\ell=1}^{K} \mu(\omega_k, \omega_\ell) \sum_{n=1}^{N} w_n \sum_{r \neq s} \mathcal{G}_n(\mathbf{x}_s, \mathbf{x}_r) Q_{r\ell} \right] \quad (5)$$

until convergence. First, the steps in (5) can be seen as convolutions with $1 \times 1$ kernels and softmax operations, thus the mean-field approximation of the CRF in (4) can be formulated as a stack of suitable CNN layers, which, when unrolled, form the recurrent structure of CRF-RNN. This formulation favors flexibility, as the compatibility transform is seen as a convolutional filter whose weights are learned, without forcefully restricting to the Potts model.

Then, multiple mean-field iterations can be implemented by repeating the above stack of layers in such a way that each iteration takes the mean-field marginals $Q_s(y_s|\mathcal{X})$ resulting from the previous iteration, as well as the unary values in their original form. This corresponds to unfolding the iterative mean-field inference, approximating the pairwise potentials of (4) through a sequence of RNN steps [24]. This equivalence enables the formulation of a “mean-field inference module,” which can be incorporated into any FCN or CNN, dedicated to computing the unary potentials. It is worth stressing that the inference algorithm in CRF-RNN is the well-known iterative mean-field approximation of denseCRF [35]. Details and proof can be found in [24].

## DPN

The approach in [27], called deep parsing network (DPN), learns a PGM model by extending a CNN – the VGG16 architecture in particular – with additional layers to model pairwise terms. In particular, unary and pairwise potentials are defined as follows ($s, r \in S, s \sim r$):

$$V_s(y_s|\mathcal{X}) = -\ln P_s^{(\text{vgg})}(y_s|\mathcal{X}), \quad V_{sr}(y_s, y_r|\mathcal{X}) = \sum_{c=1}^{C} \lambda_c \mu_{csr}(y_s, y_r) \sum_{q \sim r} d(\mathbf{x}_r, \mathbf{x}_q) P_q^{(\text{vgg})}(y_q|\mathcal{X}), \quad (6)$$

where $P_s^{(\text{vgg})}(y_s|\mathcal{X})$ denotes the output posteriors predicted on pixel $s \in S$ by VGG16. The pairwise term in (6) learns a mixture of local label contexts, with the goal of introducing local penalties that encode priors on the label assignment. $C$ denotes the number of components in the mixture, while the $\lambda_c$ coefficients ($c = 1, 2, \dots, C$) are Boolean and indicate which component is activated ($\lambda_c \in \{0, 1\}$, $\sum_{c=1}^{C} \lambda_c = 1$). The term $d(\mathbf{x}_r, \mathbf{x}_q)$ is a distance between pixel feature vectors and models a triple penalty, involving pixels $s$ and $r$, and the neighbors $q$ of $r$, implying that if $s$ and $r$ have compatible class labels, according to the compatibility function $\mu_{csr}(\cdot)$, then $s$ should be also compatible with the pixels nearby $r$. The function $\mu_{csr}(\omega_k, \omega_\ell)$ returns a labeling cost between $y_s$ and $y_r$ with
respect to their relative positions and their class labels $\omega_k, \omega_\ell$ ($k, \ell = 1, 2, \dots, K$). In particular, this compatibility function extends the one in CRF-RNN by incorporating dependence on mixture components and pixel locations.

Similar to CRF-RNN, a mean-field approximation $Q(\mathcal{Y}|\mathcal{X}) = \prod_{s \in S} Q_s(y_s|\mathcal{X})$ of the posterior $P(\mathcal{Y}|\mathcal{X})$ is introduced, which is obtained in DPN by minimizing the following functional [27]:

$$F(Q) = \sum_{s \in S} \sum_{k=1}^{K} Q_{sk} V_s(\omega_k|\mathcal{X}) + \sum_{s \in S} \sum_{k=1}^{K} \sum_{\ell=1}^{K} Q_{sk} Q_{r\ell} V_{sr}(\omega_k, \omega_\ell|\mathcal{X}) + \sum_{s \in S} \sum_{k=1}^{K} Q_{sk} \ln Q_{sk}, \quad (7)$$

where again $Q_{sk} = Q_s(\omega_k|\mathcal{X})$ ($s \in S; k = 1, 2, \dots, K$). Here, the first term is associated with the cost of each pixel label prediction, the second term characterizes the consistencies of the predictions between pixels, and the third is an entropy contribution measuring the confidence of the prediction. The following property holds [27]:

**Theorem 3.** The mean-field approximation of the CRF in (6) by the minimization of (7) can be expressed as a set of convolutional layers. ■

Specifically, minimizing (7) leads in closed form to the following iterative formulation ($s \in S; k = 1, 2, \dots, K$) [27]:

$$Q_{sk}^{\text{new}} \propto \exp \left[ -V_s(\omega_k|\mathcal{X}) - \sum_{c=1}^{C} \lambda_c \sum_{r \sim s} \sum_{\ell=1}^{L} \mu_{csr}(\omega_k, \omega_\ell) \sum_{q \sim r} d(\mathbf{x}_r, \mathbf{x}_q) Q_{r\ell} Q_{q\ell} \right], \quad (8)$$

where each $Q_{sk}$ is initialized as $P_s^{(\text{vgg})}(\omega_k|\mathcal{X})$ according to the unary prediction of VGG16. In the DPN approach, one iteration of (8) is decomposed into two steps, which are expressed in terms of properly defined convolutional layers with case-specific kernels that address the triple penalty and the compatibility function. Details can be found in [27].

Thanks to the previously stated theorem, DPN can jointly train a CNN – for instance, VGG16 – along with complex pairwise terms, using extra layers that approximate a single iteration of mean-field inference. As compared to CRF-RNN [24], DPN avoids recurrent operations for solving the PGM; instead, one iteration of the mean-field update is decomposed into dedicated convolutional layers, thus reducing computational cost.

## INTEGRATING DL AND PGM MODELS FOR SEMANTIC SEGMENTATION

In this section, we focus on recent methods in which the analogies between specific families of DL and PGM approaches are exploited to combine them into integrated frameworks for semantic segmentation [6], [28], [36]. In particular, here, we do not detail algorithms in which neural and probabilistic graphical techniques are used jointly by applying them in cascade, one after the other. While this strategy is sound, widely used, and often experimentally effective, here the focus is on mathematical models aimed at a methodological integration between the two approaches.

### FCN – Hierarchical PGM

Focusing first on the integrated DL-PGM method in [6], let us consider a hierarchical PGM over a quadtree: in this case, there is a finite set $\{S^0, S^1, \dots, S^L\}$ of distinct pixel lattices, related by a well-defined geometrical correspondence, rather than a single pixel grid as in the approaches discussed in the previous section. Specifically, each lattice $S^\ell \subset \mathbb{Z}^2$ ($\ell = 1, \dots, L$) is such that each pixel in $S^{\ell-1}$ spatially corresponds to a group of $2 \times 2$ pixels in $S^\ell$ [17] (see Figure 4). For semantic segmentation purposes, the pixel grid $S$ of the input image typically corresponds to $S^L$, whereas the other lattices $S^{L-1}, \dots, S^0$ are used to represent multiscale information: as $\ell$ decreases from $L$ to 0, the lattice $S^\ell$ is associated to coarser and coarser spatial resolution. In this case, the graph $G$, which is named quadtree, is such that all pixels of all lattices are nodes, while edges encode interactions across scales: every pixel $s \in S^\ell$ has a parent site $s^- \in S^{\ell-1}$ and four children sites $s^+ \subset S^{\ell+1}$ ($\ell = 0, 1, \dots, L$; with obvious exceptions on the border cases $\ell = 0$ and $\ell = L$; see Figure 4). Class labels are assumed to be defined not only on $S = S^L$ but also on all other nodes of $G$, i.e., on all levels of the quadtree and at all corresponding scales. A hierarchical PGM on this quadtree is characterized by Markovianity of the labels across the scales [17], [18], [41]. Additional edges, together with spatial Markovianity conditions, can also be introduced within each lattice $S^\ell$ to also represent spatial information at each individual scale [22].
The rationale of the integrated approach in [6] is that a very similar topology can be recognized in the graph representing the architecture of FCNs, or CNNs at large, thanks to the multiscale operations performed by their pooling and convolutional layers. Given the encoder-decoder structure of an FCN, the $(L-\ell)$th lattice $S^{L-\ell}$ of the quadtree is the same pixel grid shared by both the $\ell$th and the $(N-\ell)$th blocks of FCN activations in the network, where $N$ represents the last layer of the FCN ($\ell = 0, 1, \dots, L$; see Fig. 4). Accordingly, from a probabilistic viewpoint, the FCN and hierarchical PGM are merged, first, by incorporating the FCN activations in these blocks as the feature vectors $x_s (s \in S^\ell)$ of the pixels in lattice $S^\ell$ of the quadtree. Second, a “causal neighbor” relation $\lesssim$ is defined within each spatial lattice $S^\ell$ ($\ell = 0, 1, \dots, L-1$) so that $r \lesssim s$ formalizes that pixel $r$ is a neighbor of pixel $s$ and “precedes it” in some sense. Formally, the relation $\lesssim$ is defined by suitably combining a neighborhood system, as described above, and a total order relation on $S^\ell$ (details can be found in [22]). In the graph $G$, each node $s$, in addition to the aforementioned edges linking to its parent node $s^-$ in the upper level of the quadtree and its children nodes $s^+$ in the lower level, also have edges linking to all nodes $r$ such that $r \lesssim s$. Given this graph topology, the approach in [6] is based on the following Markovianity and conditional independence assumptions ($\ell = 1, 2, \dots, L$):

$$P(\mathcal{Y}^\ell | \mathcal{Y}^{\ell-1}, \mathcal{Y}^{\ell-2}, \dots, \mathcal{Y}^0) = P(\mathcal{Y}^\ell | \mathcal{Y}^{\ell-1}) \propto \prod_{s \in S^\ell} P(y_s | y_{s-}) P(y_s | y_r, r \lesssim s) \quad (9)$$

$$P(\mathcal{X}|\mathcal{Y}) = \prod_{\ell=1}^{L} \prod_{s \in S^\ell} P(x_s|y_s), \quad (10)$$

where $\mathcal{Y}^\ell = \{y_s\}_{s \in S^\ell}$ is the random process of the labels on lattice $S^\ell$, and $\mathcal{X}$ and $\mathcal{Y}$ are the random processes of the feature vectors and labels on the whole graph (hence, $\mathcal{Y} = \bigcup_{\ell=1}^L \mathcal{Y}^\ell$). Eq. (9) formulates Markovianity in the graph $G$ both across consecutive scales and within each lattice, thus jointly addressing both spatial-contextual modeling at every scale and multiscale fusion. Eq. (10) formalizes that the feature vectors $x_s$ are conditionally independent, given the labels $y_s$, a commonly accepted hypothesis with PGMs for classification. For this integrated hierarchical PGM-FCN framework, desirable mathematical properties can be proven, which make it possible to define a particularly efficient inference scheme. Specifically, the following theorem holds under mild assumptions (the proof can be found in [6], [22]).

**Theorem 4.** For the integrated hierarchical PGM-FCN framework defined by (9) and (10), the following properties hold:

1) The PGM is causal;

2) There exists a sequential formulation (with close-form analytical expressions for all processing stages) for the marginal posterior mode (MPM) criterion. ■

Here, the notion of causality in property 1 refers to the graph $G$ and its topology, as defined by the edges linking pixels (nodes) across consecutive levels of the quadtree and within each level. In general, causality in

<figure><img src="image_5.png" /><figcaption>Fig. 4. Structural analogy between an encoder-decoder based neural network (left) and a quadtree topological structure (right).</figcaption></figure>
PGMs favors the introduction of computationally efficient sequential inference algorithms. Property 2 ensures that this is precisely feasible for the MPM inference in the case of the described PGM-FCN framework. We recall that MPM is the decision rule that assigns pixel s the class label $y_s$ that maximizes its posterior $P(y_s|\mathcal{X})$ conditioned to all observations $\mathcal{X}$, and that MPM is known to be particularly appropriate in the case of a quadtree, because it allows penalizing errors as a function of the scale where they are made [18], [41]. Owing to the causal quadtree topology, MPM can be computed sequentially with closed-form updates, yielding an efficient inference algorithm that scales linearly with the number of nodes [6], [18].

## Unary-Pairwise Net

Another recent example of approach integrating PGMs and deep learning is presented in [36]. This approach, named Unary-Pairwise Net, is aimed at exploring the use of contextual information between image regions. The method includes two separate subnetworks, named unary and pairwise networks and aimed at defining the potentials of a suitable CRF.

A CNN-based CRF pairwise potential is introduced for pixels within a range box to model the semantic compatibility between two nodes and the output, for every possible value of the pair of labels of those nodes. Both the unary and pairwise potentials are formulated as multiscale CNNs for capturing semantic relations between image regions. The pairwise potentials learnt by a CNN are not restricted to a specific parametric model (e.g., the Potts model), but are non-parametric. They model the semantic compatibility between image regions with the output for every possible value of the label pair $(y_s, y_r)$, individually parameterized by CNNs. The aim is to improve a coarse-level prediction, rather than local smoothness, and thus have a different purpose compared to Potts-like pairwise potentials.

Given an input image, first a convolutional network is applied to generate a feature map at a lower spatial resolution. Then, a graph $G = (\mathcal{N}, \mathcal{E})$, where $\mathcal{N}$ and $\mathcal{E}$ denote as usual the sets of nodes and edges, is created as follows: each location in the feature map (which corresponds to a rectangular region in the input image) is a node. The edges are constructed by connecting one node to all other nodes that lie within a spatial range box. In addition to the feature vectors corresponding to the nodes of the graph, a feature vector is also associated with each edge by concatenating the feature vectors of the two nodes it connects. Recalling that $K$ is the number of classes, the unary potential at each node is the $K$-dimensional output of a shallow fully connected network applied to the node feature vectors ("Unary Net"), while the pairwise potentials are the $K \times K$ outputs of the shallow fully connected network applied to the pair of concatenated feature vectors ("Pairwise Net"). Accordingly, an energy function is defined as follows:

$$ \mathcal{U}(\mathcal{Y}|\mathcal{X}) = \sum_{s \in \mathcal{N}} V_s(y_s|\mathcal{X}, \vartheta_U) + \sum_{(s,r) \in \mathcal{E}} V_{sr}(y_s, y_r|\mathcal{X}, \vartheta_P), \quad (11) $$

where $V_s(\cdot)$ and $V_{sr}(\cdot)$ denote the outputs of the Unary and Pairwise Nets, respectively, and where $\vartheta_U$ and $\vartheta_P$ are their parameter vectors. Indeed, (11) is immediately recognized as a CRF energy with up to pairwise potentials.

Since the graph $G$ does not form a tree structure, it is generally necessary to apply an approximate inference algorithm, which is achieved through a mean field approximation via message passing [36] on the constructed graph. The training of the CRF of the Unary-Pairwise Net approach with the CNN potentials is performed in a piecewise manner, with the conditional likelihood being formulated as independent likelihoods defined on the potentials [36].

Operatively, the Unary-Pairwise Net approach produces a low-resolution prediction following CRF inference. This prediction is then up-sampled and refined using the denseCRF in [35] to obtain the final result, adding an additional mean-field inference step as post-processing. The method in [28] provides a further extension of this last method [36], aimed at hyperspectral image segmentation with 3D CNNs.

## EXAMPLES OF SEMANTIC SEGMENTATION RESULTS

To further investigate the aforementioned DL-PGM approaches, examples of experimental results are discussed in this section. In particular, the performances have been demonstrated using standard benchmark datasets, widely renowned in two major fields of image processing, i.e., remote sensing image analysis and computer vision. In particular, we shall present examples of experimental results with: (i) the ISPRS 2D Semantic Labeling Challenge Vaihingen dataset, an aerial dataset collected by the German Society for Photogrammetry, Remote Sensing, and Geoinformation (DGPF) over the city of Vaihingen, Germany; and (ii) the Pascal Visual Object Classes Challenge
<figure><img src="image_6.png" /><figcaption>Fig. 5. Example of spatially exhaustive and scarce training data in the case of the aerial ISPRS 2D Semantic Labeling Challenge Vaihingen dataset. An exhaustive label map is characterized by a semantic label on every pixel. A scarce label map only has a few labeled samples, often arranged in separate homogeneous regions, and lacks boundary information between different classes.</figcaption></figure>

2012 (VOC2012), a benchmark dataset in visual object category recognition and detection framed within the EU Network of Excellence funded under the IST Programme of the European Union. In the present tutorial article, we shall focus on high-level comparisons, besides the aforementioned remarks on the inference complexity of the considered methods. While we will not provide exhaustive implementation details here, we recall that hyper-parameter tuning guidelines as well as empirical analyses of performance sensitivity to hyper-parameter values and of computational cost are reported in the related journal publications [6], [24], [25], [27], [36].

## Remote sensing image analysis

The semantic segmentation of remote sensing images involves classifying each pixel of satellite or aerial imagery into predefined categories or classes, such as vegetation, water bodies, urban areas, roads, etc. This process allows for the extraction of detailed, spatially explicit information from remotely sensed data. In remote sensing, semantic segmentation plays a key role in a variety of applications, including land cover mapping, climate change monitoring, disaster management, and urban planning.

Scarce ground truth in remote sensing applications presents an especially significant challenge in developing accurate semantic segmentation models, particularly when it comes to training machine learning algorithms. Ground truth data—typically collected via field surveys or visual photointerpretation—is essential for labeling and validating the classes in remote sensing datasets (e.g., distinguishing between vegetation, water bodies, or urban areas). However, obtaining comprehensive and reliable ground truth can be logistically difficult, time-consuming, and expensive, especially in remote or hard-to-reach areas, or when dealing with rapidly changing landscapes such as urban environments or deforested regions. In the aforementioned ISPRS Vaihingen dataset or in similar datasets designed for benchmark competitions, the ground truth is an “ideal” spatially exhaustive one (see above Introduction), where the true label is known for all pixels in the training maps. This is generally unfeasible in real-world remote-sensing applications, where the ground truth maps are usually arranged in isolated homogeneous patches not including spatial class borders. Such a scenario exemplifies well an applicative situation in which, as discussed in the Introduction, the training set is significantly suboptimal, thus generally hindering the performances of DL-based semantic segmentation methods.

To experimentally investigate this aspect and the possible benefits of mitigating training data scarcity through the combination with structured output learning based on PGMs, the original ground truth of the Vaihingen dataset was modified to approximate the ground truths usually found in realistic remote-sensing applications, by removing entire connected components from the original exhaustive ground truth and then applying morphological erosion, until only few labeled pixels were left (see Figure 5). The results presented here refer to the experiments conducted with this scarce label information.

Examples of semantic segmentation maps obtained by some of the aforementioned approaches, when applied to the Vaihingen dataset, are presented in Figure 6. The contours of the homogeneous regions of the ground truth are superimposed both to the ground truth itself and to the segmentation maps produced by the various methods,
<figure><img src="image_7.png" /><figcaption>Fig. 6. Example of semantic segmentation results obtained using scarce training data in the case of the aerial ISPRS 2D Semantic Labeling Challenge Vaihingen dataset. The examples involve: stand-alone FCNs [2]; integration of FCN and hierarchical PGMs [6]; and methods learning a CRF through deep networks [24], [25]. Color legend for the classes: buildings (red), trees (dark green), low vegetation (green), streets (white), cars (grey).</figcaption></figure>

to make the visual comparison easier. As a baseline reference, the results of an FCN and of its combination in cascade with a CRF (first the FCN, then the CRF as post-processing) are also shown. The basic FCN in [2] and the denseCRF model in [35] are used for this purpose. In Figure 6, we focus on one of the test tiles, however, the behavior on the other test tiles of the dataset is similar. As expected, when trained with the scarce input training set, the FCN provides suboptimal results. On the contrary, the maps obtained using CRFNet and CRF-RNN – i.e., two of the techniques that characterize spatial information through deep CRFs and involve equivalence properties between neural and CRF models –, and through the integration of the FCN and of a hierarchical PGM, exhibit more homogeneous zones and sharper edges than the FCN. The classification results of CRFNet are particularly remarkable for the class “building”.

Concerning the quantitative results shown in Table I, CRFNet attains the best performances in terms of all the averaged accuracy metrics. The results are presented in terms of overall accuracy, recall, precision, and F1 score (the last three are averaged over the set of classes). Consistently with a conformal prediction approach, the results in Table I derive from 20 random calibration/test splits of non-overlapping windows pooled from the test tiles. The other considered methods also provide accurate results. CRF-RNN reaches the second highest values in terms of the F1 score and average precision, whereas FCN-Hierarchical PGM reaches the second highest values in recall and overall accuracy, together with the method employing denseCRF post-processing.

These results suggest the effectiveness of the integrated DL-PGM approaches, at least on the considered dataset. They confirm the potential of the integration between CRF and deep learning concepts, in a twofold manner. On one hand, they emphasize the opportunity to exploit deep architectures to learn the spatial modeling structure of PGM models, as pointed out by the effective performance of CRFNet and CRF-RNN. On the other hand, the results also confirm that the integration with structured output priors defined by PGMs allows mitigating the impact of spatially scarce training data on the segmentation performance of deep neural networks.

## Computer vision

The Pascal Visual Object Classes (VOC) dataset is a widely recognized benchmark in the field of computer vision, particularly for object detection, image classification, and segmentation tasks. Introduced as part of the Pascal VOC challenges, the dataset provides a diverse collection of images with rich annotations, including object bounding boxes, class labels, and segmentation masks across 20 object categories such as person, vehicle, animal,
<table><caption>TABLE I
TEST-SET RESULTS ON THE ISPRS 2D SEMANTIC LABELING CHALLENGE VAIHINGEN DATASET WITH SCARCE INPUT TRAINING SET. AVERAGED OVER 20 RANDOM CALIBRATION/TEST SPLITS OF NON-OVERLAPPING WINDOWS POOLED FROM THE 4 TEST TILES (NOTATION: MEAN ± STANDARD DEVIATION).</caption><thead><tr><th>Architecture</th><th>overall accuracy</th><th>average recall</th><th>average precision</th><th>average F1 score</th></tr></thead><tbody><tr><td><em>DL</em><br>FCN [2]</td><td>0.80 ± 0.01</td><td>0.80 ± 0.02</td><td>0.83 ± 0.01</td><td>0.81 ± 0.01</td></tr><tr><td><em>PGM as post-processing</em><br>FCN + denseCRF [35]</td><td>0.82 ± 0.01</td><td>0.78 ± 0.01</td><td>0.84 ± 0.03</td><td>0.81 ± 0.02</td></tr><tr><td><em>Integration DL + PGM</em><br>FCN-Hierarchical PGM [6]</td><td>0.82 ± 0.01</td><td>0.85 ± 0.02</td><td>0.72 ± 0.02</td><td>0.78 ± 0.02</td></tr><tr><td><em>End-to-end deep PGM</em><br>CRF-RNN [24]</td><td>0.81 ± 0.01</td><td>0.80 ± 0.02</td><td>0.84 ± 0.01</td><td>0.82 ± 0.01</td></tr><tr><td>CRFNet [25]</td><td><strong>0.85</strong> ± 0.01</td><td><strong>0.87</strong> ± 0.01</td><td><strong>0.86</strong> ± 0.02</td><td><strong>0.86</strong> ± 0.01</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>Fig. 7. Examples of results obtained on the Pascal VOC 2012 validation set. The examples involve: deep learning [34]; integration of deep learning and PGMs [27]; and methods learning a CRF through deep networks [24], [25]. Figures taken from the papers [24], [27], [34]. Color legend for the classes: person (pink) and bicycle (green).</figcaption></figure>

<table><caption>TABLE II
RESULTS ON THE PASCAL VOC 2012 VALIDATION SET. FOR [27], [34], RESULTS ARE TAKEN DIRECTLY FROM THE CORRESPONDING PUBLICATIONS, AS RETRAINING WAS NOT FEASIBLE. ALL RESULTS ARE REPORTED UNDER A CONSISTENT TRAINING STRATEGY AND WITH REGARD TO THE SAME TRAINING/TEST SPLIT.</caption><thead><tr><th></th><th>mIoU</th></tr></thead><tbody><tr><td><em>DL</em><br>FCN [2]</td><td>0.62</td></tr><tr><td>DeepLab (before CRF) [34]</td><td>0.65</td></tr><tr><td><em>PGM as post-processing</em><br>FCN + denseCRF [35]</td><td>0.64</td></tr><tr><td>DeepLab (after CRF) [34]</td><td>0.69</td></tr><tr><td><em>Integration DL + PGM</em><br>FCN-Hierarchical PGM [6]</td><td>0.65</td></tr><tr><td>Unary-Pairwise Net [36]</td><td>0.70</td></tr><tr><td><em>End-to-end deep PGM</em><br>CRF-RNN [24]</td><td>0.70</td></tr><tr><td>CRFNet [25]</td><td>0.70</td></tr><tr><td>DPN [27]</td><td>0.68</td></tr></tbody></table>
household items, etc. Its comprehensive dataset has made it a well-known resource for evaluating and comparing the performance of machine learning models, and for facilitating standardized evaluation protocols under consistent conditions. The dataset complexity, with objects in varied poses, scales, and occlusions, aims at exemplifying real-world scenarios. In this section, we present examples of experimental results obtained by some of the aforementioned DL-PGM approaches, when applied to the Pascal VOC dataset. In particular, quantitative results are presented in Table II in terms of the mean intersection over unit (mIoU), the standard average metric to evaluate performance on this data set. Here, in addition to the FCN in [2] and the processing with denseCRF, the DeepLab network, which integrates atrous convolutions and fully connected CRF, is also used as a further baseline reference [34].

Given the challenge presented by the dataset, the average performance metrics are lower than those obtained on the previous dataset; however, they emphasize even more clearly the power of approaches combining deep learning and stochastic models defined over graphs. Specifically, the best results are obtained by the Unary-Pairwise-Net [36] and by the methods defining CRF models through deep networks (namely, CRFNet [25], CRF-RNN [24], and DPN [27]). The performances of the FCNs are characterized by lower accuracies on this dataset. This is confirmed by the classification maps shown in Figure 7, which appear to be smoother and with better defined edges for the methods involving the integrated / deep CRFs [24], [25], [27], and a bit more fragmented for the purely deep learning method [34].

CONCLUSION

This paper explores the mathematical relationships between probabilistic graphical and deep learning models within the framework of image semantic segmentation. By analyzing recent advanced approaches from both methodological families and examining their integration, their complementary advantages and the potential of their synergy have been discussed. In particular, DL provides flexibility, strong representation learning, and non-parametric modeling power, while PGMs contribute structured prediction, interpretability, and the ability to incorporate spatial and multiscale dependencies. Mathematical models that integrate neural and stochastic graphical components into unified frameworks are therefore especially appealing, as they can leverage these combined strengths to achieve more robust and accurate segmentation. Moreover, for specific approaches, theoretical properties can even be proven that ensure the equivalence, under suitable assumptions, between deep neural architectures and case-specific PGM formulations. These equivalences are fruitful in both directions. On one hand, in addition to favoring structured prediction, they endow the corresponding networks with rigorous Bayesian interpretations with respect to well-defined graph topologies. On the other hand, they allow the PGM structure and parameters (e.g., CRF potentials) to be automatically learned in a non-parametric fashion through network training.

Several representative approaches discussed in this tutorial paper, such as CRF-RNN [24], CRFNet [25], DPN [27], the Unary-Pairwise Net [36], and the FCN-Hierarchical PGM framework [6], exemplify different strategies for integrating probabilistic graphical models with deep networks. CRF-RNN relies on iterative mean-field updates unrolled as recurrent layers, which enable end-to-end learning while preserving the interpretability of CRFs, though at a higher computational cost. CRFNet instead frames inference as the direct output of an FCN, which makes training straightforward and efficient. Yet, the corresponding equivalence theorem involves case-specific neighborhood assumptions. DPN leverages convolutional layers to approximate mean-field inference with higher-order potentials, thus avoiding recurrence and capturing richer dependencies, but at the expense of a more complex architecture. The Unary-Pairwise Net employs message passing over learned CNN-based potentials, combining flexibility with general graph structures while still relying on approximate inference. Finally, the FCN-Hierarchical PGM exploits quadtree structures that allow for causal sequential inference, providing efficiency and multiscale modeling capabilities.

Experimental results on well-known benchmark datasets have exemplified the benefits of these integrated DL-
PGM strategies, at least in the considered applications to remote sensing image analysis and computer vision. The
results have confirmed the capability of the presented techniques at taking benefit from spatial information and at
mitigating the training data requirements of deep networks for semantic segmentation. In particular, the considered
approaches combining DL and PGMs have demonstrated capable to produce smoother, more homogeneous clas-
sification maps with sharper edges than conventional FCNs. Moreover, the methods integrating deep learning and
stochastic models, or emulating the capabilities of stochastic models based on random fields via deep learning, are
shown to be particularly effective in limiting the impact, on segmentation performance, of scarce ground-truth data,
a common scenario in multiple image processing tasks.
Moving forward, future research is expected to focus on further refining the integration of these two methodolog-
ical families, improving computational efficiency, and expanding their applicability to emerging challenges in image
processing, computer vision, and machine learning. As the scientific community is currently moving towards the
availability of larger and larger datasets, a crucial point is the development of more scalable and efficient algorithms
that can effectively leverage these vast amounts of high-dimensional and multi-source data in real-world scenarios,
while maintaining accuracy and interpretability. From this perspective, the mathematical relations amongst deep
and probabilistic graphical approaches and the opportunity to integrate their paradigms together appear especially
promising, both to further enhance performance in current semantic segmentation applications and to provide a
more in-depth insight into these machine learning methodologies, paving the way for more robust and interpretable
semantic segmentation techniques.

APPENDIX

Sketch of proof of Theorem 1. Let us focus on pixel $s \in S$ and class $\omega_k$ ($k = 1, 2, \dots, M$). First, the probability $P_{sk}^{(\text{fcn})}(\mathcal{X})$ of $\omega_k$ predicted by the network on pixel $s$ is obtained on the output of the cascade of a softmax and a convolutional layers. Given the definition of discrete 2D convolution on the lattice $S$ and suitably relating the support of the filters with the neighborhood system (details can be found in [25]), the predicted probability can be expressed, in log-scale, as:

$$
\ln P_{sk}^{(\text{fcn})}(\mathcal{X}) = h_0(\omega_k) f_s(\omega_k | \mathcal{X}) + \sum_{r \sim s} h_{r-s}(\omega_k) f_s(\omega_k | \mathcal{X}) + [\dots], \quad (12)
$$

where [...] indicates additive terms that do not depend on the class $\omega_k$.

Second, let us consider the CRF model with up to nonzero pairwise potentials defined by (2). The local posterior energy function is given by:

$$
\mathcal{U}_s(y_s|y_r, \forall r \sim s, \mathcal{X}) = V_s(y_s|\mathcal{X}) + \sum_{r \sim s} V_{sr}(y_s, y_r|\mathcal{X}) \quad (13)
$$

and is, up to an additive constant, the negative logarithm of the local posterior distribution determined by the CRF, i.e.,

$$
\begin{align}
\ln P^{(\text{crf})}(y_s = \omega_k | y_r, \forall r \sim s, \mathcal{X}) &= -\mathcal{U}_s(\omega_k | y_r, \forall r \sim s, \mathcal{X}) + [\dots] = \nonumber \\
&= h_0(\omega_k) f_s(\omega_k | \mathcal{X}) + \sum_{r \sim s} h_{r-s}(\omega_k) f_r(\omega_k | \mathcal{X}) \delta(\omega_k, y_r) + [\dots], \tag{14}
\end{align}
$$

where we have plugged the expressions of the potentials in (2). If $y_r = \omega_k$ for all $r \sim s$, i.e., if the neighborhood of pixel $s$ shares the same class label $\omega_k$, then $\delta(\omega_k, y_r) = 1$ for all $r \sim s$. Therefore (12) and (14) coincide up to an additive term that does not depend on the class $\omega_k$:

$$
\ln P_{sk}^{(\text{fcn})}(\mathcal{X}) = \ln P^{(\text{crf})}(y_s = \omega_k | y_r, \forall r \sim s, \mathcal{X}) + [\dots], \quad (15)
$$

Due to the obvious sum-to-1 constraint with respect to $k$, (15) proves the theorem.

<figure><img src="image_9.png" /></figure>

REFERENCES

[1] C. M. Bishop and H. Bishop, *Deep Learning: Foundations and Concepts*, 2023.

[2] E. Shelhamer, J. Long, and T. Darrell, “Fully convolutional networks for semantic segmentation,” *IEEE Transactions on Pattern Analysis and Machine Intelligence*, vol. 39, no. 4, pp. 640–651, 2017.

[3] H. L. Van Trees, K. L. Bell, and Z. Tian, *Part 1: Detection, estimation, and filtering theory*. Wiley, 2013.

[4] V. N. Vapnik, *Statistical learning theory*. Wiley Interscience, 1998.

[5] N. Cristianini and J. Shawe-Taylor, *An Introduction to Support Vector Machines*. Cambridge University Press, 2000.

[6] M. Pastorino, G. Moser, S. B. Serpico, and J. Zerubia, “Semantic segmentation of remote-sensing images through fully convolutional neural networks and hierarchical probabilistic graphical models,” *IEEE Trans. Geosci. Remote Sens.*, vol. 60, no. 5407116, pp. 1–16, 2022.

[7] P.-H. Conze, G. Andrade-Miranda, V. K. Singh, V. Jaouen, and D. Visvikis, “Current and emerging trends in medical image segmentation with deep learning,” *IEEE Transactions on Radiation and Plasma Medical Sciences*, vol. 7, no. 6, pp. 545–569, 2023.

[8] M. Awais, M. Naseer, S. Khan, R. M. Anwer, H. Cholakkal, M. Shah, M.-H. Yang, and F. S. Khan, “Foundation models defining a new era in vision: A survey and outlook,” *IEEE Transactions on Pattern Analysis and Machine Intelligence*, vol. 47, no. 4, pp. 2245–2264, 2025.
[9] S. Geman and D. Geman, "Stochastic relaxation, Gibbs distributions, and the Bayesian restoration of images," *IEEE Trans. Pattern Anal. Mach. Intell.*, vol. PAMI-6, no. 6, pp. 721-741, 1984.

[10] D. Zheng, X. Zhang, K. Ma, and C. Bao, "Learn from unpaired data for image restoration: A variational Bayes approach," *IEEE Trans. Pattern Anal. Mach. Intell.*, pp. 1-15, 2022.

[11] J. Shen, P.-C. Su, S.-C. S. Cheung, and J. Zhao, "Virtual mirror rendering with stationary RGB-D cameras and stored 3-D background," *IEEE Trans. Image Process.*, vol. 22, no. 9, pp. 3433-3448, 2013.

[12] D. Liu, G. Hua, and T. Chen, "A hierarchical visual model for video object summarization," *IEEE Trans. Pattern Anal. Mach. Intell.*, vol. 32, no. 12, pp. 2178-2190, 2010.

[13] J. van de Ven, F. Ramos, and G. D. Tipaldi, "An integrated probabilistic model for scan-matching, moving object detection and motion estimation," in *Int. Conf. Robot. Autom. (ICRA)*, 2010, pp. 887-894.

[14] Y. Wang, J. Lin, Q. Cai, Y. Pan, T. Yao, H. Chao, and T. Mei, "A low rank promoting prior for unsupervised contrastive learning," *IEEE Trans. Pattern Anal. Mach. Intell.*, pp. 1-14, 2022.

[15] N. Friedman, D. Geiger, and M. Goldszmidt, "Bayesian network classifiers," *Mach. Learn.*, vol. 29, pp. 131-163, 1997.

[16] P. A. Devijver, "Hidden Markov mesh random field models in image analysis," *J. Appl. Stat.*, vol. 20, no. 5-6, pp. 187-227, 1993.

[17] Z. Kato and J. Zerubia, "Markov random fields in image segmentation," *Found. Trends Signal Process.*, vol. 5, no. 1-2, pp. 1-155, 2012.

[18] J. Laferté, P. Pérez, and F. Heitz, "Discrete Markov image modeling and inference on the quadtree," *IEEE Trans. Image Process.*, vol. 9, no. 3, pp. 390-404, 2000.

[19] A. Blake, P. Kohli, and C. Rother, *Markov Random Fields for Vision and Image Processing*. The MIT Press, 2011.

[20] J. D. Lafferty, A. McCallum, and F. C. N. Pereira, "Conditional random fields: Probabilistic models for segmenting and labeling sequence data," in *Proc. Int. Conf. Mach. Learn. (ICML)*. San Francisco, CA, USA: Morgan Kaufmann Publishers Inc., 2001, p. 282-289.

[21] I. Hedhli, G. Moser, S. B. Serpico, and J. Zerubia, "Classification of multisensor and multiresolution remote sensing images through hierarchical Markov random fields," *IEEE Geosci. Remote Sens. Lett.*, vol. 14, no. 12, pp. 2448-2452, 2017.

[22] M. Pastorino, A. Montaldo, L. Fronda, I. Hedhli, G. Moser, S. B. Serpico, and J. Zerubia, "Multisensor and multiresolution remote sensing image classification through a causal hierarchical Markov framework and decision tree ensembles," *Remote Sens.*, vol. 13, no. 5, p. 849, 2021.

[23] C. Sutton, A. McCallum, and F. Pereira, "An introduction to conditional random fields," *Found. Trends Mach. Learn.*, vol. 4, no. 4, pp. 267-373, 2011.

[24] S. Zheng, S. Jayasumana, B. Romera-Paredes, V. Vineet, Z. Su, D. Du, C. Huang, and P. H. S. Torr, "Conditional random fields as recurrent neural networks," in *IEEE Int. Conf. Comp. Vis. (ICCV)*, 2015, pp. 1529-1537.

[25] M. Pastorino, G. Moser, S. B. Serpico, and J. Zerubia, "CRFNet: A deep convolutional network to learn the potentials of a CRF for the semantic segmentation of remote sensing images," *IEEE Transactions on Geoscience and Remote Sensing*, vol. 62, pp. 1-19, 2024.

[26] A. Arnab, S. Zheng, S. Jayasumana, B. Romera-Paredes, M. Larsson, A. Kirillov, B. Savchynskyy, C. Rother, F. Kahl, and P. H. Torr, "Conditional random fields meet deep neural networks for semantic segmentation: Combining probabilistic graphical models with deep learning for structured prediction," *IEEE Signal Process. Mag.*, vol. 35, no. 1, pp. 37-52, 2018.

[27] Z. Liu, X. Li, P. Luo, C. C. Loy, and X. Tang, "Deep learning Markov random field for semantic segmentation," *IEEE Transactions on Pattern Analysis and Machine Intelligence*, vol. 40, no. 8, pp. 1814-1828, 2018.

[28] F. I. Alam, J. Zhou, A. W.-C. Liew, X. Jia, J. Chanussot, and Y. Gao, "Conditional random field and deep feature learning for hyperspectral image classification," *IEEE Trans. Geosci. Remote Sens.*, vol. 57, no. 3, pp. 1612-1628, 2019.

[29] Z. Zhang, T. Zhao, H. Gay, W. Zhang, and B. Sun, "ARPM-net: A novel CNN-based adversarial method with Markov random field enhancement for prostate and organs at risk segmentation in pelvic CT images," *Med Phys*, vol. 48, no. 1, pp. 227-237, 2021.

[30] Y. Yang, Y. Zhuang, F. Bi, H. Shi, and Y. Xie, "M-FCN: Effective fully convolutional network-based airplane detection framework," *IEEE Geosci. Remote Sens. Lett.*, vol. 14, no. 8, pp. 1293-1297, 2017.

[31] S. Wang, H. Wang, S. She, Y. Zhang, Q. Qiu, and Z. Xiao, "Swin-T-NFC CRFs: An encoder-decoder neural model for high-precision UAV positioning via point cloud super resolution and image semantic segmentation," *Computer Communications*, vol. 197, pp. 52-60, 2023.

[32] Y. Konforti, A. Shpigler, B. Lerner, and A. Bar-Hillel, "Sign: Statistical inference graphs based on probabilistic network activity interpretation," *IEEE Transactions on Pattern Analysis and Machine Intelligence*, vol. 45, no. 3, pp. 3783-3797, 2023.

[33] F. Coeurdoux, N. Dobigeon, and P. Chainais, "Plug-and-play split Gibbs sampler: Embedding deep generative priors in Bayesian inference," *IEEE Transactions on Image Processing*, vol. 33, pp. 3496-3507, 2024.

[34] L.-C. Chen, G. Papandreou, I. Kokkinos, K. Murphy, and A. L. Yuille, "DeepLab: Semantic image segmentation with deep convolutional nets, atrous convolution, and fully connected CRFs," *IEEE Trans. Pattern Anal. Mach. Intell.*, vol. 40, no. 4, pp. 834-848, 2018.

[35] P. Krähenbühl and V. Koltun, "Efficient inference in fully connected CRFs with Gaussian edge potentials," in *Proc. Neural Inf. Process. Syst.*, pp. 109-117, 2011.

[36] G. Lin, C. Shen, A. van den Hengel, and I. Reid, "Exploring context with deep structured models for semantic segmentation," *IEEE Transactions on Pattern Analysis and Machine Intelligence*, vol. 40, no. 6, pp. 1352-1366, 2018.

[37] Y. Bengio, "Learning deep architectures for AI," *Foundations and Trends in Machine Learning*, vol. 2, pp. 1-55, 2009.

[38] D. H. Ackley, G. E. Hinton, and T. J. Sejnowski, "A learning algorithm for Boltzmann machines," *Cognitive Science*, vol. 9, no. 1, pp. 147-169, 1985.

[39] L.-C. Chen, A. Schwing, A. Yuille, and R. Urtasun, "Learning deep structured models," in *Proc. Int. Conf. Mach. Learn. (ICML)*, ser. Proc. Mach. Learn. Res., F. Bach and D. Blei, Eds., vol. 37. Lille, France: PMLR, 07-09 Jul 2015, pp. 1785-1794.

[40] J. Xiang, G. Xu, C. Ma, and J. Hou, "End-to-end learning deep CRF models for multi-object tracking deep CRF models," *IEEE Trans. Circuits Syst. Video Technol.*, vol. 31, no. 1, pp. 275-288, 2021.

[41] S. Theodoridis, *Machine Learning: A Bayesian and Optimization Perspective*, 2nd ed. Academic Press, 2020.
BIOGRAPHIES

**Martina Pastorino** (Member, IEEE) received her M.Sc. degree from the University of Genoa, Italy and IMT Atlantique, Brest, France and her joint Ph.D. degree in Science and Technology for Electronic and Telecommunication Engineering from the University of Genoa, Italy and Inria d'Université Côte d'Azur, Sophia-Antipolis, France. She is an Assistant Professor at the University of Genoa. She received the Best Student Paper Award at IGARSS 2021, the Prix d'excellence d'Université Côte d'Azur in 2021, the IEEE GRSS Italian Chapter Best Master Thesis and Best PhD Thesis Prizes in 2022 and 2025, respectively, and the Best PhD Thesis Prize from the Ecole Doctorale STIC, Université Côte d'Azur in 2024. Her research activity is focused on the combination of stochastic models and deep learning techniques for remote sensing image analysis.

**Gabriele Moser** (Fellow, IEEE) received his Laurea (M.Sc. equivalent) degree and the Ph.D. degree from the University of Genoa, Italy. He is a Full Professor of Telecommunications at the University of Genoa. His research activity is focused on machine learning, data fusion, and image processing methodologies for remote sensing and energy applications. He was Associate Editor (2011-2015) and then Area Editor of Pattern Recognition Letters (2015-2018), has been Associate Editor of IEEE Geoscience and Remote Sensing Letters since 2008, and was Guest Co-Editor of the Sept. 2015 special issue of IEEE Geoscience and Remote Sensing Magazine. He served as Chair or Co-Chair of the Image Analysis and Data Fusion Technical Committee of the IEEE Geoscience and Remote Sensing Society (GRSS; 2013-2017), was Publication Co-Chair of IGARSS 2015, Technical Program Co-Chair of the IEEE GRSS EarthVision workshop at CVPR 2015, and Co-Organizer of the second edition of EarthVision at CVPR 2017. He received the Best Paper Award at the 2010 IEEE Workshop on Hyperspectral Image and Signal Processing (WHISPERS) and the Interactive Symposium Paper Award at IGARSS 2016. Together with Prof. Zerubia, he was editor of the book *Mathematical Models for Remote Sensing Image Processing* (Springer, 2019).

**Sebastiano B. Serpico** (Life Fellow, IEEE) received the Laurea (M.Sc. equivalent) degree and the Ph.D. degree from the University of Genoa, Italy. He is a Full Professor of Telecommunications at the Polytechnic School of the University of Genoa, where he is the coordinator of the research group on Signal Processing and Recognition Methods and Systems of the Department of Electrical, Electronic, Telecommunications Engineering, and Naval Architecture. His research interests are mainly in pattern recognition for remote sensing and biomedical image analysis. He was the Chairman of the Institute of Advanced Studies in Information and Communication Technologies (ISICT) from 2003 to 2019, and has been the Chairman of the Advanced Education School of the University of Genoa since 2022. He has been the project manager of numerous research projects and an evaluator of project proposals for various programs of the European Union, Italian Space Agency, Italian Ministry of Education and Research, etc. He received the 2019 IEEE GRSS Education Award, the Interactive Symposium Paper Award at IGARSS 2016, and the Best Paper Award at the IEEE WHISPERS 2010. He was an Associate Editor of the IEEE transactions on Geoscience and Remote Sensing from 2001 to 2022. He was Co-Chair of the IEEE International Geoscience and Remote Sensing Symposium in 2015 (Milan, Italy).

**Josiane Zerubia** (Life Fellow, IEEE) received her M.Sc. degree from the École nationale supérieure d'ingénieurs électriciens de Grenoble (ENSIEG), Grenoble, and her PhD from the University of Nice Sophia Antipolis, France. She is a Research Director at INRIA, Sophia Antipolis (DR Exceptional Class since 2023). She was head of the PASTIS laboratory (INRIA, 1995-1997) and of the ARIANA research group (INRIA/CNRS/University of Nice, 1998-2011), which worked on inverse problems in remote sensing and biological imaging. She headed the AYIN research group (INRIA-SAM; 2012-2016) dedicated to models of spatio-temporal structures for high-resolution image processing with focus on remote sensing and skincare imaging. She has been heading the AYANA exploratory research group since 2020. AYANA is an interdisciplinary project using knowledge in stochastic modeling, image processing, artificial intelligence, remote sensing and embedded electronics/computing. She was professor (PR1) at SUPAERO (ISAE) in Toulouse from 1999 to 2020. She received a Doctor Honoris Causa degree from the University of Szeged, Hungary, in 2020, and three times the Excellence Award from University of Nice. She is a Fellow of IEEE, EURASIP, and IAPR. She was IEEE SPS Distinguished Lecturer (2016-2017), Associate Editor (1998-2002) and then Area Editor of IEEE Transactions on Image Processing (2003-2006), Guest Co-Editor of a special issue of IEEE Transactions on Pattern Analysis and Machine Intelligence (2003), member of the Senior editorial board of IEEE Signal Processing Magazine (2018-2022), member of the editorial board of International Journal of Computer Vision (2004-2013), member-at-large of the Board of Governors of the IEEE SPS (2002-2004), member of the Awards Board of the IEEE SPS (2020-2022) and senior member of the Women in Signal Processing Committee of the IEEE SPS (2024-2025). Her main research interest is in image processing using probabilistic models. She also works on parameter estimation, statistical learning, optimization techniques, and artificial intelligence. Together with Prof. Moser, she was editor of the book *Mathematical Models for Remote Sensing Image Processing* (Springer, 2019).